make_opf_feasible

MAKE_OPF_FEASIBLE Attempts to relax constraints to make an OPF feasible.
   [RESULTS, CHGS] = MAKE_OPF_FEASIBLE(MPC, MPOPT) 

   Attempts to automate the process of finding a feasible OPF solution when
   starting with an infeasible case.

   MPC - initial (possibly infeasible) MATPOWER case struct
   MPOPT - (optional) MATPOWER options struct.

   Work-in-progress:   CHGS are currently not returned.
                       At this point, the code *is* the documentation.

   See also RUNOPF.

0001 function [r, chgs] = make_opf_feasible(mpc, mpopt)
0002 %MAKE_OPF_FEASIBLE Attempts to relax constraints to make an OPF feasible.
0003 %   [RESULTS, CHGS] = MAKE_OPF_FEASIBLE(MPC, MPOPT)
0004 %
0005 %   Attempts to automate the process of finding a feasible OPF solution when
0006 %   starting with an infeasible case.
0007 %
0008 %   MPC - initial (possibly infeasible) MATPOWER case struct
0009 %   MPOPT - (optional) MATPOWER options struct.
0010 %
0011 %   Work-in-progress:   CHGS are currently not returned.
0012 %                       At this point, the code *is* the documentation.
0013 %
0014 %   See also RUNOPF.
0015 
0016 %   MATPOWER
0017 %   $Id: make_opf_feasible.m 2471 2014-12-16 15:32:35Z ray $
0018 %   by Ray Zimmerman, PSERC Cornell
0019 %   Copyright (c) 2014 by Power System Engineering Research Center (PSERC)
0020 %
0021 %   This file is part of MATPOWER.
0022 %   See http://www.pserc.cornell.edu/matpower/ for more info.
0023 %
0024 %   MATPOWER is free software: you can redistribute it and/or modify
0025 %   it under the terms of the GNU General Public License as published
0026 %   by the Free Software Foundation, either version 3 of the License,
0027 %   or (at your option) any later version.
0028 %
0029 %   MATPOWER is distributed in the hope that it will be useful,
0030 %   but WITHOUT ANY WARRANTY; without even the implied warranty of
0031 %   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
0032 %   GNU General Public License for more details.
0033 %
0034 %   You should have received a copy of the GNU General Public License
0035 %   along with MATPOWER. If not, see <http://www.gnu.org/licenses/>.
0036 %
0037 %   Additional permission under GNU GPL version 3 section 7
0038 %
0039 %   If you modify MATPOWER, or any covered work, to interface with
0040 %   other modules (such as MATLAB code and MEX-files) available in a
0041 %   MATLAB(R) or comparable environment containing parts covered
0042 %   under other licensing terms, the licensors of MATPOWER grant
0043 %   you additional permission to convey the resulting work.
0044 
0045 %% define constants
0046 [PQ, PV, REF, NONE, BUS_I, BUS_TYPE, PD, QD, GS, BS, BUS_AREA, VM, ...
0047     VA, BASE_KV, ZONE, VMAX, VMIN, LAM_P, LAM_Q, MU_VMAX, MU_VMIN] = idx_bus;
0048 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0049     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0050     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0051 [GEN_BUS, PG, QG, QMAX, QMIN, VG, MBASE, GEN_STATUS, PMAX, PMIN, ...
0052     MU_PMAX, MU_PMIN, MU_QMAX, MU_QMIN, PC1, PC2, QC1MIN, QC1MAX, ...
0053     QC2MIN, QC2MAX, RAMP_AGC, RAMP_10, RAMP_30, RAMP_Q, APF] = idx_gen;
0054 
0055 %%-----  process inputs  -----
0056 if nargin < 2
0057     mpopt = mpoption();
0058 end
0059 
0060 %% initialize
0061 msg = '';
0062 done = 0;
0063 rc = 0;     %% return code
0064 
0065 %%-----  on-line generation capacity <= fixed load  -----
0066 ig = find(~isload(mpc.gen) & mpc.gen(:, GEN_STATUS) > 0);
0067 total_Pmax = sum(mpc.gen(ig, PMAX));
0068 total_Pd   = total_load(mpc.bus, [], 'all');
0069 if total_Pmax <= total_Pd
0070     msg = sprintf('%stotal fixed load (%g MW) > total on-line generation capacity (%g MW)\n', ...
0071         msg, total_Pd, total_Pmax);
0072     r = mpc;
0073     r.success = 0;
0074     r.f = NaN;
0075     rc = -2;        %% inadequate on-line generation capacity
0076     done = 1;
0077 end
0078 
0079 %%-----  check for infeasible limits  -----
0080 if ~done
0081     %% generator P limits
0082     g = find(mpc.gen(:, PMIN) > mpc.gen(:, PMAX) & mpc.gen(:, GEN_STATUS) > 0);
0083     if ~isempty(g)
0084         for k = 1:length(g)
0085             msg = sprintf('%sPmax (%g) < Pmin (%g) for generator %d at bus %d, ', ...
0086                 msg, mpc.gen(g(k), PMAX), mpc.gen(g(k), PMIN), g(k), mpc.gen(g(k), GEN_BUS));
0087             %% decide which limit to move based on where Pg is
0088             if mpc.gen(g(k), PG) >= mpc.gen(g(k), PMIN)
0089                 mpc.gen(g(k), PMAX) = mpc.gen(g(k), PMIN);      %% move Pmax
0090                 msg = sprintf('%smove Pmax\n', msg);
0091             elseif mpc.gen(g(k), PG) <= mpc.gen(g(k), PMAX)
0092                 mpc.gen(g(k), PMAX) = mpc.gen(g(k), PMIN);      %% move Pmin
0093                 msg = sprintf('%smove Pmax\n', msg);
0094             else
0095                 tmp = mpc.gen(g(k), PMAX);                      %% swap 'em
0096                 mpc.gen(g(k), PMAX) = mpc.gen(g(k), PMIN);
0097                 mpc.gen(g(k), PMIN) = tmp;
0098                 msg = sprintf('%sswap Pmin and Pmax\n', msg);
0099             end
0100         end
0101         rc = 2;         %% fixed generation capacity limits
0102     end
0103     
0104     %% generator Q limits
0105 %     rc = 3;         %% fixed reactive generation capacity limits
0106 
0107     %% voltage angle limits
0108 %     rc = 4;         %% fixed voltage angle limits
0109 
0110     %% voltage magnitude limits
0111 %     rc = 5;         %% fixed voltage magnitude limits
0112 
0113     %% negative branch flow limits
0114 %     rc = 6;         %% fixed negative branch flow limits
0115 
0116     %% run initial (infeasible?) case
0117     r = runopf(mpc, mpopt);
0118     if r.success
0119         done = 1;
0120     end
0121 end
0122 
0123 %%-----  attempt with short term branch ratings  -----
0124 if ~done
0125     %% get branch limits
0126     rate_a = mpc.branch(:, RATE_A);
0127     rate_b = mpc.branch(:, RATE_B);
0128     rate_c = mpc.branch(:, RATE_C);
0129     rate_a(rate_a == 0) = Inf;
0130     rate_b(rate_b == 0) = Inf;
0131     rate_c(rate_c == 0) = Inf;
0132 
0133     %% set short term limits
0134     rating = max([rate_a rate_b], [], 2);
0135     rating(isinf(rating)) = 0;
0136     mpc.branch(:, RATE_A) = rating;
0137     r = runopf(mpc, mpopt);
0138     if r.success
0139         done = 1;
0140         rc = 7;         %% using short-term branch ratings
0141         msg = sprintf('%susing short-term branch ratings\n', msg);
0142     end
0143 end
0144 
0145 %%-----  attempt with emergency branch ratings  -----
0146 if ~done
0147     %% set emergency limits
0148     rating = max([rate_a rate_b rate_c], [], 2);
0149     rating(isinf(rating)) = 0;
0150     mpc.branch(:, RATE_A) = rating;
0151     r = runopf(mpc, mpopt);
0152 
0153     if r.success
0154         done = 1;
0155         rc = 8;         %% using emergency branch ratings
0156         msg = sprintf('%susing emergency branch ratings\n', msg);
0157     end
0158 end
0159 
0160 %%-----  attempt without branch flow limits  -----
0161 if ~done
0162     %% try with no line limits
0163     mpc.branch(:, RATE_A) = 0;
0164     r = runopf(mpc, mpopt);
0165                 
0166     if ~r.success
0167         done = 1;
0168         msg = sprintf('%sstill infeasible without line limits\n', msg);
0169     else
0170         %% the following lines should be eliminated when I do further
0171         %% relaxations
0172         msg = sprintf('%susing no branch limits\n', msg);
0173         rc = 9;         %% using no branch ratings
0174     end
0175 end
0176 
0177 %%-----  try again, relaxing only violated limits  -----
0178 if ~done
0179 
0180 end
0181 
0182 %% include msg in output
0183 r.msg = msg;
0184 r.rc  = rc;
0185 if mpopt.verbose
0186     fprintf('%s', r.msg);
0187 end