toggle_softlims

TOGGLE_SOFTLIMS Relax DC optimal power flow branch limits.
   MPC = TOGGLE_SOFTLIMS(MPC, 'on')
   MPC = TOGGLE_SOFTLIMS(MPC, 'off')
   T_F = TOGGLE_SOFTLIMS(MPC, 'status')

   Enables, disables or checks the status of a set of OPF userfcn
   callbacks to implement relaxed branch flow limits a DC optimal flow model.

   These callbacks expect to find a 'softlims' field in the input MPC,
   where MPC.softlims is a struct with the following fields:
       idx     (optional) n x 1, index vector for branches whos flow
               limits are to be relaxed, default is to use all on-line
               branches with non-zero limits specified in RATE_A
       cost    (optional) n x 1, linear marginal cost per MW of exceeding
               RATE_A. Can optionally be a scalar, in which case it is
               applied to all soft limits. Default if not specified is
               $1000/MW.

   The 'int2ext' callback also packages up results and stores them in
   the following output fields of results.softlims:
       overload - nl x 1, amount of overload of each line in MW
       ovl_cost - nl x 1, total cost of overload in $/hr

   The shadow prices on the soft limit constraints are also returned in the
   MU_SF and MU_ST columns of the branch matrix.
   Note: These shadow prices are equal to the corresponding hard limit
       shadow prices when the soft limits are not violated. When violated,
       the shadow price on a soft limit constraint is equal to the
       user-specified soft limit violation cost.

   See also ADD_USERFCN, REMOVE_USERFCN, RUN_USERFCN, T_CASE30_USERFCNS.

0001 function mpc = toggle_softlims(mpc, on_off)
0002 %TOGGLE_SOFTLIMS Relax DC optimal power flow branch limits.
0003 %   MPC = TOGGLE_SOFTLIMS(MPC, 'on')
0004 %   MPC = TOGGLE_SOFTLIMS(MPC, 'off')
0005 %   T_F = TOGGLE_SOFTLIMS(MPC, 'status')
0006 %
0007 %   Enables, disables or checks the status of a set of OPF userfcn
0008 %   callbacks to implement relaxed branch flow limits a DC optimal flow model.
0009 %
0010 %   These callbacks expect to find a 'softlims' field in the input MPC,
0011 %   where MPC.softlims is a struct with the following fields:
0012 %       idx     (optional) n x 1, index vector for branches whos flow
0013 %               limits are to be relaxed, default is to use all on-line
0014 %               branches with non-zero limits specified in RATE_A
0015 %       cost    (optional) n x 1, linear marginal cost per MW of exceeding
0016 %               RATE_A. Can optionally be a scalar, in which case it is
0017 %               applied to all soft limits. Default if not specified is
0018 %               $1000/MW.
0019 %
0020 %   The 'int2ext' callback also packages up results and stores them in
0021 %   the following output fields of results.softlims:
0022 %       overload - nl x 1, amount of overload of each line in MW
0023 %       ovl_cost - nl x 1, total cost of overload in $/hr
0024 %
0025 %   The shadow prices on the soft limit constraints are also returned in the
0026 %   MU_SF and MU_ST columns of the branch matrix.
0027 %   Note: These shadow prices are equal to the corresponding hard limit
0028 %       shadow prices when the soft limits are not violated. When violated,
0029 %       the shadow price on a soft limit constraint is equal to the
0030 %       user-specified soft limit violation cost.
0031 %
0032 %   See also ADD_USERFCN, REMOVE_USERFCN, RUN_USERFCN, T_CASE30_USERFCNS.
0033 
0034 %   To do for future versions:
0035 %       - make softlims input field optional, convert all lines if missing
0036 %       Inputs:
0037 %       cost    n x 3, linear marginal cost per MW of exceeding each of
0038 %               RATE_A, RATE_B and RATE_C. Columns 2 and 3 are optional.
0039 %       brkpts  n x npts, allow to specify arbitrary breakpoints at which
0040 %               cost increases, defined as percentages above RATE_A.
0041 %       base_flow   n x 1, arbitrary baseline (other than RATE_A)
0042 
0043 %   MATPOWER
0044 %   Copyright (c) 2009-2015 by Power System Engineering Research Center (PSERC)
0045 %   by Ray Zimmerman, PSERC Cornell
0046 %
0047 %   $Id: toggle_softlims.m 2644 2015-03-11 19:34:22Z ray $
0048 %
0049 %   This file is part of MATPOWER.
0050 %   Covered by the 3-clause BSD License (see LICENSE file for details).
0051 %   See http://www.pserc.cornell.edu/matpower/ for more info.
0052 
0053 if strcmp(upper(on_off), 'ON')
0054     %% check for proper softlims inputs
0055     %% (inputs are optional, defaults handled in ext2int callback)
0056     
0057     %% add callback functions
0058     %% note: assumes all necessary data included in 1st arg (mpc, om, results)
0059     %%       so, no additional explicit args are needed
0060     mpc = add_userfcn(mpc, 'ext2int', @userfcn_softlims_ext2int);
0061     mpc = add_userfcn(mpc, 'formulation', @userfcn_softlims_formulation);
0062     mpc = add_userfcn(mpc, 'int2ext', @userfcn_softlims_int2ext);
0063     mpc = add_userfcn(mpc, 'printpf', @userfcn_softlims_printpf);
0064     mpc = add_userfcn(mpc, 'savecase', @userfcn_softlims_savecase);
0065     mpc.userfcn.status.softlims = 1;
0066 elseif strcmp(upper(on_off), 'OFF')
0067     mpc = remove_userfcn(mpc, 'savecase', @userfcn_softlims_savecase);
0068     mpc = remove_userfcn(mpc, 'printpf', @userfcn_softlims_printpf);
0069     mpc = remove_userfcn(mpc, 'int2ext', @userfcn_softlims_int2ext);
0070     mpc = remove_userfcn(mpc, 'formulation', @userfcn_softlims_formulation);
0071     mpc = remove_userfcn(mpc, 'ext2int', @userfcn_softlims_ext2int);
0072     mpc.userfcn.status.softlims = 0;
0073 elseif strcmp(upper(on_off), 'STATUS')
0074     if isfield(mpc, 'userfcn') && isfield(mpc.userfcn, 'status') && ...
0075             isfield(mpc.userfcn.status, 'softlims')
0076         mpc = mpc.userfcn.status.softlims;
0077     else
0078         mpc = 0;
0079     end
0080 else
0081     error('toggle_softlims: 2nd argument must be ''on'', ''off'' or ''status''');
0082 end
0083 
0084 
0085 %%-----  ext2int  ------------------------------------------------------
0086 function mpc = userfcn_softlims_ext2int(mpc, args)
0087 %
0088 %   mpc = userfcn_softlims_ext2int(mpc, args)
0089 %
0090 %   This is the 'ext2int' stage userfcn callback that prepares the input
0091 %   data for the formulation stage. It expects to find a 'softlims' field in
0092 %   mpc as described above. The optional args are not currently used.
0093 
0094 %% define named indices into data matrices
0095 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0096     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0097     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0098 
0099 %% check for proper softlims inputs
0100 default_cost = 1000;    %% used if cost is not specified
0101 if isfield(mpc, 'softlims')
0102     if ~isfield(mpc.softlims, 'idx')
0103         mpc.softlims.idx = [];
0104     end
0105     if ~isfield(mpc.softlims, 'cost')
0106         mpc.softlims.cost = default_cost;
0107     end
0108 else
0109     mpc.softlims.idx = [];
0110     mpc.softlims.cost = default_cost;
0111 end
0112 
0113 %% initialize some things
0114 s = mpc.softlims;
0115 o = mpc.order;
0116 nl0 = size(o.ext.branch, 1);    %% original number of branches
0117 nl = size(mpc.branch, 1);       %% number of on-line branches
0118 
0119 %% save softlims struct for external indexing
0120 mpc.order.ext.softlims = s;
0121 
0122 %%-----  convert stuff to internal indexing  -----
0123 s = softlims_defaults(s, o.ext.branch);     %% get defaults
0124 e2i = zeros(nl0, 1);
0125 e2i(o.branch.status.on) = (1:nl)';  %% ext->int branch index mapping
0126 s.idx = e2i(s.idx);
0127 k = find(s.idx == 0);   %% find idxes corresponding to off-line branches
0128 s.idx(k) = [];          %% delete them
0129 s.cost(k, :) = [];
0130 k = find(mpc.branch(s.idx, RATE_A) <= 0);   %% find branches w/o limits
0131 s.idx(k) = [];          %% delete them
0132 s.cost(k, :) = [];
0133 
0134 %%-----  remove hard limits on branches with soft limits  -----
0135 s.Pfmax = mpc.branch(s.idx, RATE_A) / mpc.baseMVA;  %% save limit first
0136 mpc.branch(s.idx, RATE_A) = 0;                      %% then remove it
0137 
0138 mpc.softlims = s;
0139 mpc.order.int.softlims = s;
0140 
0141 
0142 %%-----  formulation  --------------------------------------------------
0143 function om = userfcn_softlims_formulation(om, args)
0144 %
0145 %   om = userfcn_softlims_formulation(om, args)
0146 %
0147 %   This is the 'formulation' stage userfcn callback that defines the
0148 %   user costs and constraints for interface flow limits. It expects to
0149 %   find a 'softlims' field in the mpc stored in om, as described above. The
0150 %   optional args are not currently used.
0151 
0152 %% define named indices into data matrices
0153 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0154     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0155     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0156 
0157 %% initialize some things
0158 mpc = get_mpc(om);
0159 [baseMVA, bus, branch] = deal(mpc.baseMVA, mpc.bus, mpc.branch);
0160 s = mpc.softlims;
0161 
0162 %% form B matrices for DC model
0163 [B, Bf, Pbusinj, Pfinj] = makeBdc(baseMVA, bus, branch);
0164 n = size(Bf, 2);        %% dim of theta
0165 
0166 %% form constraints (flv is flow limit violation variable)
0167 %%   -Bf * Va - flv <=  Pfinj + Pfmax
0168 %%    Bf * Va - flv <= -Pfinj + Pfmax
0169 ns = length(s.idx);     %% number of soft limits
0170 I = speye(ns, ns);
0171 As = [-Bf(s.idx, :) -I; Bf(s.idx, :) -I];
0172 ls = -Inf(2*ns, 1);
0173 us = [   Pfinj(s.idx) + s.Pfmax;
0174         -Pfinj(s.idx) + s.Pfmax ];
0175 
0176 %% costs on flv variable
0177 Cw = s.cost(:, 1) * mpc.baseMVA;
0178 
0179 %% add vars, costs, constraints
0180 om = add_vars(om, 'flv', ns, zeros(ns, 1), zeros(ns, 1), Inf(ns, 1));
0181 om = add_costs(om, 'vc', struct('N', I, 'Cw', Cw), {'flv'});
0182 om = add_constraints(om, 'softlims',  As, ls, us, {'Va', 'flv'});   %% 2*ns
0183 
0184 
0185 %%-----  int2ext  ------------------------------------------------------
0186 function results = userfcn_softlims_int2ext(results, args)
0187 %
0188 %   results = userfcn_softlims_int2ext(results, args)
0189 %
0190 %   This is the 'int2ext' stage userfcn callback that converts everything
0191 %   back to external indexing and packages up the results. It expects to
0192 %   find a 'softlims' field in the results struct as described for mpc above.
0193 %   It also expects the results to contain solved branch flows and linear
0194 %   constraints named 'softlims' which are used to populate output fields
0195 %   in results.softlims. The optional args are not currently used.
0196 
0197 %% define named indices into data matrices
0198 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0199     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0200     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0201 
0202 %% get internal softlims struct
0203 s = results.softlims;
0204 
0205 %%-----  convert stuff back to external indexing  -----
0206 o = results.order;
0207 nl0 = size(o.ext.branch, 1);    %% original number of branches
0208 nl = size(results.branch, 1);   %% number of on-line branches
0209 o.branch.status.on;
0210 results.softlims = results.order.ext.softlims;
0211 
0212 %%-----  restore hard limits  -----
0213 results.branch(s.idx, RATE_A) = s.Pfmax * results.baseMVA;
0214 
0215 %%-----  results post-processing  -----
0216 %% get shadow prices
0217 n = size(results.lin.mu.u.softlims, 1) / 2;
0218 results.branch(s.idx, MU_ST) = results.lin.mu.u.softlims(1:n) / results.baseMVA;
0219 results.branch(s.idx, MU_SF) = results.lin.mu.u.softlims(n+1:end) / results.baseMVA;
0220 %% get overloads and overload costs
0221 results.softlims.overload = zeros(nl0, 1);
0222 k = find(results.branch(:, RATE_A) & ...
0223          abs(results.branch(:, PF)) > results.branch(:, RATE_A) );
0224 results.softlims.overload(o.branch.status.on(k)) = ...
0225         abs(results.branch(k, PF)) - results.branch(k, RATE_A);
0226 results.softlims.ovl_cost = zeros(size(results.softlims.overload));
0227 results.softlims.ovl_cost(o.branch.status.on(s.idx)) = ...
0228     results.softlims.overload(o.branch.status.on(s.idx)) .* s.cost(:, 1);
0229 
0230 
0231 %%-----  printpf  ------------------------------------------------------
0232 function results = userfcn_softlims_printpf(results, fd, mpopt, args)
0233 %
0234 %   results = userfcn_softlims_printpf(results, fd, mpopt, args)
0235 %
0236 %   This is the 'printpf' stage userfcn callback that pretty-prints the
0237 %   results. It expects a results struct, a file descriptor and a MATPOWER
0238 %   options struct. The optional args are not currently used.
0239 
0240 %% define named indices into data matrices
0241 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0242     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0243     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0244 
0245 %%-----  print results  -----
0246 ptol = 1e-6;        %% tolerance for displaying shadow prices
0247 isOPF           = isfield(results, 'f') && ~isempty(results.f);
0248 SUPPRESS        = mpopt.out.suppress_detail;
0249 if SUPPRESS == -1
0250     if size(results.bus, 1) > 500
0251         SUPPRESS = 1;
0252     else
0253         SUPPRESS = 0;
0254     end
0255 end
0256 OUT_ALL         = mpopt.out.all;
0257 OUT_FORCE       = mpopt.out.force;
0258 OUT_BRANCH      = OUT_ALL == 1 || (OUT_ALL == -1 && ~SUPPRESS && mpopt.out.branch);
0259 
0260 if isOPF && OUT_BRANCH && (results.success || OUT_FORCE)
0261     s = softlims_defaults(results.softlims, results.branch);   %% get defaults
0262     k = find(s.overload(s.idx) | sum(results.branch(s.idx, MU_SF:MU_ST), 2) > ptol);
0263 
0264     fprintf(fd, '\n================================================================================');
0265     fprintf(fd, '\n|     Soft Flow Limits                                                         |');
0266     fprintf(fd, '\n================================================================================');
0267     fprintf(fd, '\nBrnch   From   To     Flow      Limit    Overload     mu');
0268     fprintf(fd, '\n  #     Bus    Bus    (MW)      (MW)       (MW)     ($/MW)');
0269     fprintf(fd, '\n-----  -----  -----  --------  --------  --------  ---------');
0270     fprintf(fd, '\n%4d%7d%7d%10.2f%10.2f%10.2f%11.3f', ...
0271             [   s.idx(k), results.branch(s.idx(k), [F_BUS, T_BUS]), ...
0272                 results.branch(s.idx(k), [PF, RATE_A]), ...
0273                 s.overload(s.idx(k)), ...
0274                 sum(results.branch(s.idx(k), [MU_SF:MU_ST]), 2) ...
0275             ]');
0276     fprintf(fd, '\n                                         --------');
0277     fprintf(fd, '\n                                Total:%10.2f', ...
0278             sum(s.overload(s.idx(k))));
0279     fprintf(fd, '\n');
0280 end
0281 
0282 
0283 %%-----  savecase  -----------------------------------------------------
0284 function mpc = userfcn_softlims_savecase(mpc, fd, prefix, args)
0285 %
0286 %   mpc = userfcn_softlims_savecase(mpc, fd, mpopt, args)
0287 %
0288 %   This is the 'savecase' stage userfcn callback that prints the M-file
0289 %   code to save the 'softlims' field in the case file. It expects a
0290 %   MATPOWER case struct (mpc), a file descriptor and variable prefix
0291 %   (usually 'mpc.'). The optional args are not currently used.
0292 
0293 if isfield(mpc, 'softlims')
0294     s = mpc.softlims;
0295 
0296     fprintf(fd, '\n%%%%-----  Soft Flow Limit Data  -----%%%%\n');
0297 
0298     if isfield(s, 'idx')
0299         fprintf(fd, '%%%% branch indexes\n');
0300         fprintf(fd, '%%\tbranchidx\n');
0301         if isempty(s.idx)
0302             fprintf(fd, '%ssoftlims.idx = [];\n\n', prefix);
0303         else
0304             fprintf(fd, '%ssoftlims.idx = [\n', prefix);
0305             fprintf(fd, '\t%d;\n', s.idx);
0306             fprintf(fd, '];\n\n');
0307         end
0308     end
0309 
0310     fprintf(fd, '%%%% violation cost coefficients\n');
0311     fprintf(fd, '%%\trate_a_cost\n');
0312     fprintf(fd, '%ssoftlims.cost = [\n', prefix);
0313     fprintf(fd, '\t%g;\n', s.cost);
0314     fprintf(fd, '];\n');
0315 
0316     %% save output fields for solved case
0317     if isfield(mpc.softlims, 'overload')
0318         fprintf(fd, '\n%%%% overloads\n');
0319         fprintf(fd, '%%\toverload\n');
0320         fprintf(fd, '%ssoftlims.overload = [\n', prefix);
0321         fprintf(fd, '\t%g;\n', s.overload);
0322         fprintf(fd, '];\n');
0323 
0324         fprintf(fd, '\n%%%% overload costs\n');
0325         fprintf(fd, '%%\toverload_costs\n');
0326         fprintf(fd, '%ssoftlims.ovl_cost = [\n', prefix);
0327         fprintf(fd, '\t%g;\n', s.ovl_cost);
0328         fprintf(fd, '];\n');
0329     end
0330 end
0331 
0332 %%-----  softlims_defaults  --------------------------------------------
0333 function s = softlims_defaults(s, branch)
0334 %
0335 %   s = softlims_defaults(s, branch)
0336 %
0337 %   Takes a softlims struct that could have an empty 'idx' field or a
0338 %   scalar 'cost' field and fills them out with the defaults, where the
0339 %   default for 'idx' includes all on-line branches with non-zero RATE_A,
0340 %   and the default for the cost is to apply the scalar to each soft limit
0341 %   violation.
0342 
0343 %% define named indices into data matrices
0344 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0345     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0346     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0347 
0348 if isempty(s.idx)
0349     s.idx = find(branch(:, BR_STATUS) > 0 & branch(:, RATE_A) > 0);
0350 end
0351 if length(s.cost) == 1 && length(s.idx) > 1
0352     s.cost = s.cost * ones(size(s.idx));
0353 end