Description of opf

0001 function [results, success, raw] = opf_execute(om, mpopt)
0002 %OPF_EXECUTE  Executes the OPF specified by an OPF model object.
0003 %   [RESULTS, SUCCESS, RAW] = OPF_EXECUTE(OM, MPOPT)
0004 %
0005 %   RESULTS are returned with internal indexing, all equipment
0006 %   in-service, etc.
0007 %
0008 %   See also OPF, OPF_SETUP.
0009 
0010 %   MATPOWER
0011 %   $Id: opf_execute.m,v 1.7 2011/06/17 20:28:44 cvs Exp $
0012 %   by Ray Zimmerman, PSERC Cornell
0013 %   Copyright (c) 2009-2010 by Power System Engineering Research Center (PSERC)
0014 %
0015 %   This file is part of MATPOWER.
0016 %   See http://www.pserc.cornell.edu/matpower/ for more info.
0017 %
0018 %   MATPOWER is free software: you can redistribute it and/or modify
0019 %   it under the terms of the GNU General Public License as published
0020 %   by the Free Software Foundation, either version 3 of the License,
0021 %   or (at your option) any later version.
0022 %
0023 %   MATPOWER is distributed in the hope that it will be useful,
0024 %   but WITHOUT ANY WARRANTY; without even the implied warranty of
0025 %   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
0026 %   GNU General Public License for more details.
0027 %
0028 %   You should have received a copy of the GNU General Public License
0029 %   along with MATPOWER. If not, see <http://www.gnu.org/licenses/>.
0030 %
0031 %   Additional permission under GNU GPL version 3 section 7
0032 %
0033 %   If you modify MATPOWER, or any covered work, to interface with
0034 %   other modules (such as MATLAB code and MEX-files) available in a
0035 %   MATLAB(R) or comparable environment containing parts covered
0036 %   under other licensing terms, the licensors of MATPOWER grant
0037 %   you additional permission to convey the resulting work.
0038 
0039 %% define named indices into data matrices
0040 [PQ, PV, REF, NONE, BUS_I, BUS_TYPE, PD, QD, GS, BS, BUS_AREA, VM, ...
0041     VA, BASE_KV, ZONE, VMAX, VMIN, LAM_P, LAM_Q, MU_VMAX, MU_VMIN] = idx_bus;
0042 [GEN_BUS, PG, QG, QMAX, QMIN, VG, MBASE, GEN_STATUS, PMAX, PMIN, ...
0043     MU_PMAX, MU_PMIN, MU_QMAX, MU_QMIN, PC1, PC2, QC1MIN, QC1MAX, ...
0044     QC2MIN, QC2MAX, RAMP_AGC, RAMP_10, RAMP_30, RAMP_Q, APF] = idx_gen;
0045 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0046     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0047     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0048 
0049 %%-----  setup  -----
0050 %% options
0051 dc  = mpopt(10);        %% PF_DC        : 1 = DC OPF, 0 = AC OPF
0052 alg = mpopt(11);        %% OPF_ALG
0053 verbose = mpopt(31);    %% VERBOSE
0054 
0055 %% build user-defined costs
0056 om = build_cost_params(om);
0057 
0058 %% get indexing
0059 [vv, ll, nn] = get_idx(om);
0060 
0061 if verbose > 0
0062     v = mpver('all');
0063     fprintf('\nMATPOWER Version %s, %s', v.Version, v.Date);
0064 end
0065 
0066 %%-----  run DC OPF solver  -----
0067 if dc
0068   if verbose > 0
0069     fprintf(' -- DC Optimal Power Flow\n');
0070   end
0071   [results, success, raw] = dcopf_solver(om, mpopt);
0072 else
0073   %%-----  run AC OPF solver  -----
0074   if verbose > 0
0075     fprintf(' -- AC Optimal Power Flow\n');
0076   end
0077 
0078   %% if OPF_ALG not set, choose best available option
0079   if alg == 0
0080     if have_fcn('pdipmopf')
0081       alg = 540;                %% PDIPM
0082     else
0083       alg = 560;                %% MIPS
0084     end
0085   end
0086 
0087   %% update deprecated algorithm codes to new, generalized formulation equivalents
0088   if alg == 100 || alg == 200        %% CONSTR
0089     alg = 300;
0090   elseif alg == 120 || alg == 220    %% dense LP
0091     alg = 320;
0092   elseif alg == 140 || alg == 240    %% sparse (relaxed) LP
0093     alg = 340;
0094   elseif alg == 160 || alg == 260    %% sparse (full) LP
0095     alg = 360;
0096   end
0097   mpopt(11) = alg;
0098 
0099   %% run specific AC OPF solver
0100   if alg == 560 || alg == 565                   %% MIPS
0101     if have_fcn('anon_fcns')
0102       solver = @mipsopf_solver;
0103     else
0104       solver = @mips6opf_solver;
0105     end
0106     [results, success, raw] = feval(solver, om, mpopt);
0107   elseif alg == 580                             %% IPOPT
0108     if ~have_fcn('ipopt')
0109       error('opf_execute: OPF_ALG %d requires IPOPT (see https://projects.coin-or.org/Ipopt/)', alg);
0110     end
0111     [results, success, raw] = ipoptopf_solver(om, mpopt);
0112   elseif alg == 540 || alg == 545 || alg == 550 %% PDIPM_OPF, SCPDIPM_OPF, or TRALM_OPF
0113     if alg == 540                               %% PDIPM_OPF
0114       if ~have_fcn('pdipmopf')
0115         error('opf_execute: OPF_ALG %d requires PDIPMOPF (see http://www.pserc.cornell.edu/tspopf/)', alg);
0116       end
0117     elseif alg == 545                           %% SCPDIPM_OPF
0118       if ~have_fcn('scpdipmopf')
0119         error('opf_execute: OPF_ALG %d requires SCPDIPMOPF (see http://www.pserc.cornell.edu/tspopf/)', alg);
0120       end
0121     elseif alg == 550                           %% TRALM_OPF
0122       if ~have_fcn('tralmopf')
0123         error('opf_execute: OPF_ALG %d requires TRALM (see http://www.pserc.cornell.edu/tspopf/)', alg);
0124       end
0125     end
0126     [results, success, raw] = tspopf_solver(om, mpopt);
0127   elseif alg == 500                             %% MINOPF
0128     if ~have_fcn('minopf')
0129       error('opf_execute: OPF_ALG %d requires MINOPF (see http://www.pserc.cornell.edu/minopf/)', alg);
0130     end
0131     [results, success, raw] = mopf_solver(om, mpopt);
0132   elseif alg == 520                             %% FMINCON
0133     if ~have_fcn('fmincon')
0134       error('opf_execute: OPF_ALG %d requires FMINCON (Optimization Toolbox 2.x or later)', alg);
0135     end
0136     if have_fcn('anon_fcns')
0137       solver = @fmincopf_solver;
0138     else
0139       solver = @fmincopf6_solver;
0140     end
0141     [results, success, raw] = feval(solver, om, mpopt);
0142   elseif alg == 600                             %% KNITRO
0143     if ~have_fcn('knitro')
0144       error('opf_execute: OPF_ALG %d requires KNITRO (see http://www.ziena.com/)', alg);
0145     end
0146     [results, success, raw] = ktropf_solver(om, mpopt);
0147   elseif alg == 300                             %% CONSTR
0148     if ~have_fcn('constr')
0149       error('opf_execute: OPF_ALG %d requires CONSTR (Optimization Toolbox 1.x)', alg);
0150     end
0151     [results, success, raw] = copf_solver(om, mpopt);
0152   elseif alg == 320 || alg == 340 || alg == 360 %% LP
0153     [results, success, raw] = lpopf_solver(om, mpopt);
0154   else
0155     error('opf_execute: OPF_ALG %d is not a valid algorithm code', alg);
0156   end
0157 end
0158 if ~isfield(raw, 'output') || ~isfield(raw.output, 'alg') || isempty(raw.output.alg)
0159     raw.output.alg = alg;
0160 end
0161 
0162 if success
0163   if ~dc
0164     %% copy bus voltages back to gen matrix
0165     results.gen(:, VG) = results.bus(results.gen(:, GEN_BUS), VM);
0166 
0167     %% gen PQ capability curve multipliers
0168     if ll.N.PQh > 0 || ll.N.PQl > 0
0169       mu_PQh = results.mu.lin.l(ll.i1.PQh:ll.iN.PQh) - results.mu.lin.u(ll.i1.PQh:ll.iN.PQh);
0170       mu_PQl = results.mu.lin.l(ll.i1.PQl:ll.iN.PQl) - results.mu.lin.u(ll.i1.PQl:ll.iN.PQl);
0171       Apqdata = userdata(om, 'Apqdata');
0172       results.gen = update_mupq(results.baseMVA, results.gen, mu_PQh, mu_PQl, Apqdata);
0173     end
0174 
0175     %% compute g, dg, f, df, d2f if requested by RETURN_RAW_DER = 1
0176     if mpopt(52)
0177       %% move from results to raw if using v4.0 of MINOPF or TSPOPF
0178       if isfield(results, 'dg')
0179         raw.dg = results.dg;
0180         raw.g = results.g;
0181       end
0182       %% compute g, dg, unless already done by post-v4.0 MINOPF or TSPOPF
0183       if ~isfield(raw, 'dg')
0184         mpc = get_mpc(om);
0185         [Ybus, Yf, Yt] = makeYbus(mpc.baseMVA, mpc.bus, mpc.branch);
0186         [g, geq, dg, dgeq] = opf_consfcn(results.x, om, Ybus, Yf, Yt, mpopt);
0187         raw.g = [ geq; g];
0188         raw.dg = [ dgeq'; dg'];   %% true Jacobian organization
0189       end
0190       %% compute df, d2f
0191       [f, df, d2f] = opf_costfcn(results.x, om);
0192       raw.df = df;
0193       raw.d2f = d2f;
0194     end
0195   end
0196 
0197   %% delete g and dg fieldsfrom results if using v4.0 of MINOPF or TSPOPF
0198   if isfield(results, 'dg')
0199     rmfield(results, 'dg');
0200     rmfield(results, 'g');
0201   end
0202 
0203   %% angle limit constraint multipliers
0204   if ll.N.ang > 0
0205     iang = userdata(om, 'iang');
0206     results.branch(iang, MU_ANGMIN) = results.mu.lin.l(ll.i1.ang:ll.iN.ang) * pi/180;
0207     results.branch(iang, MU_ANGMAX) = results.mu.lin.u(ll.i1.ang:ll.iN.ang) * pi/180;
0208   end
0209 else
0210   %% assign empty g, dg, f, df, d2f if requested by RETURN_RAW_DER = 1
0211   if ~dc && mpopt(52)
0212     raw.dg = [];
0213     raw.g = [];
0214     raw.df = [];
0215     raw.d2f = [];
0216   end
0217 end
0218 
0219 %% assign values and limit shadow prices for variables
0220 om_var_order = get(om, 'var', 'order');
0221 for k = 1:length(om_var_order)
0222   name = om_var_order{k};
0223   if getN(om, 'var', name)
0224     idx = vv.i1.(name):vv.iN.(name);
0225     results.var.val.(name) = results.x(idx);
0226     results.var.mu.l.(name) = results.mu.var.l(idx);
0227     results.var.mu.u.(name) = results.mu.var.u(idx);
0228   end
0229 end
0230 
0231 %% assign shadow prices for linear constraints
0232 om_lin_order = get(om, 'lin', 'order');
0233 for k = 1:length(om_lin_order)
0234   name = om_lin_order{k};
0235   if getN(om, 'lin', name)
0236     idx = ll.i1.(name):ll.iN.(name);
0237     results.lin.mu.l.(name) = results.mu.lin.l(idx);
0238     results.lin.mu.u.(name) = results.mu.lin.u(idx);
0239   end
0240 end
0241 
0242 %% assign shadow prices for nonlinear constraints
0243 if ~dc
0244   om_nln_order = get(om, 'nln', 'order');
0245   for k = 1:length(om_nln_order)
0246     name = om_nln_order{k};
0247     if getN(om, 'nln', name)
0248       idx = nn.i1.(name):nn.iN.(name);
0249       results.nln.mu.l.(name) = results.mu.nln.l(idx);
0250       results.nln.mu.u.(name) = results.mu.nln.u(idx);
0251     end
0252   end
0253 end
0254 
0255 %% assign values for components of user cost
0256 om_cost_order = get(om, 'cost', 'order');
0257 for k = 1:length(om_cost_order)
0258   name = om_cost_order{k};
0259   if getN(om, 'cost', name)
0260     results.cost.(name) = compute_cost(om, results.x, name);
0261   end
0262 end
0263 
0264 %% if single-block PWL costs were converted to POLY, insert dummy y into x
0265 %% Note: The "y" portion of x will be nonsense, but everything should at
0266 %%       least be in the expected locations.
0267 pwl1 = userdata(om, 'pwl1');
0268 if ~isempty(pwl1) && alg ~= 545 && alg ~= 550
0269   %% get indexing
0270   vv = get_idx(om);
0271   if dc
0272     nx = vv.iN.Pg;
0273   else
0274     nx = vv.iN.Qg;
0275   end
0276   y = zeros(length(pwl1), 1);
0277   raw.xr = [ raw.xr(1:nx); y; raw.xr(nx+1:end)];
0278   results.x = [ results.x(1:nx); y; results.x(nx+1:end)];
0279 end
0280
opf_execute

PURPOSE

SYNOPSIS

DESCRIPTION

CROSS-REFERENCE INFORMATION

SOURCE CODE
