Home > matpower7.1 > lib > t > t_vdep_load.m

t_vdep_load

PURPOSE ^

T_VDEP_LOAD Test voltage dependent ZIP load model for PF, CPF, OPF.

SYNOPSIS ^

function t_vdep_load(quiet)

DESCRIPTION ^

T_VDEP_LOAD    Test voltage dependent ZIP load model for PF, CPF, OPF.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function t_vdep_load(quiet)
0002 %T_VDEP_LOAD    Test voltage dependent ZIP load model for PF, CPF, OPF.
0003 
0004 %   MATPOWER
0005 %   Copyright (c) 2009-2020, Power Systems Engineering Research Center (PSERC)
0006 %   by Ray Zimmerman, PSERC Cornell
0007 %
0008 %   This file is part of MATPOWER.
0009 %   Covered by the 3-clause BSD License (see LICENSE file for details).
0010 %   See https://matpower.org for more info.
0011 
0012 if nargin < 1
0013     quiet = 0;
0014 end
0015 
0016 pfalgs = {'NR', 'FDXB', 'FDBX'};
0017 cpfparm = {'natural', 'arc-len', 'pseudo arc-len'};
0018 
0019 t_begin(33+39*length(pfalgs)+39*2*length(cpfparm), quiet);
0020 
0021 if quiet
0022     verbose = 0;
0023 else
0024     verbose = 0;
0025 end
0026 verbose = 0;
0027 
0028 casefile = 't_case9_opf';
0029 mpopt = mpoption;
0030 mpopt = mpoption(mpopt, 'out.gen', 1);
0031 mpopt = mpoption(mpopt, 'pf.tol', 1e-9);
0032 mpopt = mpoption(mpopt, 'opf.ac.solver', 'MIPS');
0033 mpopt = mpoption(mpopt, 'opf.violation', 1e-6, 'mips.gradtol', 1e-9, ...
0034         'mips.comptol', 1e-8, 'mips.costtol', 1e-8);
0035 mpopt = mpoption(mpopt, 'verbose', verbose);
0036 if ~verbose
0037     mpopt = mpoption(mpopt, 'out.all', 0);
0038 end
0039 
0040 %% define named indices into data matrices
0041 [PQ, PV, REF, NONE, BUS_I, BUS_TYPE, PD, QD, GS, BS, BUS_AREA, VM, ...
0042     VA, BASE_KV, ZONE, VMAX, VMIN, LAM_P, LAM_Q, MU_VMAX, MU_VMIN] = idx_bus;
0043 [GEN_BUS, PG, QG, QMAX, QMIN, VG, MBASE, GEN_STATUS, PMAX, PMIN, ...
0044     MU_PMAX, MU_PMIN, MU_QMAX, MU_QMIN, PC1, PC2, QC1MIN, QC1MAX, ...
0045     QC2MIN, QC2MAX, RAMP_AGC, RAMP_10, RAMP_30, RAMP_Q, APF] = idx_gen;
0046 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0047     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0048     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0049 
0050 %% load base case file
0051 mpc = loadcase(casefile);
0052 nb = size(mpc.bus, 1);
0053 
0054 t = 'PF - base case : ';
0055 mpopt0 = mpopt;
0056 r0 = runpf(mpc, mpopt);
0057 t_ok(r0.success, [t 'success']);
0058 
0059 for k = 1:length(pfalgs)
0060     mpopt = mpoption(mpopt0, 'pf.alg', pfalgs{k});
0061 
0062     t = sprintf('PF (%s) - constant power only : ', pfalgs{k});
0063     zip_w = [1 0 0];
0064     mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads', struct('pw', zip_w, 'qw', zip_w));
0065     r = runpf(mpc, mpopt);
0066     t_ok(r.success, [t 'success']);
0067     t_is(r.bus, r0.bus, 6, [t 'bus']);
0068     t_is(r.gen, r0.gen, 6, [t 'gen']);
0069     t_is(r.branch, r0.branch, 6, [t 'branch']);
0070 
0071     t = sprintf('PF (%s) - constant current only : ', pfalgs{k});
0072     zip_w = [0 1 0];
0073     mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads', struct('pw', zip_w, 'qw', zip_w));
0074     r = runpf(mpc, mpopt);
0075     [tPd, tQd] = total_load(r, 'bus', [], mpopt);
0076     % [tPd, tQd] = total_load(r, 'bus');
0077     mpc1 = mpc;
0078     mpc1.bus(:, PD) = mpc1.bus(:, PD) .* r.bus(:, VM);
0079     mpc1.bus(:, QD) = mpc1.bus(:, QD) .* r.bus(:, VM);
0080     r1 = runpf(mpc1, mpopt0);
0081     t_ok(r1.success, [t 'success']);
0082     t_ok(r.success, [t 'success']);
0083     t_is(r.bus(:, [VM VA]), r1.bus(:, [VM VA]), 6, [t 'bus voltages']);
0084     t_is(tPd, r1.bus(:, PD), 6, [t 'bus P loads']);
0085     t_is(tPd, [0; 0; 0; 0; 87.937521; 0; 98.674071; 0; 120.041555], 6, [t 'bus P loads']);
0086     t_is(tQd, r1.bus(:, QD), 6, [t 'bus Q loads']);
0087     t_is(tQd, [0; 0; 0; 0; 29.312507; 0; 34.535925; 0; 48.016622], 6, [t 'bus Q loads']);
0088     t_is(r.gen(:, [PG QG]), r1.gen(:, [PG QG]), 6, [t 'gen dispatches']);
0089     t_is(r.branch, r1.branch, 6, [t 'branch']);
0090 
0091     t = sprintf('PF (%s) - constant impedance only : ', pfalgs{k});
0092     mpc1 = mpc;
0093     mpc1.bus(:, GS) =  mpc1.bus(:, PD);
0094     mpc1.bus(:, BS) = -mpc1.bus(:, QD);
0095     mpc1.bus(:, PD) = 0;
0096     mpc1.bus(:, QD) = 0;
0097     r1 = runpf(mpc1, mpopt0);
0098     t_ok(r1.success, [t 'success']);
0099 
0100     zip_w = [0 0 1];
0101     mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads', struct('pw', zip_w, 'qw', []));
0102     r = runpf(mpc, mpopt);
0103     t_ok(r.success, [t 'success']);
0104     t_is(r.bus(:, [VM VA]), r1.bus(:, [VM VA]), 6, [t 'bus voltages']);
0105     t_is(r.gen(:, [PG QG]), r1.gen(:, [PG QG]), 6, [t 'gen dispatches']);
0106     t_is(r.branch, r1.branch, 6, [t 'branch']);
0107 
0108     t = sprintf('PF (%s) - combo ZIP loads : ', pfalgs{k});
0109     zip_w = [0.1 0.4 0.5];
0110     mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads.pw', zip_w);
0111     r = runpf(mpc, mpopt);
0112     [tPd, tQd] = total_load(r, 'bus', [], mpopt);
0113     mpc1 = mpc;
0114     Vm = r.bus(:, VM);
0115     scale = [Vm.^0 Vm Vm.^2] * zip_w';
0116     mpc1.bus(:, PD) = mpc1.bus(:, PD) .* scale;
0117     mpc1.bus(:, QD) = mpc1.bus(:, QD) .* scale;
0118     r1 = runpf(mpc1, mpopt0);
0119     t_ok(r1.success, [t 'success']);
0120 
0121     t_ok(r.success, [t 'success']);
0122     t_is(r.bus(:, [VM VA]), r1.bus(:, [VM VA]), 6, [t 'bus voltages']);
0123     t_is(tPd, r1.bus(:, PD), 6, [t 'bus P loads']);
0124     t_is(tPd, [0; 0; 0; 0; 87.205063; 0; 98.205247; 0; 118.314510], 6, [t 'bus P loads']);
0125     t_is(tQd, r1.bus(:, QD), 6, [t 'bus Q loads']);
0126     t_is(tQd, [0; 0; 0; 0; 29.068354; 0; 34.371836; 0; 47.325804], 6, [t 'bus Q loads']);
0127     t_is(r.gen(:, [PG QG]), r1.gen(:, [PG QG]), 6, [t 'gen dispatches']);
0128     t_is(r.branch, r1.branch, 6, [t 'branch']);
0129 
0130     t = sprintf('PF (%s) - combo ZIP loads + Q lims : ', pfalgs{k});
0131     zip_w = [0.1 0.4 0.5];
0132     mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads.pw', zip_w);
0133 
0134 %     mpopt2 = mpoption(mpopt, 'pf.enforce_q_lims', 1, 'verbose', 2, 'out.all', -1, 'out.gen', 1);
0135     mpopt2 = mpoption(mpopt, 'pf.enforce_q_lims', 1);
0136     mpc2 = mpc;
0137     mpc2.gen(2:3, QMAX) = 1;
0138     mpc2.gen(2:3, QMIN) = -1;
0139     r = runpf(mpc2, mpopt2);
0140     [tPd, tQd] = total_load(r, 'bus', [], mpopt2);
0141     mpc1 = mpc2;
0142     Vm = r.bus(:, VM);
0143     scale = [Vm.^0 Vm Vm.^2] * zip_w';
0144     mpc1.bus(:, PD) = mpc1.bus(:, PD) .* scale;
0145     mpc1.bus(:, QD) = mpc1.bus(:, QD) .* scale;
0146     mpopt1 = mpoption(mpopt0, 'pf.enforce_q_lims', 1);
0147     r1 = runpf(mpc1, mpopt1);
0148     t_ok(r1.success, [t 'success']);
0149 
0150     t_ok(r.success, [t 'success']);
0151     t_is(r.gen(2:3, QG), [-1; 1], 12, [t 'Qg(2:3)']);
0152     t_is(r.branch(4, QF), -1, 6, [t 'Qf(4)']);
0153     t_is(r.branch(7, QT), 1, 6, [t 'Qt(7)']);
0154     t_is(r.bus(:, [VM VA]), r1.bus(:, [VM VA]), 6, [t 'bus voltages']);
0155     t_is(tPd, r1.bus(:, PD), 6, [t 'bus P loads']);
0156     t_is(tPd, [0; 0; 0; 0; 85.9288588; 0; 95.2041458; 0; 115.8912983], 6, [t 'bus P loads']);
0157     t_is(tQd, r1.bus(:, QD), 6, [t 'bus Q loads']);
0158     t_is(tQd, [0; 0; 0; 0; 28.6429529; 0; 33.321451; 0; 46.3565193], 6, [t 'bus Q loads']);
0159     t_is(r.gen(:, [PG QG]), r1.gen(:, [PG QG]), 6, [t 'gen dispatches']);
0160     t_is(r.branch, r1.branch, 6, [t 'branch']);
0161 end
0162 
0163 t = 'OPF - base case : ';
0164 % mpopt0 = mpoption(mpopt0, 'verbose', 2);
0165 mpopt = mpopt0;
0166 r0 = runopf(mpc, mpopt);
0167 t_ok(r0.success, [t 'success']);
0168 
0169 t = 'OPF - constant power only : ';
0170 zip_w = [1 0 0];
0171 mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads.pw', zip_w);
0172 r = runopf(mpc, mpopt);
0173 t_ok(r.success, [t 'success']);
0174 t_is(r.f, r0.f, 6, [t 'f']);
0175 t_is(r.bus, r0.bus, 6, [t 'bus']);
0176 t_is(r.gen, r0.gen, 6, [t 'gen']);
0177 t_is(r.branch, r0.branch, 6, [t 'branch']);
0178 
0179 t = 'OPF - constant current only : ';
0180 zip_w = [0 1 0];
0181 mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads.pw', zip_w);
0182 r = runopf(mpc, mpopt);
0183 [tPd, tQd] = total_load(r, 'bus', [], mpopt);
0184 % [tPd, tQd] = total_load(r, 'bus');
0185 mpc1 = mpc;
0186 mpc1.bus(:, PD) = mpc1.bus(:, PD) .* r.bus(:, VM);
0187 mpc1.bus(:, QD) = mpc1.bus(:, QD) .* r.bus(:, VM);
0188 k = find(mpc1.bus(:, PD));  %% buses with non-zero loads
0189 mpc1.bus(k, VMAX) = r.bus(k, VM);
0190 mpc1.bus(k, VMIN) = r.bus(k, VM);
0191 r1 = runopf(mpc1, mpopt0);
0192 t_ok(r1.success, [t 'success']);
0193 t_ok(r.success, [t 'success']);
0194 t_is(r.f, r1.f, 4, [t 'f']);
0195 t_is(r.bus(:, [VM VA]), r1.bus(:, [VM VA]), 6, [t 'bus voltages']);
0196 t_is(tPd, r1.bus(:, PD), 6, [t 'bus P loads']);
0197 t_is(tPd, [0; 0; 0; 0; 81.708978; 0; 90.914185; 0; 112.500000], 6, [t 'bus P loads']);
0198 t_is(tQd, r1.bus(:, QD), 6, [t 'bus Q loads']);
0199 t_is(tQd, [0; 0; 0; 0; 27.236326; 0; 31.819965; 0; 45.000000], 6, [t 'bus Q loads']);
0200 t_is(r.gen(:, [PG QG]), r1.gen(:, [PG QG]), 6, [t 'gen dispatches']);
0201 t_is(r.branch(:, 1:ANGMAX), r1.branch(:, 1:ANGMAX), 6, [t 'branch']);
0202 
0203 t = 'OPF - constant impedance only : ';
0204 mpc1 = mpc;
0205 mpc1.bus(:, GS) =  mpc1.bus(:, PD);
0206 mpc1.bus(:, BS) = -mpc1.bus(:, QD);
0207 mpc1.bus(:, PD) = 0;
0208 mpc1.bus(:, QD) = 0;
0209 r1 = runopf(mpc1, mpopt0);
0210 t_ok(r1.success, [t 'success']);
0211 
0212 zip_w = [0 0 1];
0213 mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads.pw', zip_w);
0214 r = runopf(mpc, mpopt);
0215 t_ok(r.success, [t 'success']);
0216 t_is(r.f, r1.f, 6, [t 'f']);
0217 t_is(r.bus(:, [VM VA]), r1.bus(:, [VM VA]), 6, [t 'bus voltages']);
0218 t_is(r.gen(:, [PG QG]), r1.gen(:, [PG QG]), 6, [t 'gen dispatches']);
0219 t_is(r.branch, r1.branch, 6, [t 'branch']);
0220 
0221 t = 'OPF - combo ZIP loads : ';
0222 zip_w = [0.1 0.4 0.5];
0223 mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads.pw', zip_w);
0224 r = runopf(mpc, mpopt);
0225 [tPd, tQd] = total_load(r, 'bus', [], mpopt);
0226 mpc1 = mpc;
0227 Vm = r.bus(:, VM);
0228 scale = [Vm.^0 Vm Vm.^2] * zip_w';
0229 mpc1.bus(:, PD) = mpc1.bus(:, PD) .* scale;
0230 mpc1.bus(:, QD) = mpc1.bus(:, QD) .* scale;
0231 k = find(mpc1.bus(:, PD));  %% buses with non-zero loads
0232 mpc1.bus(k, VMAX) = r.bus(k, VM);
0233 mpc1.bus(k, VMIN) = r.bus(k, VM);
0234 r1 = runopf(mpc1, mpopt0);
0235 t_ok(r1.success, [t 'success']);
0236 
0237 t_ok(r.success, [t 'success']);
0238 t_is(r.f, r1.f, 4, [t 'f']);
0239 t_is(r.bus(:, [VM VA]), r1.bus(:, [VM VA]), 6, [t 'bus voltages']);
0240 t_is(tPd, r1.bus(:, PD), 6, [t 'bus P loads']);
0241 t_is(tPd, [0; 0; 0; 0; 78.909380; 0; 87.300375; 0; 108.125000], 6, [t 'bus P loads']);
0242 t_is(tQd, r1.bus(:, QD), 6, [t 'bus Q loads']);
0243 t_is(tQd, [0; 0; 0; 0; 26.303127; 0; 30.555131; 0; 43.250000], 6, [t 'bus Q loads']);
0244 t_is(r.gen(:, [PG QG]), r1.gen(:, [PG QG]), 6, [t 'gen dispatches']);
0245 t_is(r.branch(:, 1:ANGMAX), r1.branch(:, 1:ANGMAX), 6, [t 'branch']);
0246 
0247 
0248 mpopt =  mpopt0;
0249     plot_nose_curve = 0;
0250     mpopt = mpoption(mpopt, 'cpf.stop_at', 1);
0251     mpopt = mpoption(mpopt, 'cpf.step', 0.05);
0252     %mpopt = mpoption(mpopt, 'cpf.adapt_step', 1);
0253     %mpopt = mpoption(mpopt, 'cpf.adapt_step_tol', 2e-5);
0254     mpopt = mpoption(mpopt, 'cpf.plot.level', plot_nose_curve);
0255     %mpopt = mpoption(mpopt, 'cpf.plot.bus', 9);
0256 mpopt0 = mpopt;
0257 
0258     %% set up base and target cases
0259     mpcb = loadcase(casefile);
0260     mpct = mpcb;
0261     factor = 2.5 * 0.7;
0262     mpct.gen(:, [PG QG]) = mpct.gen(:, [PG QG]) * factor;
0263     mpct.bus(:, [PD QD]) = mpct.bus(:, [PD QD]) * factor;
0264 
0265     iterations = [
0266         20 20 20 20;
0267         21 22 22 22;
0268         21 22 22 22;
0269     ];
0270 
0271 for k = 1:length(cpfparm)
0272     mpopt0 = mpoption(mpopt0, 'cpf.parameterization', k);
0273 
0274     t = sprintf('CPF (%s) - base case : ', cpfparm{k});
0275     r0 = runcpf(mpcb, mpct, mpopt0);
0276     t_ok(r0.success, [t 'success']);
0277     t_is(r0.cpf.iterations, iterations(k, 1), 12, [t 'iterations']);
0278     t_is(r0.cpf.max_lam, 1, 12, [t 'max_lam']);
0279 
0280     mpopt = mpopt0;
0281     t = sprintf('CPF (%s) - constant power only : ', cpfparm{k});
0282     zip_w = [1 0 0];
0283     mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads.pw', zip_w);
0284     r = runcpf(mpcb, mpct, mpopt);
0285     t_ok(r.success, [t 'success']);
0286     t_is(r.cpf.iterations, r0.cpf.iterations, 12, [t 'iterations']);
0287     t_is(r.cpf.max_lam, r0.cpf.max_lam, 12, [t 'max_lam']);
0288     t_is(r.bus, r0.bus, 6, [t 'bus']);
0289     t_is(r.gen, r0.gen, 6, [t 'gen']);
0290     t_is(r.branch, r0.branch, 6, [t 'branch']);
0291 
0292     t = sprintf('CPF (%s) - constant current only : ', cpfparm{k});
0293     zip_w = [0 1 0];
0294     mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads.pw', zip_w);
0295     rb = runpf(mpcb, mpopt);
0296     r = runcpf(mpcb, mpct, mpopt);
0297     [tPd, tQd] = total_load(r, 'bus', [], mpopt);
0298     mpc1b = mpcb;
0299     mpc1b.bus(:, PD) = mpc1b.bus(:, PD) .* rb.bus(:, VM);
0300     mpc1b.bus(:, QD) = mpc1b.bus(:, QD) .* rb.bus(:, VM);
0301     mpc1t = mpct;
0302     mpc1t.bus(:, PD) = mpc1t.bus(:, PD) .* r.bus(:, VM);
0303     mpc1t.bus(:, QD) = mpc1t.bus(:, QD) .* r.bus(:, VM);
0304     r1 = runcpf(mpc1b, mpc1t, mpopt0);
0305     t_ok(r1.success, [t 'success']);
0306     t_is(r1.cpf.iterations, iterations(k,2), 12, [t 'iterations']);
0307     t_ok(r.success, [t 'success']);
0308     t_is(r.bus(:, [VM VA]), r1.bus(:, [VM VA]), 6, [t 'bus voltages']);
0309     t_is(tPd, r1.bus(:, PD), 6, [t 'bus P loads']);
0310     t_is(tPd, [0; 0; 0; 0; 143.494121; 0; 163.580859; 0; 191.911932], 6, [t 'bus P loads']);
0311     t_is(tQd, r1.bus(:, QD), 6, [t 'bus Q loads']);
0312     t_is(tQd, [0; 0; 0; 0; 47.831374; 0; 57.253301; 0; 76.764773], 6, [t 'bus Q loads']);
0313     t_is(r.gen(:, [PG QG]), r1.gen(:, [PG QG]), 6, [t 'gen dispatches']);
0314     t_is(r.branch, r1.branch, 6, [t 'branch']);
0315 
0316     t = sprintf('CPF (%s) - constant impedance only : ', cpfparm{k});
0317     zip_w = [0 0 1];
0318     mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads.pw', zip_w);
0319     rb = runpf(mpcb, mpopt);
0320     r = runcpf(mpcb, mpct, mpopt);
0321     [tPd, tQd] = total_load(r, 'bus', [], mpopt);
0322     mpc1b = mpcb;
0323     mpc1b.bus(:, PD) = mpc1b.bus(:, PD) .* rb.bus(:, VM).^2;
0324     mpc1b.bus(:, QD) = mpc1b.bus(:, QD) .* rb.bus(:, VM).^2;
0325     mpc1t = mpct;
0326     mpc1t.bus(:, PD) = mpc1t.bus(:, PD) .* r.bus(:, VM).^2;
0327     mpc1t.bus(:, QD) = mpc1t.bus(:, QD) .* r.bus(:, VM).^2;
0328     r1 = runcpf(mpc1b, mpc1t, mpopt0);
0329     t_ok(r1.success, [t 'success']);
0330     t_is(r1.cpf.iterations, iterations(k,3), 12, [t 'iterations']);
0331     t_ok(r.success, [t 'success']);
0332     t_is(r.bus(:, [VM VA]), r1.bus(:, [VM VA]), 6, [t 'bus voltages']);
0333     t_is(tPd, r1.bus(:, PD), 6, [t 'bus P loads']);
0334     t_is(tPd, [0; 0; 0; 0; 132.977098; 0; 154.969557; 0; 172.896177], 6, [t 'bus P loads']);
0335     t_is(tQd, r1.bus(:, QD), 6, [t 'bus Q loads']);
0336     t_is(tQd, [0; 0; 0; 0; 44.325699; 0; 54.239345; 0; 69.158471], 6, [t 'bus Q loads']);
0337     t_is(r.gen(:, [PG QG]), r1.gen(:, [PG QG]), 6, [t 'gen dispatches']);
0338     t_is(r.branch, r1.branch, 6, [t 'branch']);
0339 
0340     t = sprintf('CPF (%s) - combo ZIP loads : ', cpfparm{k});
0341     zip_w = [0.1 0.4 0.5];
0342     mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads.pw', zip_w);
0343     rb = runpf(mpcb, mpopt);
0344     r = runcpf(mpcb, mpct, mpopt);
0345     [tPd, tQd] = total_load(r, 'bus', [], mpopt);
0346     Vm = abs(rb.bus(:, VM));
0347     scaleb = [Vm.^0 Vm Vm.^2] * zip_w';
0348     Vm = abs(r.bus(:, VM));
0349     scalet = [Vm.^0 Vm Vm.^2] * zip_w';
0350     mpc1b = mpcb;
0351     mpc1b.bus(:, PD) = mpc1b.bus(:, PD) .* scaleb;
0352     mpc1b.bus(:, QD) = mpc1b.bus(:, QD) .* scaleb;
0353     mpc1t = mpct;
0354     mpc1t.bus(:, PD) = mpc1t.bus(:, PD) .* scalet;
0355     mpc1t.bus(:, QD) = mpc1t.bus(:, QD) .* scalet;
0356     r1 = runcpf(mpc1b, mpc1t, mpopt0);
0357     t_ok(r1.success, [t 'success']);
0358     t_is(r1.cpf.iterations, iterations(k,4), 12, [t 'iterations']);
0359     t_ok(r.success, [t 'success']);
0360     t_is(r.bus(:, [VM VA]), r1.bus(:, [VM VA]), 6, [t 'bus voltages']);
0361     t_is(tPd, r1.bus(:, PD), 6, [t 'bus P loads']);
0362     t_is(tPd, [0; 0; 0; 0; 139.202168; 0; 160.005294; 0; 184.209518], 6, [t 'bus P loads']);
0363     t_is(tQd, r1.bus(:, QD), 6, [t 'bus Q loads']);
0364     t_is(tQd, [0; 0; 0; 0; 46.400723; 0; 56.001853; 0; 73.683807], 6, [t 'bus Q loads']);
0365     t_is(r.gen(:, [PG QG]), r1.gen(:, [PG QG]), 6, [t 'gen dispatches']);
0366     t_is(r.branch, r1.branch, 6, [t 'branch']);
0367 end
0368 
0369 
0370 %% set up base and target cases
0371     mpcb = loadcase(casefile);
0372     mpct = mpcb;
0373     factor = 2.5 * 0.6;
0374     mpct.gen(:, [PG QG]) = mpct.gen(:, [PG QG]) * factor;
0375     mpct.bus(:, [PD QD]) = mpct.bus(:, [PD QD]) * factor;
0376     mpct.gen(:, QMAX) = 75;
0377     mpcb.gen(:, QMAX) = 75;
0378 
0379     %% with Q limits enforced
0380     iterations = [
0381         22 20 20 20;
0382         22 21 21 21;
0383         22 21 21 21;
0384     ];
0385 
0386 for k = 1:length(cpfparm)
0387     mpopt0 = mpoption(mpopt0, 'cpf.parameterization', k, 'cpf.enforce_q_lims', 1);
0388 
0389     t = sprintf('CPF (%s) w/Q lims - base case : ', cpfparm{k});
0390     r0 = runcpf(mpcb, mpct, mpopt0);
0391     t_ok(r0.success, [t 'success']);
0392     t_is(r0.cpf.iterations, iterations(k, 1), 12, [t 'iterations']);
0393     t_is(r0.cpf.max_lam, 1, 12, [t 'max_lam']);
0394 
0395     mpopt = mpopt0;
0396     t = sprintf('CPF (%s) w/Q lims - constant power only : ', cpfparm{k});
0397     zip_w = [1 0 0];
0398     mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads.pw', zip_w);
0399     r = runcpf(mpcb, mpct, mpopt);
0400     t_ok(r.success, [t 'success']);
0401     t_is(r.cpf.iterations, r0.cpf.iterations, 12, [t 'iterations']);
0402     t_is(r.cpf.max_lam, r0.cpf.max_lam, 12, [t 'max_lam']);
0403     t_is(r.bus, r0.bus, 6, [t 'bus']);
0404     t_is(r.gen, r0.gen, 6, [t 'gen']);
0405     t_is(r.branch, r0.branch, 6, [t 'branch']);
0406 
0407     t = sprintf('CPF (%s) w/Q lims - constant current only : ', cpfparm{k});
0408     zip_w = [0 1 0];
0409     mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads.pw', zip_w);
0410     rb = runpf(mpcb, mpopt);
0411     r = runcpf(mpcb, mpct, mpopt);
0412     [tPd, tQd] = total_load(r, 'bus', [], mpopt);
0413     mpc1b = mpcb;
0414     mpc1b.bus(:, PD) = mpc1b.bus(:, PD) .* rb.bus(:, VM);
0415     mpc1b.bus(:, QD) = mpc1b.bus(:, QD) .* rb.bus(:, VM);
0416     mpc1t = mpct;
0417     mpc1t.bus(:, PD) = mpc1t.bus(:, PD) .* r.bus(:, VM);
0418     mpc1t.bus(:, QD) = mpc1t.bus(:, QD) .* r.bus(:, VM);
0419     r1 = runcpf(mpc1b, mpc1t, mpopt0);
0420     t_ok(r1.success, [t 'success']);
0421     t_is(r1.cpf.iterations, iterations(k,2), 12, [t 'iterations']);
0422     t_ok(r.success, [t 'success']);
0423     t_is(r.bus(:, [VM VA]), r1.bus(:, [VM VA]), 6, [t 'bus voltages']);
0424     t_is(tPd, r1.bus(:, PD), 6, [t 'bus P loads']);
0425     t_is(tPd, [0; 0; 0; 0; 126.246737; 0; 143.114682; 0; 170.231775], 6, [t 'bus P loads']);
0426     t_is(tQd, r1.bus(:, QD), 6, [t 'bus Q loads']);
0427     t_is(tQd, [0; 0; 0; 0; 42.082245; 0; 50.090138; 0; 68.092710], 6, [t 'bus Q loads']);
0428     t_is(r.gen(:, [PG QG]), r1.gen(:, [PG QG]), 6, [t 'gen dispatches']);
0429     t_is(r.branch, r1.branch, 6, [t 'branch']);
0430 
0431     t = sprintf('CPF (%s) w/Q lims - constant impedance only : ', cpfparm{k});
0432     zip_w = [0 0 1];
0433     mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads.pw', zip_w);
0434     rb = runpf(mpcb, mpopt);
0435     r = runcpf(mpcb, mpct, mpopt);
0436     [tPd, tQd] = total_load(r, 'bus', [], mpopt);
0437     mpc1b = mpcb;
0438     mpc1b.bus(:, PD) = mpc1b.bus(:, PD) .* rb.bus(:, VM).^2;
0439     mpc1b.bus(:, QD) = mpc1b.bus(:, QD) .* rb.bus(:, VM).^2;
0440     mpc1t = mpct;
0441     mpc1t.bus(:, PD) = mpc1t.bus(:, PD) .* r.bus(:, VM).^2;
0442     mpc1t.bus(:, QD) = mpc1t.bus(:, QD) .* r.bus(:, VM).^2;
0443     r1 = runcpf(mpc1b, mpc1t, mpopt0);
0444     t_ok(r1.success, [t 'success']);
0445     t_is(r1.cpf.iterations, iterations(k,3), 12, [t 'iterations']);
0446     t_ok(r.success, [t 'success']);
0447     t_is(r.bus(:, [VM VA]), r1.bus(:, [VM VA]), 6, [t 'bus voltages']);
0448     t_is(tPd, r1.bus(:, PD), 6, [t 'bus P loads']);
0449     t_is(tPd, [0; 0; 0; 0; 119.354069; 0; 137.673893; 0; 157.193899], 6, [t 'bus P loads']);
0450     t_is(tQd, r1.bus(:, QD), 6, [t 'bus Q loads']);
0451     t_is(tQd, [0; 0; 0; 0; 39.784689; 0; 48.185862; 0; 62.877559], 6, [t 'bus Q loads']);
0452     t_is(r.gen(:, [PG QG]), r1.gen(:, [PG QG]), 6, [t 'gen dispatches']);
0453     t_is(r.branch, r1.branch, 6, [t 'branch']);
0454 
0455     t = sprintf('CPF (%s) w/Q lims - combo ZIP loads : ', cpfparm{k});
0456     zip_w = [0.1 0.4 0.5];
0457     mpopt = mpoption(mpopt, 'exp.sys_wide_zip_loads.pw', zip_w);
0458     rb = runpf(mpcb, mpopt);
0459     r = runcpf(mpcb, mpct, mpopt);
0460     [tPd, tQd] = total_load(r, 'bus', [], mpopt);
0461     Vm = abs(rb.bus(:, VM));
0462     scaleb = [Vm.^0 Vm Vm.^2] * zip_w';
0463     Vm = abs(r.bus(:, VM));
0464     scalet = [Vm.^0 Vm Vm.^2] * zip_w';
0465     mpc1b = mpcb;
0466     mpc1b.bus(:, PD) = mpc1b.bus(:, PD) .* scaleb;
0467     mpc1b.bus(:, QD) = mpc1b.bus(:, QD) .* scaleb;
0468     mpc1t = mpct;
0469     mpc1t.bus(:, PD) = mpc1t.bus(:, PD) .* scalet;
0470     mpc1t.bus(:, QD) = mpc1t.bus(:, QD) .* scalet;
0471     r1 = runcpf(mpc1b, mpc1t, mpopt0);
0472     t_ok(r1.success, [t 'success']);
0473     t_is(r1.cpf.iterations, iterations(k,4), 12, [t 'iterations']);
0474     t_ok(r.success, [t 'success']);
0475     t_is(r.bus(:, [VM VA]), r1.bus(:, [VM VA]), 6, [t 'bus voltages']);
0476     t_is(tPd, r1.bus(:, PD), 6, [t 'bus P loads']);
0477     t_is(tPd, [0; 0; 0; 0; 123.418124; 0; 140.853612; 0; 164.913503], 6, [t 'bus P loads']);
0478     t_is(tQd, r1.bus(:, QD), 6, [t 'bus Q loads']);
0479     t_is(tQd, [0; 0; 0; 0; 41.1393746; 0; 49.298764; 0; 65.9654012], 6, [t 'bus Q loads']);
0480     t_is(r.gen(:, [PG QG]), r1.gen(:, [PG QG]), 6, [t 'gen dispatches']);
0481     t_is(r.branch, r1.branch, 6, [t 'branch']);
0482 end
0483 
0484 t_end;

Generated on Fri 09-Oct-2020 11:21:31 by m2html © 2005