0001 function [results, success, raw] = opf_execute(om, mpopt)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031
0032
0033
0034
0035
0036
0037
0038
0039
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
0050
0051 dc = strcmp(upper(mpopt.model), 'DC');
0052 alg = upper(mpopt.opf.ac.solver);
0053 sdp = strcmp(alg, 'SDPOPF');
0054
0055
0056 om = build_cost_params(om);
0057
0058
0059 [vv, ll, nn] = get_idx(om);
0060
0061 if mpopt.verbose > 0
0062 v = mpver('all');
0063 fprintf('\nMATPOWER Version %s, %s', v.Version, v.Date);
0064 end
0065
0066
0067 if dc
0068 if mpopt.verbose > 0
0069 fprintf(' -- DC Optimal Power Flow\n');
0070 end
0071 [results, success, raw] = dcopf_solver(om, mpopt);
0072 else
0073
0074 if mpopt.verbose > 0
0075 fprintf(' -- AC Optimal Power Flow\n');
0076 end
0077
0078
0079 if strcmp(alg, 'DEFAULT')
0080 if have_fcn('pdipmopf')
0081 alg = 'PDIPM';
0082 else
0083 alg = 'MIPS';
0084 end
0085 mpopt = mpoption(mpopt, 'opf.ac.solver', alg);
0086 end
0087
0088
0089 switch alg
0090 case 'MIPS'
0091 [results, success, raw] = mipsopf_solver(om, mpopt);
0092 case 'IPOPT'
0093 if ~have_fcn('ipopt')
0094 error('opf_execute: MPOPT.opf.ac.solver = ''%s'' requires IPOPT (see https://projects.coin-or.org/Ipopt/)', alg);
0095 end
0096 [results, success, raw] = ipoptopf_solver(om, mpopt);
0097 case 'PDIPM'
0098 if mpopt.pdipm.step_control
0099 if ~have_fcn('scpdipmopf')
0100 error('opf_execute: MPOPT.opf.ac.solver = ''%s'' requires SCPDIPMOPF (see http://www.pserc.cornell.edu/tspopf/)', alg);
0101 end
0102 else
0103 if ~have_fcn('pdipmopf')
0104 error('opf_execute: MPOPT.opf.ac.solver = ''%s'' requires PDIPMOPF (see http://www.pserc.cornell.edu/tspopf/)', alg);
0105 end
0106 end
0107 [results, success, raw] = tspopf_solver(om, mpopt);
0108 case 'TRALM'
0109 if ~have_fcn('tralmopf')
0110 error('opf_execute: MPOPT.opf.ac.solver = ''%s'' requires TRALM (see http://www.pserc.cornell.edu/tspopf/)', alg);
0111 end
0112 [results, success, raw] = tspopf_solver(om, mpopt);
0113 case 'MINOPF'
0114 if ~have_fcn('minopf')
0115 error('opf_execute: MPOPT.opf.ac.solver = ''%s'' requires MINOPF (see http://www.pserc.cornell.edu/minopf/)', alg);
0116 end
0117 [results, success, raw] = mopf_solver(om, mpopt);
0118 case 'FMINCON'
0119 if ~have_fcn('fmincon')
0120 error('opf_execute: MPOPT.opf.ac.solver = ''%s'' requires FMINCON (Optimization Toolbox 2.x or later)', alg);
0121 end
0122 [results, success, raw] = fmincopf_solver(om, mpopt);
0123 case 'KNITRO'
0124 if ~have_fcn('knitro')
0125 error('opf_execute: MPOPT.opf.ac.solver = ''%s'' requires KNITRO (see http://www.ziena.com/)', alg);
0126 end
0127 [results, success, raw] = ktropf_solver(om, mpopt);
0128 case 'SDPOPF'
0129 if ~have_fcn('yalmip')
0130 error('opf_execute: MPOPT.opf.ac.solver = ''%s'' requires YALMIP (see http://users.isy.liu.se/johanl/yalmip/)', alg);
0131 end
0132 [results, success, raw] = sdpopf_solver(om, mpopt);
0133 otherwise
0134 error('opf_execute: MPOPT.opf.ac.solver = ''%s'' is not a valid AC OPF solver selection', alg);
0135 end
0136 end
0137 if ~isfield(raw, 'output') || ~isfield(raw.output, 'alg') || isempty(raw.output.alg)
0138 raw.output.alg = alg;
0139 end
0140
0141 if success
0142 if ~dc && ~sdp
0143
0144 results.gen(:, VG) = results.bus(results.gen(:, GEN_BUS), VM);
0145
0146
0147 if ll.N.PQh > 0 || ll.N.PQl > 0
0148 mu_PQh = results.mu.lin.l(ll.i1.PQh:ll.iN.PQh) - results.mu.lin.u(ll.i1.PQh:ll.iN.PQh);
0149 mu_PQl = results.mu.lin.l(ll.i1.PQl:ll.iN.PQl) - results.mu.lin.u(ll.i1.PQl:ll.iN.PQl);
0150 Apqdata = userdata(om, 'Apqdata');
0151 results.gen = update_mupq(results.baseMVA, results.gen, mu_PQh, mu_PQl, Apqdata);
0152 end
0153
0154
0155 if mpopt.opf.return_raw_der
0156
0157 if isfield(results, 'dg')
0158 raw.dg = results.dg;
0159 raw.g = results.g;
0160 end
0161
0162 if ~isfield(raw, 'dg')
0163 mpc = get_mpc(om);
0164 [Ybus, Yf, Yt] = makeYbus(mpc.baseMVA, mpc.bus, mpc.branch);
0165 [g, geq, dg, dgeq] = opf_consfcn(results.x, om, Ybus, Yf, Yt, mpopt);
0166 raw.g = [ geq; g];
0167 raw.dg = [ dgeq'; dg'];
0168 end
0169
0170 [f, df, d2f] = opf_costfcn(results.x, om);
0171 raw.df = df;
0172 raw.d2f = d2f;
0173 end
0174 end
0175
0176
0177 if isfield(results, 'dg')
0178 rmfield(results, 'dg');
0179 rmfield(results, 'g');
0180 end
0181
0182
0183 if ~sdp && ll.N.ang > 0
0184 iang = userdata(om, 'iang');
0185 results.branch(iang, MU_ANGMIN) = results.mu.lin.l(ll.i1.ang:ll.iN.ang) * pi/180;
0186 results.branch(iang, MU_ANGMAX) = results.mu.lin.u(ll.i1.ang:ll.iN.ang) * pi/180;
0187 end
0188 else
0189
0190 if ~dc && mpopt.opf.return_raw_der
0191 raw.dg = [];
0192 raw.g = [];
0193 raw.df = [];
0194 raw.d2f = [];
0195 end
0196 end
0197
0198 if ~sdp
0199
0200 om_var_order = get(om, 'var', 'order');
0201 for k = 1:length(om_var_order)
0202 name = om_var_order(k).name;
0203 if getN(om, 'var', name)
0204 idx = vv.i1.(name):vv.iN.(name);
0205 results.var.val.(name) = results.x(idx);
0206 results.var.mu.l.(name) = results.mu.var.l(idx);
0207 results.var.mu.u.(name) = results.mu.var.u(idx);
0208 end
0209 end
0210
0211
0212 om_lin_order = get(om, 'lin', 'order');
0213 for k = 1:length(om_lin_order)
0214 name = om_lin_order(k).name;
0215 if getN(om, 'lin', name)
0216 idx = ll.i1.(name):ll.iN.(name);
0217 results.lin.mu.l.(name) = results.mu.lin.l(idx);
0218 results.lin.mu.u.(name) = results.mu.lin.u(idx);
0219 end
0220 end
0221
0222
0223 if ~dc
0224 om_nln_order = get(om, 'nln', 'order');
0225 for k = 1:length(om_nln_order)
0226 name = om_nln_order(k).name;
0227 if getN(om, 'nln', name)
0228 idx = nn.i1.(name):nn.iN.(name);
0229 results.nln.mu.l.(name) = results.mu.nln.l(idx);
0230 results.nln.mu.u.(name) = results.mu.nln.u(idx);
0231 end
0232 end
0233 end
0234
0235
0236 om_cost_order = get(om, 'cost', 'order');
0237 for k = 1:length(om_cost_order)
0238 name = om_cost_order(k).name;
0239 if getN(om, 'cost', name)
0240 results.cost.(name) = compute_cost(om, results.x, name);
0241 end
0242 end
0243
0244
0245
0246
0247 pwl1 = userdata(om, 'pwl1');
0248 if ~isempty(pwl1) && ~strcmp(alg, 'TRALM') && ~(strcmp(alg, 'PDIPM') && mpopt.pdipm.step_control)
0249
0250 vv = get_idx(om);
0251 if dc
0252 nx = vv.iN.Pg;
0253 else
0254 nx = vv.iN.Qg;
0255 end
0256 y = zeros(length(pwl1), 1);
0257 raw.xr = [ raw.xr(1:nx); y; raw.xr(nx+1:end)];
0258 results.x = [ results.x(1:nx); y; results.x(nx+1:end)];
0259 end
0260 end