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 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 2nd and 3rd columns specify
0020 %               the lower and upper limits on the (DC model) flow
0021 %               across the interface, respectively. Normally, the lower
0022 %               limit is negative, indicating a flow in the opposite
0023 %               direction.
0024 %
0025 %   The 'int2ext' callback also packages up results and stores them in
0026 %   the following output fields of results.if:
0027 %       P       - nif x 1, actual flow across each interface in MW
0028 %       mu.l    - nif x 1, shadow price on lower flow limit, ($/MW)
0029 %       mu.u    - nif x 1, shadow price on upper flow limit, ($/MW)
0030 %
0031 %   See also ADD_USERFCN, REMOVE_USERFCN, RUN_USERFCN, T_CASE30_USERFCNS.
0032 
0033 %   MATPOWER
0034 %   $Id: toggle_iflims.m 2430 2014-12-02 16:25:08Z ray $
0035 %   by Ray Zimmerman, PSERC Cornell
0036 %   Copyright (c) 2009-2010 by Power System Engineering Research Center (PSERC)
0037 %
0038 %   This file is part of MATPOWER.
0039 %   See http://www.pserc.cornell.edu/matpower/ for more info.
0040 %
0041 %   MATPOWER is free software: you can redistribute it and/or modify
0042 %   it under the terms of the GNU General Public License as published
0043 %   by the Free Software Foundation, either version 3 of the License,
0044 %   or (at your option) any later version.
0045 %
0046 %   MATPOWER is distributed in the hope that it will be useful,
0047 %   but WITHOUT ANY WARRANTY; without even the implied warranty of
0048 %   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
0049 %   GNU General Public License for more details.
0050 %
0051 %   You should have received a copy of the GNU General Public License
0052 %   along with MATPOWER. If not, see <http://www.gnu.org/licenses/>.
0053 %
0054 %   Additional permission under GNU GPL version 3 section 7
0055 %
0056 %   If you modify MATPOWER, or any covered work, to interface with
0057 %   other modules (such as MATLAB code and MEX-files) available in a
0058 %   MATLAB(R) or comparable environment containing parts covered
0059 %   under other licensing terms, the licensors of MATPOWER grant
0060 %   you additional permission to convey the resulting work.
0061 
0062 if strcmp(upper(on_off), 'ON')
0063     %% check for proper reserve inputs
0064     if ~isfield(mpc, 'if') || ~isstruct(mpc.if) || ...
0065             ~isfield(mpc.if, 'map') || ...
0066             ~isfield(mpc.if, 'lims')
0067         error('toggle_iflims: case must contain an ''if'' field, a struct defining ''map'' and ''lims''');
0068     end
0069     
0070     %% add callback functions
0071     %% note: assumes all necessary data included in 1st arg (mpc, om, results)
0072     %%       so, no additional explicit args are needed
0073     mpc = add_userfcn(mpc, 'ext2int', @userfcn_iflims_ext2int);
0074     mpc = add_userfcn(mpc, 'formulation', @userfcn_iflims_formulation);
0075     mpc = add_userfcn(mpc, 'int2ext', @userfcn_iflims_int2ext);
0076     mpc = add_userfcn(mpc, 'printpf', @userfcn_iflims_printpf);
0077     mpc = add_userfcn(mpc, 'savecase', @userfcn_iflims_savecase);
0078     mpc.userfcn.status.iflims = 1;
0079 elseif strcmp(upper(on_off), 'OFF')
0080     mpc = remove_userfcn(mpc, 'savecase', @userfcn_iflims_savecase);
0081     mpc = remove_userfcn(mpc, 'printpf', @userfcn_iflims_printpf);
0082     mpc = remove_userfcn(mpc, 'int2ext', @userfcn_iflims_int2ext);
0083     mpc = remove_userfcn(mpc, 'formulation', @userfcn_iflims_formulation);
0084     mpc = remove_userfcn(mpc, 'ext2int', @userfcn_iflims_ext2int);
0085     mpc.userfcn.status.iflims = 0;
0086 elseif strcmp(upper(on_off), 'STATUS')
0087     if isfield(mpc, 'userfcn') && isfield(mpc.userfcn, 'status') && ...
0088             isfield(mpc.userfcn.status, 'iflims')
0089         mpc = mpc.userfcn.status.iflims;
0090     else
0091         mpc = 0;
0092     end
0093 else
0094     error('toggle_iflims: 2nd argument must be ''on'', ''off'' or ''status''');
0095 end
0096 
0097 
0098 %%-----  ext2int  ------------------------------------------------------
0099 function mpc = userfcn_iflims_ext2int(mpc, args)
0100 %
0101 %   mpc = userfcn_iflims_ext2int(mpc, args)
0102 %
0103 %   This is the 'ext2int' stage userfcn callback that prepares the input
0104 %   data for the formulation stage. It expects to find an 'if' field in
0105 %   mpc as described above. The optional args are not currently used.
0106 
0107 %% initialize some things
0108 ifmap = mpc.if.map;
0109 o = mpc.order;
0110 nl0 = size(o.ext.branch, 1);    %% original number of branches
0111 nl = size(mpc.branch, 1);       %% number of on-line branches
0112 
0113 %% save if.map for external indexing
0114 mpc.order.ext.ifmap = ifmap;
0115 
0116 %%-----  convert stuff to internal indexing  -----
0117 e2i = zeros(nl0, 1);
0118 e2i(o.branch.status.on) = (1:nl)';  %% ext->int branch index mapping
0119 d = sign(ifmap(:, 2));
0120 br = abs(ifmap(:, 2));
0121 ifmap(:, 2) = d .* e2i(br);
0122 ifmap(ifmap(:, 2) == 0, :) = [];    %% delete branches that are out
0123 
0124 mpc.if.map = ifmap;
0125 
0126 
0127 %%-----  formulation  --------------------------------------------------
0128 function om = userfcn_iflims_formulation(om, args)
0129 %
0130 %   om = userfcn_iflims_formulation(om, args)
0131 %
0132 %   This is the 'formulation' stage userfcn callback that defines the
0133 %   user costs and constraints for interface flow limits. It expects to
0134 %   find an 'if' field in the mpc stored in om, as described above. The
0135 %   optional args are not currently used.
0136 
0137 %% define named indices into data matrices
0138 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0139     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0140     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0141 
0142 %% initialize some things
0143 mpc = get_mpc(om);
0144 [baseMVA, bus, branch] = deal(mpc.baseMVA, mpc.bus, mpc.branch);
0145 ifmap = mpc.if.map;
0146 iflims = mpc.if.lims;
0147 
0148 %% form B matrices for DC model
0149 [B, Bf, Pbusinj, Pfinj] = makeBdc(baseMVA, bus, branch);
0150 n = size(Bf, 2);                    %% dim of theta
0151 
0152 %% form constraints
0153 ifidx = unique(iflims(:, 1));   %% interface number list
0154 nifs = length(ifidx);           %% number of interfaces
0155 Aif = sparse(nifs, n);
0156 lif = zeros(nifs, 1);
0157 uif = zeros(nifs, 1);
0158 for k = 1:nifs
0159     %% extract branch indices
0160     br = ifmap(ifmap(:, 1) == ifidx(k), 2);
0161     if isempty(br)
0162         error('userfcn_iflims_formulation: interface %d has no in-service branches', k);
0163     end
0164     d = sign(br);
0165     br = abs(br);
0166     Ak = sparse(1, n);              %% Ak = sum( d(i) * Bf(i, :) )
0167     bk = 0;                         %% bk = sum( d(i) * Pfinj(i) )
0168     for i = 1:length(br)
0169         Ak = Ak + d(i) * Bf(br(i), :);
0170         bk = bk + d(i) * Pfinj(br(i));
0171     end
0172     Aif(k, :) = Ak;
0173     lif(k) = iflims(k, 2) / baseMVA - bk;
0174     uif(k) = iflims(k, 3) / baseMVA - bk;
0175 end
0176 
0177 %% add interface constraint
0178 om = add_constraints(om, 'iflims',  Aif, lif, uif, {'Va'});      %% nifs
0179 
0180 
0181 %%-----  int2ext  ------------------------------------------------------
0182 function results = userfcn_iflims_int2ext(results, args)
0183 %
0184 %   results = userfcn_iflims_int2ext(results, args)
0185 %
0186 %   This is the 'int2ext' stage userfcn callback that converts everything
0187 %   back to external indexing and packages up the results. It expects to
0188 %   find an 'if' field in the results struct as described for mpc above.
0189 %   It also expects the results to contain solved branch flows and linear
0190 %   constraints named 'iflims' which are used to populate output fields
0191 %   in results.if. The optional args are not currently used.
0192 
0193 %% define named indices into data matrices
0194 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0195     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0196     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0197 
0198 %% get internal ifmap
0199 ifmap = results.if.map;
0200 iflims = results.if.lims;
0201 
0202 %%-----  convert stuff back to external indexing  -----
0203 results.if.map = results.order.ext.ifmap;
0204 
0205 %%-----  results post-processing  -----
0206 ifidx = unique(iflims(:, 1));   %% interface number list
0207 nifs = length(ifidx);           %% number of interfaces
0208 results.if.P = zeros(nifs, 1);
0209 for k = 1:nifs
0210     %% extract branch indices
0211     br = ifmap(ifmap(:, 1) ==  ifidx(k), 2);
0212     d = sign(br);
0213     br = abs(br);
0214     results.if.P(k) = sum( d .* results.branch(br, PF) );
0215 end
0216 results.if.mu.l = results.lin.mu.l.iflims;
0217 results.if.mu.u = results.lin.mu.u.iflims;
0218 
0219 
0220 %%-----  printpf  ------------------------------------------------------
0221 function results = userfcn_iflims_printpf(results, fd, mpopt, args)
0222 %
0223 %   results = userfcn_iflims_printpf(results, fd, mpopt, args)
0224 %
0225 %   This is the 'printpf' stage userfcn callback that pretty-prints the
0226 %   results. It expects a results struct, a file descriptor and a MATPOWER
0227 %   options struct. The optional args are not currently used.
0228 
0229 %% define named indices into data matrices
0230 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0231     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0232     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0233 
0234 %%-----  print results  -----
0235 ptol = 1e-6;        %% tolerance for displaying shadow prices
0236 isOPF           = isfield(results, 'f') && ~isempty(results.f);
0237 SUPPRESS        = mpopt.out.suppress_detail;
0238 if SUPPRESS == -1
0239     if size(results.bus, 1) > 500
0240         SUPPRESS = 1;
0241     else
0242         SUPPRESS = 0;
0243     end
0244 end
0245 OUT_ALL         = mpopt.out.all;
0246 OUT_FORCE       = mpopt.out.force;
0247 OUT_BRANCH      = OUT_ALL == 1 || (OUT_ALL == -1 && ~SUPPRESS && mpopt.out.branch);
0248 
0249 if isOPF && OUT_BRANCH && (results.success || OUT_FORCE)
0250     iflims = results.if.lims;
0251     fprintf(fd, '\n================================================================================');
0252     fprintf(fd, '\n|     Interface Flow Limits                                                    |');
0253     fprintf(fd, '\n================================================================================');
0254     fprintf(fd, '\n Interface  Shadow Prc  Lower Lim      Flow      Upper Lim   Shadow Prc');
0255     fprintf(fd, '\n     #        ($/MW)       (MW)        (MW)        (MW)       ($/MW)   ');
0256     fprintf(fd, '\n----------  ----------  ----------  ----------  ----------  -----------');
0257     ifidx = unique(iflims(:, 1));   %% interface number list
0258     nifs = length(ifidx);           %% number of interfaces
0259     for k = 1:nifs
0260         fprintf(fd, '\n%6d ', iflims(k, 1));
0261         if results.if.mu.l(k) > ptol
0262             fprintf(fd, '%14.3f', results.if.mu.l(k));
0263         else
0264             fprintf(fd, '          -   ');
0265         end
0266         fprintf(fd, '%12.2f%12.2f%12.2f', iflims(k, 2), results.if.P(k), iflims(k, 3));
0267         if results.if.mu.u(k) > ptol
0268             fprintf(fd, '%13.3f', results.if.mu.u(k));
0269         else
0270             fprintf(fd, '         -     ');
0271         end
0272     end
0273     fprintf(fd, '\n');
0274 end
0275 
0276 
0277 %%-----  savecase  -----------------------------------------------------
0278 function mpc = userfcn_iflims_savecase(mpc, fd, prefix, args)
0279 %
0280 %   mpc = userfcn_iflims_savecase(mpc, fd, mpopt, args)
0281 %
0282 %   This is the 'savecase' stage userfcn callback that prints the M-file
0283 %   code to save the 'if' field in the case file. It expects a
0284 %   MATPOWER case struct (mpc), a file descriptor and variable prefix
0285 %   (usually 'mpc.'). The optional args are not currently used.
0286 
0287 ifmap = mpc.if.map;
0288 iflims = mpc.if.lims;
0289 
0290 fprintf(fd, '\n%%%%-----  Interface Flow Limit Data  -----%%%%\n');
0291 fprintf(fd, '%%%% interface<->branch map data\n');
0292 fprintf(fd, '%%\tifnum\tbranchidx (negative defines opposite direction)\n');
0293 fprintf(fd, '%sif.map = [\n', prefix);
0294 fprintf(fd, '\t%d\t%d;\n', ifmap');
0295 fprintf(fd, '];\n');
0296 
0297 fprintf(fd, '\n%%%% interface flow limit data (based on DC model)\n');
0298 fprintf(fd, '%%%% (lower limit should be negative for opposite direction)\n');
0299 fprintf(fd, '%%\tifnum\tlower\tupper\n');
0300 fprintf(fd, '%sif.lims = [\n', prefix);
0301 fprintf(fd, '\t%d\t%g\t%g;\n', iflims');
0302 fprintf(fd, '];\n');
0303 
0304 %% save output fields for solved case
0305 if isfield(mpc.if, 'P')
0306     if exist('serialize', 'file') == 2
0307         fprintf(fd, '\n%%%% solved values\n');
0308         fprintf(fd, '%sif.P = %s\n', prefix, serialize(mpc.if.P));
0309         fprintf(fd, '%sif.mu.l = %s\n', prefix, serialize(mpc.if.mu.l));
0310         fprintf(fd, '%sif.mu.u = %s\n', prefix, serialize(mpc.if.mu.u));
0311     else
0312         url = 'http://www.mathworks.com/matlabcentral/fileexchange/12063';
0313         warning('MATPOWER:serialize', ...
0314             '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);
0315     end
0316 end