0001 function om = opf_setup(mpc, 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
0041 dc = strcmp(upper(mpopt.model), 'DC');
0042 alg = upper(mpopt.opf.ac.solver);
0043
0044
0045 [PQ, PV, REF, NONE, BUS_I, BUS_TYPE, PD, QD, GS, BS, BUS_AREA, VM, ...
0046 VA, BASE_KV, ZONE, VMAX, VMIN, LAM_P, LAM_Q, MU_VMAX, MU_VMIN] = idx_bus;
0047 [GEN_BUS, PG, QG, QMAX, QMIN, VG, MBASE, GEN_STATUS, PMAX, PMIN, ...
0048 MU_PMAX, MU_PMIN, MU_QMAX, MU_QMIN, PC1, PC2, QC1MIN, QC1MAX, ...
0049 QC2MIN, QC2MAX, RAMP_AGC, RAMP_10, RAMP_30, RAMP_Q, APF] = idx_gen;
0050 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0051 TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0052 ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0053 [PW_LINEAR, POLYNOMIAL, MODEL, STARTUP, SHUTDOWN, NCOST, COST] = idx_cost;
0054
0055
0056 nb = size(mpc.bus, 1);
0057 nl = size(mpc.branch, 1);
0058 ng = size(mpc.gen, 1);
0059 if isfield(mpc, 'A')
0060 nusr = size(mpc.A, 1);
0061 else
0062 nusr = 0;
0063 end
0064 if isfield(mpc, 'N')
0065 nw = size(mpc.N, 1);
0066 else
0067 nw = 0;
0068 end
0069
0070 if dc
0071
0072 mpc.gencost = pqcost(mpc.gencost, ng);
0073
0074
0075 if nusr || nw
0076 acc = [nb+(1:nb) 2*nb+ng+(1:ng)];
0077 if nusr && size(mpc.A, 2) >= 2*nb + 2*ng
0078
0079 if any(any(mpc.A(:, acc)))
0080 error('opf_setup: attempting to solve DC OPF with user constraints on Vm or Qg');
0081 end
0082 mpc.A(:, acc) = [];
0083 end
0084 if nw && size(mpc.N, 2) >= 2*nb + 2*ng
0085
0086 if any(any(mpc.N(:, acc)))
0087 [ii, jj] = find(mpc.N(:, acc));
0088 ii = unique(ii);
0089 if any(mpc.Cw(ii)) || (isfield(mpc, 'H') && ~isempty(mpc.H) && ...
0090 any(any(mpc.H(:, ii))))
0091 error('opf_setup: attempting to solve DC OPF with user costs on Vm or Qg');
0092 end
0093 end
0094 mpc.N(:, acc) = [];
0095 end
0096 end
0097 end
0098
0099
0100 pwl1 = find(mpc.gencost(:, MODEL) == PW_LINEAR & mpc.gencost(:, NCOST) == 2);
0101
0102 if ~isempty(pwl1)
0103 x0 = mpc.gencost(pwl1, COST);
0104 y0 = mpc.gencost(pwl1, COST+1);
0105 x1 = mpc.gencost(pwl1, COST+2);
0106 y1 = mpc.gencost(pwl1, COST+3);
0107 m = (y1 - y0) ./ (x1 - x0);
0108 b = y0 - m .* x0;
0109 mpc.gencost(pwl1, MODEL) = POLYNOMIAL;
0110 mpc.gencost(pwl1, NCOST) = 2;
0111 mpc.gencost(pwl1, COST:COST+1) = [m b];
0112 end
0113
0114
0115 [baseMVA, bus, gen, branch, gencost, Au, lbu, ubu, mpopt, ...
0116 N, fparm, H, Cw, z0, zl, zu, userfcn] = opf_args(mpc, mpopt);
0117
0118
0119 refs = find(bus(:, BUS_TYPE) == REF);
0120 if length(refs) > 1 && mpopt.verbose > 0
0121 errstr = ['\nopf_setup: Warning: Multiple reference buses.\n', ...
0122 ' For a system with islands, a reference bus in each island\n', ...
0123 ' may help convergence, but in a fully connected system such\n', ...
0124 ' a situation is probably not reasonable.\n\n' ];
0125 fprintf(errstr);
0126 end
0127
0128
0129 Va = bus(:, VA) * (pi/180);
0130 Vm = bus(:, VM);
0131 Pg = gen(:, PG) / baseMVA;
0132 Qg = gen(:, QG) / baseMVA;
0133 Pmin = gen(:, PMIN) / baseMVA;
0134 Pmax = gen(:, PMAX) / baseMVA;
0135 Qmin = gen(:, QMIN) / baseMVA;
0136 Qmax = gen(:, QMAX) / baseMVA;
0137
0138 if dc
0139
0140 nv = 0;
0141 nq = 0;
0142 q1 = [];
0143
0144
0145 [B, Bf, Pbusinj, Pfinj] = makeBdc(baseMVA, bus, branch);
0146 neg_Cg = sparse(gen(:, GEN_BUS), 1:ng, -1, nb, ng);
0147 Amis = [B neg_Cg];
0148 bmis = -(bus(:, PD) + bus(:, GS)) / baseMVA - Pbusinj;
0149
0150
0151 il = find(branch(:, RATE_A) ~= 0 & branch(:, RATE_A) < 1e10);
0152 nl2 = length(il);
0153 if nl2
0154 upf = branch(il, RATE_A) / baseMVA - Pfinj(il);
0155 upt = branch(il, RATE_A) / baseMVA + Pfinj(il);
0156 else
0157 upf = [];
0158 upt = [];
0159 end
0160
0161 user_vars = {'Va', 'Pg'};
0162 ycon_vars = {'Pg', 'y'};
0163 else
0164
0165 nv = nb;
0166 nq = ng;
0167 q1 = 1+ng;
0168
0169
0170 [Avl, lvl, uvl] = makeAvl(baseMVA, gen);
0171
0172
0173 [Apqh, ubpqh, Apql, ubpql, Apqdata] = makeApq(baseMVA, gen);
0174
0175 user_vars = {'Va', 'Vm', 'Pg', 'Qg'};
0176 ycon_vars = {'Pg', 'Qg', 'y'};
0177 end
0178
0179
0180 Vau = Inf * ones(nb, 1);
0181 Val = -Vau;
0182 Vau(refs) = Va(refs);
0183 Val(refs) = Va(refs);
0184
0185
0186 [Aang, lang, uang, iang] = makeAang(baseMVA, branch, nb, mpopt);
0187
0188
0189 if (strcmp(alg, 'PDIPM') && mpopt.pdipm.step_control) || strcmp(alg, 'TRALM')
0190
0191 ny = 0;
0192 Ay = sparse(0, ng+nq);
0193 by =[];
0194 else
0195 ipwl = find(gencost(:, MODEL) == PW_LINEAR);
0196 ny = size(ipwl, 1);
0197 [Ay, by] = makeAy(baseMVA, ng, gencost, 1, q1, 1+ng+nq);
0198 end
0199 if any(gencost(:, MODEL) ~= POLYNOMIAL & gencost(:, MODEL) ~= PW_LINEAR)
0200 error('opf_setup: some generator cost rows have invalid MODEL value');
0201 end
0202
0203
0204
0205 nx = nb+nv + ng+nq;
0206 if nusr
0207 nz = size(mpc.A, 2) - nx;
0208 if nz < 0
0209 error('opf_setup: user supplied A matrix must have at least %d columns.', nx);
0210 end
0211 else
0212 nz = 0;
0213 if nw
0214 if size(mpc.N, 2) ~= nx;
0215 error('opf_setup: user supplied N matrix must have %d columns.', nx);
0216 end
0217 end
0218 end
0219
0220
0221 om = opf_model(mpc);
0222 if ~isempty(pwl1)
0223 om = userdata(om, 'pwl1', pwl1);
0224 end
0225 if dc
0226 om = userdata(om, 'Bf', Bf);
0227 om = userdata(om, 'Pfinj', Pfinj);
0228 om = userdata(om, 'iang', iang);
0229 om = add_vars(om, 'Va', nb, Va, Val, Vau);
0230 om = add_vars(om, 'Pg', ng, Pg, Pmin, Pmax);
0231 om = add_constraints(om, 'Pmis', Amis, bmis, bmis, {'Va', 'Pg'});
0232 om = add_constraints(om, 'Pf', Bf(il,:), -upt, upf, {'Va'});
0233 om = add_constraints(om, 'ang', Aang, lang, uang, {'Va'});
0234 else
0235 om = userdata(om, 'Apqdata', Apqdata);
0236 om = userdata(om, 'iang', iang);
0237 om = add_vars(om, 'Va', nb, Va, Val, Vau);
0238 om = add_vars(om, 'Vm', nb, Vm, bus(:, VMIN), bus(:, VMAX));
0239 om = add_vars(om, 'Pg', ng, Pg, Pmin, Pmax);
0240 om = add_vars(om, 'Qg', ng, Qg, Qmin, Qmax);
0241 om = add_constraints(om, 'Pmis', nb, 'nonlinear');
0242 om = add_constraints(om, 'Qmis', nb, 'nonlinear');
0243 om = add_constraints(om, 'Sf', nl, 'nonlinear');
0244 om = add_constraints(om, 'St', nl, 'nonlinear');
0245 om = add_constraints(om, 'PQh', Apqh, [], ubpqh, {'Pg', 'Qg'});
0246 om = add_constraints(om, 'PQl', Apql, [], ubpql, {'Pg', 'Qg'});
0247 om = add_constraints(om, 'vl', Avl, lvl, uvl, {'Pg', 'Qg'});
0248 om = add_constraints(om, 'ang', Aang, lang, uang, {'Va'});
0249 end
0250
0251
0252 if ny > 0
0253 om = add_vars(om, 'y', ny);
0254 om = add_constraints(om, 'ycon', Ay, [], by, ycon_vars);
0255 end
0256
0257
0258 if nz > 0
0259 om = add_vars(om, 'z', nz, z0, zl, zu);
0260 user_vars{end+1} = 'z';
0261 end
0262 if nusr
0263 om = add_constraints(om, 'usr', mpc.A, lbu, ubu, user_vars);
0264 end
0265 if nw
0266 user_cost.N = mpc.N;
0267 user_cost.Cw = Cw;
0268 if ~isempty(fparm)
0269 user_cost.dd = fparm(:, 1);
0270 user_cost.rh = fparm(:, 2);
0271 user_cost.kk = fparm(:, 3);
0272 user_cost.mm = fparm(:, 4);
0273 end
0274 if ~isempty(H)
0275 user_cost.H = H;
0276 end
0277 om = add_costs(om, 'usr', user_cost, user_vars);
0278 end
0279
0280
0281 om = run_userfcn(userfcn, 'formulation', om);