toggle_iflims

TOGGLE_IFLIMS Enable, disable or check status of set of interface flow limits.
   MPC = TOGGLE_IFLIMS(MPC, 'on')
   MPC = TOGGLE_IFLIMS(MPC, 'off')
   T_F = TOGGLE_IFLIMS(MPC, 'status')

   Enables, disables or checks the status of a set of OPF userfcn
   callbacks to implement interface flow limits based on a DC flow model.

   These callbacks expect to find an 'if' field in the input MPC, where
   MPC.if is a struct with the following fields:
       map     n x 2, defines each interface in terms of a set of 
               branch indices and directions. Interface I is defined
               by the set of rows whose 1st col is equal to I. The
               2nd column is a branch index multiplied by 1 or -1
               respectively for lines whose orientation is the same
               as or opposite to that of the interface.
       lims    nif x 3, defines the DC model flow limits in MW
               for specified interfaces. The first column is the index
               of the interface, the 2nd and 3rd columns specify
               the lower and upper limits on the (DC model) flow
               across the interface, respectively. Normally, the lower
               limit is negative, indicating a flow in the opposite
               direction.

   The 'int2ext' callback also packages up results and stores them in
   the following output fields of results.if:
       P       - nif x 1, actual flow across each interface in MW
       mu.l    - nif x 1, shadow price on lower flow limit, ($/MW)
       mu.u    - nif x 1, shadow price on upper flow limit, ($/MW)

   See also ADD_USERFCN, REMOVE_USERFCN, RUN_USERFCN, T_CASE30_USERFCNS.

0001 function mpc = toggle_iflims(mpc, on_off)
0002 %TOGGLE_IFLIMS Enable, disable or check status of set of interface flow limits.
0003 %   MPC = TOGGLE_IFLIMS(MPC, 'on')
0004 %   MPC = TOGGLE_IFLIMS(MPC, 'off')
0005 %   T_F = TOGGLE_IFLIMS(MPC, 'status')
0006 %
0007 %   Enables, disables or checks the status of a set of OPF userfcn
0008 %   callbacks to implement interface flow limits based on a DC flow model.
0009 %
0010 %   These callbacks expect to find an 'if' field in the input MPC, where
0011 %   MPC.if is a struct with the following fields:
0012 %       map     n x 2, defines each interface in terms of a set of
0013 %               branch indices and directions. Interface I is defined
0014 %               by the set of rows whose 1st col is equal to I. The
0015 %               2nd column is a branch index multiplied by 1 or -1
0016 %               respectively for lines whose orientation is the same
0017 %               as or opposite to that of the interface.
0018 %       lims    nif x 3, defines the DC model flow limits in MW
0019 %               for specified interfaces. The first column is the index
0020 %               of the interface, the 2nd and 3rd columns specify
0021 %               the lower and upper limits on the (DC model) flow
0022 %               across the interface, respectively. Normally, the lower
0023 %               limit is negative, indicating a flow in the opposite
0024 %               direction.
0025 %
0026 %   The 'int2ext' callback also packages up results and stores them in
0027 %   the following output fields of results.if:
0028 %       P       - nif x 1, actual flow across each interface in MW
0029 %       mu.l    - nif x 1, shadow price on lower flow limit, ($/MW)
0030 %       mu.u    - nif x 1, shadow price on upper flow limit, ($/MW)
0031 %
0032 %   See also ADD_USERFCN, REMOVE_USERFCN, RUN_USERFCN, T_CASE30_USERFCNS.
0033 
0034 %   MATPOWER
0035 %   Copyright (c) 2009-2016, Power Systems Engineering Research Center (PSERC)
0036 %   by Ray Zimmerman, PSERC Cornell
0037 %
0038 %   This file is part of MATPOWER.
0039 %   Covered by the 3-clause BSD License (see LICENSE file for details).
0040 %   See https://matpower.org for more info.
0041 
0042 if strcmp(upper(on_off), 'ON')
0043     %% check for proper reserve inputs
0044     if ~isfield(mpc, 'if') || ~isstruct(mpc.if) || ...
0045             ~isfield(mpc.if, 'map') || ...
0046             ~isfield(mpc.if, 'lims')
0047         error('toggle_iflims: case must contain an ''if'' field, a struct defining ''map'' and ''lims''');
0048     end
0049     
0050     %% add callback functions
0051     %% note: assumes all necessary data included in 1st arg (mpc, om, results)
0052     %%       so, no additional explicit args are needed
0053     mpc = add_userfcn(mpc, 'ext2int', @userfcn_iflims_ext2int);
0054     mpc = add_userfcn(mpc, 'formulation', @userfcn_iflims_formulation);
0055     mpc = add_userfcn(mpc, 'int2ext', @userfcn_iflims_int2ext);
0056     mpc = add_userfcn(mpc, 'printpf', @userfcn_iflims_printpf);
0057     mpc = add_userfcn(mpc, 'savecase', @userfcn_iflims_savecase);
0058     mpc.userfcn.status.iflims = 1;
0059 elseif strcmp(upper(on_off), 'OFF')
0060     mpc = remove_userfcn(mpc, 'savecase', @userfcn_iflims_savecase);
0061     mpc = remove_userfcn(mpc, 'printpf', @userfcn_iflims_printpf);
0062     mpc = remove_userfcn(mpc, 'int2ext', @userfcn_iflims_int2ext);
0063     mpc = remove_userfcn(mpc, 'formulation', @userfcn_iflims_formulation);
0064     mpc = remove_userfcn(mpc, 'ext2int', @userfcn_iflims_ext2int);
0065     mpc.userfcn.status.iflims = 0;
0066 elseif strcmp(upper(on_off), 'STATUS')
0067     if isfield(mpc, 'userfcn') && isfield(mpc.userfcn, 'status') && ...
0068             isfield(mpc.userfcn.status, 'iflims')
0069         mpc = mpc.userfcn.status.iflims;
0070     else
0071         mpc = 0;
0072     end
0073 else
0074     error('toggle_iflims: 2nd argument must be ''on'', ''off'' or ''status''');
0075 end
0076 
0077 
0078 %%-----  ext2int  ------------------------------------------------------
0079 function mpc = userfcn_iflims_ext2int(mpc, mpopt, args)
0080 %
0081 %   mpc = userfcn_iflims_ext2int(mpc, mpopt, args)
0082 %
0083 %   This is the 'ext2int' stage userfcn callback that prepares the input
0084 %   data for the formulation stage. It expects to find an 'if' field in
0085 %   mpc as described above. The optional args are not currently used.
0086 
0087 %% initialize some things
0088 ifmap = mpc.if.map;
0089 o = mpc.order;
0090 nl0 = size(o.ext.branch, 1);    %% original number of branches
0091 nl = size(mpc.branch, 1);       %% number of on-line branches
0092 
0093 %% save if.map for external indexing
0094 mpc.order.ext.ifmap = ifmap;
0095 
0096 %%-----  convert stuff to internal indexing  -----
0097 e2i = zeros(nl0, 1);
0098 e2i(o.branch.status.on) = (1:nl)';  %% ext->int branch index mapping
0099 d = sign(ifmap(:, 2));
0100 br = abs(ifmap(:, 2));
0101 ifmap(:, 2) = d .* e2i(br);
0102 ifmap(ifmap(:, 2) == 0, :) = [];    %% delete branches that are out
0103 
0104 mpc.if.map = ifmap;
0105 
0106 
0107 %%-----  formulation  --------------------------------------------------
0108 function om = userfcn_iflims_formulation(om, mpopt, args)
0109 %
0110 %   om = userfcn_iflims_formulation(om, mpopt, args)
0111 %
0112 %   This is the 'formulation' stage userfcn callback that defines the
0113 %   user costs and constraints for interface flow limits. It expects to
0114 %   find an 'if' field in the mpc stored in om, as described above. The
0115 %   optional args are not currently used.
0116 
0117 %% define named indices into data matrices
0118 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0119     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0120     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0121 
0122 %% initialize some things
0123 mpc = om.get_mpc();
0124 [baseMVA, bus, branch] = deal(mpc.baseMVA, mpc.bus, mpc.branch);
0125 ifmap = mpc.if.map;
0126 iflims = mpc.if.lims;
0127 
0128 %% form B matrices for DC model
0129 [B, Bf, Pbusinj, Pfinj] = makeBdc(baseMVA, bus, branch);
0130 n = size(Bf, 2);                    %% dim of theta
0131 
0132 %% form constraints
0133 ifidx = unique(iflims(:, 1));   %% interface number list
0134 nifs = length(ifidx);           %% number of interfaces
0135 Aif = sparse(nifs, n);
0136 lif = zeros(nifs, 1);
0137 uif = zeros(nifs, 1);
0138 for k = 1:nifs
0139     %% extract branch indices
0140     br = ifmap(ifmap(:, 1) == ifidx(k), 2);
0141     if isempty(br)
0142         error('userfcn_iflims_formulation: interface %d has no in-service branches', k);
0143     end
0144     d = sign(br);
0145     br = abs(br);
0146     Ak = sparse(1, n);              %% Ak = sum( d(i) * Bf(i, :) )
0147     bk = 0;                         %% bk = sum( d(i) * Pfinj(i) )
0148     for i = 1:length(br)
0149         Ak = Ak + d(i) * Bf(br(i), :);
0150         bk = bk + d(i) * Pfinj(br(i));
0151     end
0152     Aif(k, :) = Ak;
0153     lif(k) = iflims(k, 2) / baseMVA - bk;
0154     uif(k) = iflims(k, 3) / baseMVA - bk;
0155 end
0156 
0157 %% add interface constraint
0158 om.add_lin_constraint('iflims',  Aif, lif, uif, {'Va'});   %% nifs
0159 
0160 
0161 %%-----  int2ext  ------------------------------------------------------
0162 function results = userfcn_iflims_int2ext(results, mpopt, args)
0163 %
0164 %   results = userfcn_iflims_int2ext(results, mpopt, args)
0165 %
0166 %   This is the 'int2ext' stage userfcn callback that converts everything
0167 %   back to external indexing and packages up the results. It expects to
0168 %   find an 'if' field in the results struct as described for mpc above.
0169 %   It also expects the results to contain solved branch flows and linear
0170 %   constraints named 'iflims' which are used to populate output fields
0171 %   in results.if. The optional args are not currently used.
0172 
0173 %% define named indices into data matrices
0174 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0175     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0176     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0177 
0178 %% get internal ifmap
0179 ifmap = results.if.map;
0180 iflims = results.if.lims;
0181 
0182 %%-----  convert stuff back to external indexing  -----
0183 results.if.map = results.order.ext.ifmap;
0184 
0185 %%-----  results post-processing  -----
0186 ifidx = unique(iflims(:, 1));   %% interface number list
0187 nifs = length(ifidx);           %% number of interfaces
0188 results.if.P = zeros(nifs, 1);
0189 for k = 1:nifs
0190     %% extract branch indices
0191     br = ifmap(ifmap(:, 1) ==  ifidx(k), 2);
0192     d = sign(br);
0193     br = abs(br);
0194     results.if.P(k) = sum( d .* results.branch(br, PF) );
0195 end
0196 results.if.mu.l = results.lin.mu.l.iflims / results.baseMVA;
0197 results.if.mu.u = results.lin.mu.u.iflims / results.baseMVA;
0198 
0199 
0200 %%-----  printpf  ------------------------------------------------------
0201 function results = userfcn_iflims_printpf(results, fd, mpopt, args)
0202 %
0203 %   results = userfcn_iflims_printpf(results, fd, mpopt, args)
0204 %
0205 %   This is the 'printpf' stage userfcn callback that pretty-prints the
0206 %   results. It expects a results struct, a file descriptor and a MATPOWER
0207 %   options struct. The optional args are not currently used.
0208 
0209 %% define named indices into data matrices
0210 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0211     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0212     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0213 
0214 %%-----  print results  -----
0215 ptol = 1e-6;        %% tolerance for displaying shadow prices
0216 isOPF           = isfield(results, 'f') && ~isempty(results.f);
0217 SUPPRESS        = mpopt.out.suppress_detail;
0218 if SUPPRESS == -1
0219     if size(results.bus, 1) > 500
0220         SUPPRESS = 1;
0221     else
0222         SUPPRESS = 0;
0223     end
0224 end
0225 OUT_ALL         = mpopt.out.all;
0226 OUT_FORCE       = mpopt.out.force;
0227 OUT_BRANCH      = OUT_ALL == 1 || (OUT_ALL == -1 && ~SUPPRESS && mpopt.out.branch);
0228 
0229 if isOPF && OUT_BRANCH && (results.success || OUT_FORCE)
0230     iflims = results.if.lims;
0231     fprintf(fd, '\n================================================================================');
0232     fprintf(fd, '\n|     Interface Flow Limits                                                    |');
0233     fprintf(fd, '\n================================================================================');
0234     fprintf(fd, '\n Interface  Shadow Prc  Lower Lim      Flow      Upper Lim   Shadow Prc');
0235     fprintf(fd, '\n     #        ($/MW)       (MW)        (MW)        (MW)       ($/MW)   ');
0236     fprintf(fd, '\n----------  ----------  ----------  ----------  ----------  -----------');
0237     ifidx = unique(iflims(:, 1));   %% interface number list
0238     nifs = length(ifidx);           %% number of interfaces
0239     for k = 1:nifs
0240         fprintf(fd, '\n%6d ', iflims(k, 1));
0241         if results.if.mu.l(k) > ptol
0242             fprintf(fd, '%14.3f', results.if.mu.l(k));
0243         else
0244             fprintf(fd, '          -   ');
0245         end
0246         fprintf(fd, '%12.2f%12.2f%12.2f', iflims(k, 2), results.if.P(k), iflims(k, 3));
0247         if results.if.mu.u(k) > ptol
0248             fprintf(fd, '%13.3f', results.if.mu.u(k));
0249         else
0250             fprintf(fd, '         -     ');
0251         end
0252     end
0253     fprintf(fd, '\n');
0254 end
0255 
0256 
0257 %%-----  savecase  -----------------------------------------------------
0258 function mpc = userfcn_iflims_savecase(mpc, fd, prefix, args)
0259 %
0260 %   mpc = userfcn_iflims_savecase(mpc, fd, prefix, args)
0261 %
0262 %   This is the 'savecase' stage userfcn callback that prints the M-file
0263 %   code to save the 'if' field in the case file. It expects a
0264 %   MATPOWER case struct (mpc), a file descriptor and variable prefix
0265 %   (usually 'mpc.'). The optional args are not currently used.
0266 
0267 ifmap = mpc.if.map;
0268 iflims = mpc.if.lims;
0269 
0270 fprintf(fd, '\n%%%%-----  Interface Flow Limit Data  -----%%%%\n');
0271 fprintf(fd, '%%%% interface<->branch map data\n');
0272 fprintf(fd, '%%\tifnum\tbranchidx (negative defines opposite direction)\n');
0273 fprintf(fd, '%sif.map = [\n', prefix);
0274 fprintf(fd, '\t%d\t%d;\n', ifmap');
0275 fprintf(fd, '];\n');
0276 
0277 fprintf(fd, '\n%%%% interface flow limit data (based on DC model)\n');
0278 fprintf(fd, '%%%% (lower limit should be negative for opposite direction)\n');
0279 fprintf(fd, '%%\tifnum\tlower\tupper\n');
0280 fprintf(fd, '%sif.lims = [\n', prefix);
0281 fprintf(fd, '\t%d\t%g\t%g;\n', iflims');
0282 fprintf(fd, '];\n');
0283 
0284 %% save output fields for solved case
0285 if isfield(mpc.if, 'P')
0286     if exist('serialize', 'file') == 2
0287         fprintf(fd, '\n%%%% solved values\n');
0288         fprintf(fd, '%sif.P = %s\n', prefix, serialize(mpc.if.P));
0289         fprintf(fd, '%sif.mu.l = %s\n', prefix, serialize(mpc.if.mu.l));
0290         fprintf(fd, '%sif.mu.u = %s\n', prefix, serialize(mpc.if.mu.u));
0291     else
0292         url = 'https://www.mathworks.com/matlabcentral/fileexchange/12063-serialize';
0293         warning('MATPOWER:serialize', ...
0294             'userfcn_iflims_savecase: Cannot save the ''iflims'' output fields without the ''serialize'' function, which is available as a free download from:\n<%s>\n\n', url);
0295     end
0296 end