0001 function t_pf_ac(quiet)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012 if nargin < 1
0013 quiet = 0;
0014 end
0015
0016 AC_alg = {'NR', 'NR-SP', 'NR-SC', 'NR-SH', 'NR-IP', 'NR-IC', 'NR-IH', 'FDXB', 'FDBX', 'GS'};
0017 AC_name = {
0018 'Newton (default, power-polar)',
0019 'Newton (power-polar)',
0020 'Newton (power-cartesian)',
0021 'Newton (power-hybrid)',
0022 'Newton (current-polar)',
0023 'Newton (current-cartesian)',
0024 'Newton (current-hybrid)',
0025 'Fast Decoupled (XB)',
0026 'Fast Decoupled (BX)',
0027 'Gauss-Seidel'
0028 };
0029
0030
0031
0032 t_begin(length(AC_alg)*48, quiet);
0033
0034 casefile = 't_case9_pf';
0035 if quiet
0036 verbose = 0;
0037 else
0038 verbose = 1;
0039 end
0040 if have_feature('octave')
0041 if have_feature('octave', 'vnum') >= 4
0042 file_in_path_warn_id = 'Octave:data-file-in-path';
0043 else
0044 file_in_path_warn_id = 'Octave:load-file-in-path';
0045 end
0046 s1 = warning('query', file_in_path_warn_id);
0047 warning('off', file_in_path_warn_id);
0048 end
0049 mpopt0 = mpoption('out.all', 0, 'pf.tol', 1e-9, 'verbose', 0);
0050 mpopt0 = mpoption(mpopt0, 'verbose', verbose);
0051
0052
0053 [PQ, PV, REF, NONE, BUS_I, BUS_TYPE, PD, QD, GS, BS, BUS_AREA, VM, ...
0054 VA, BASE_KV, ZONE, VMAX, VMIN, LAM_P, LAM_Q, MU_VMAX, MU_VMIN] = idx_bus;
0055 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0056 TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0057 ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0058 [GEN_BUS, PG, QG, QMAX, QMIN, VG, MBASE, GEN_STATUS, PMAX, PMIN, ...
0059 MU_PMAX, MU_PMIN, MU_QMAX, MU_QMIN, PC1, PC2, QC1MIN, QC1MAX, ...
0060 QC2MIN, QC2MAX, RAMP_AGC, RAMP_10, RAMP_30, RAMP_Q, APF] = idx_gen;
0061
0062
0063 mpc0 = loadcase(casefile);
0064 mpc0.gen(1, PG) = 60;
0065 mpc0.gen(1, [PMIN PMAX QMIN QMAX PG QG]) = mpc0.gen(1, [PMIN PMAX QMIN QMAX PG QG]) / 2;
0066 mpc0.gen = [mpc0.gen(1, :); mpc0.gen];
0067 mpc1 = mpc0;
0068 mpc = mpc0;
0069 nb = size(mpc.bus, 1);
0070 mpc1.bus(:, BUS_I) = mpc1.bus(:, BUS_I) + nb;
0071 mpc1.branch(:, F_BUS) = mpc1.branch(:, F_BUS) + nb;
0072 mpc1.branch(:, T_BUS) = mpc1.branch(:, T_BUS) + nb;
0073 mpc1.gen(:, GEN_BUS) = mpc1.gen(:, GEN_BUS) + nb;
0074 mpc.bus = [mpc.bus; mpc1.bus];
0075 mpc.branch = [mpc.branch; mpc1.branch];
0076 mpc.gen = [mpc.gen; mpc1.gen];
0077 mpc1 = mpc;
0078
0079
0080
0081 load soln9_pf;
0082 soln_vg = load('soln9_pf_vg');
0083
0084
0085 for k = 1:length(AC_alg)
0086 t = sprintf('AC PF - %s : ', AC_name{k});
0087 mpopt = mpoption(mpopt0, 'pf.alg', AC_alg{k});
0088 [baseMVA, bus, gen, branch, success, et] = runpf(casefile, mpopt);
0089 t_ok(success, [t 'success']);
0090 t_is(bus, bus_soln, 6, [t 'bus']);
0091 t_is(gen, gen_soln, 6, [t 'gen']);
0092 t_is(branch, branch_soln, 6, [t 'branch']);
0093
0094 r = runpf(casefile, mpopt);
0095 t_ok(r.success, [t 'success']);
0096 t_is(r.bus, bus_soln, 6, [t 'bus']);
0097 t_is(r.gen, gen_soln, 6, [t 'gen']);
0098 t_is(r.branch, branch_soln, 6, [t 'branch']);
0099
0100
0101 t = sprintf('%s - Vg ~= 1 : ', AC_alg{k});
0102 mpopt = mpoption(mpopt, 'verbose', 0);
0103 mpc = loadcase(casefile);
0104 mpc.gen(:, VG) = [1.04; 1.025; 1.025];
0105 r = runpf(mpc, mpopt);
0106 t_ok(r.success, [t 'success']);
0107 t_is(r.bus, soln_vg.bus_soln, 6, [t 'bus']);
0108 t_is(r.gen, soln_vg.gen_soln, 6, [t 'gen']);
0109 t_is(r.branch, soln_vg.branch_soln, 6, [t 'branch']);
0110
0111
0112 t = sprintf('%s - check Qg : ', AC_alg{k});
0113 mpc = loadcase(casefile);
0114 mpc.gen(1, [QMIN QMAX]) = [20 20];
0115 r = runpf(mpc, mpopt);
0116 t_ok(r.success, [t 'success']);
0117 t_is(r.gen(1, QG), 24.07, 2, [t 'single gen, Qmin = Qmax']);
0118
0119 mpc.gen = [mpc.gen(1, :); mpc.gen];
0120 mpc.gen(1, [QMIN QMAX]) = [10 10];
0121 mpc.gen(2, [QMIN QMAX]) = [0 50];
0122 r = runpf(mpc, mpopt);
0123 t_ok(r.success, [t 'success']);
0124 t_is(r.gen(1:2, QG), [10; 14.07], 2, [t '2 gens, Qmin = Qmax for one']);
0125
0126 mpc.gen(1, [QMIN QMAX]) = [10 10];
0127 mpc.gen(2, [QMIN QMAX]) = [-50 -50];
0128 r = runpf(mpc, mpopt);
0129 t_ok(r.success, [t 'success']);
0130 t_is(r.gen(1:2, QG), [42.03; -17.97], 2, [t '2 gens, Qmin = Qmax for both']);
0131
0132 mpc.gen(1, [QMIN QMAX]) = [0 50];
0133 mpc.gen(2, [QMIN QMAX]) = [0 100];
0134 r = runpf(mpc, mpopt);
0135 t_ok(r.success, [t 'success']);
0136 t_is(r.gen(1:2, QG), [8.02; 16.05], 2, [t '2 gens, proportional']);
0137
0138 mpc.gen(1, [QMIN QMAX]) = [-50 0];
0139 mpc.gen(2, [QMIN QMAX]) = [50 150];
0140 r = runpf(mpc, mpopt);
0141 t_ok(r.success, [t 'success']);
0142 t_is(r.gen(1:2, QG), [-50+8.02; 50+16.05], 2, [t '2 gens, proportional']);
0143
0144 mpc.gen(1, [QMIN QMAX]) = [-50 Inf];
0145 mpc.gen(2, [QMIN QMAX]) = [50 150];
0146 r = runpf(mpc, mpopt);
0147 t_ok(r.success, [t 'success']);
0148 t_is(r.gen(1:2, QG), [-31.61; 55.68], 2, [t '2 gens, one infinite range']);
0149
0150 mpc.gen(1, [QMIN QMAX]) = [-50 Inf];
0151 mpc.gen(2, [QMIN QMAX]) = [50 Inf];
0152 r = runpf(mpc, mpopt);
0153 t_ok(r.success, [t 'success']);
0154 t_is(r.gen(1:2, QG), [-33.12; 57.18], 2, [t '2 gens, both infinite range']);
0155
0156 mpc.gen(1, [QMIN QMAX]) = [-50 Inf];
0157 mpc.gen(2, [QMIN QMAX]) = [-Inf 150];
0158 r = runpf(mpc, mpopt);
0159 t_ok(r.success, [t 'success']);
0160 t_is(r.gen(1:2, QG), [76.07; -52], 2, [t '2 gens, both infinite range']);
0161
0162 mpc.gen(1, [QMIN QMAX]) = [-Inf Inf];
0163 mpc.gen(2, [QMIN QMAX]) = [-Inf Inf];
0164 r = runpf(mpc, mpopt);
0165 t_ok(r.success, [t 'success']);
0166 t_is(r.gen(1:2, QG), [12.03; 12.03], 2, [t '2 gens, both infinite range']);
0167
0168 t = sprintf('%s - reactive generation allocation : ', AC_alg{k});
0169 mpc = loadcase(casefile);
0170
0171
0172 mpc.gen = [
0173 1 0 0 300 -300 1 100 1 250 10 0 0 0 0 0 0 0 0 0 0 0;
0174 2 54 0 0 -5 1 100 1 300 10 0 0 0 0 0 0 0 0 0 0 0;
0175 2 54 0 5 -5 1 100 1 300 10 0 0 0 0 0 0 0 0 0 0 0;
0176 2 55 0 25 10 1 100 1 300 10 0 0 0 0 0 0 0 0 0 0 0;
0177 30 25 1 300 -300 1 100 1 270 10 0 0 0 0 0 0 0 0 0 0 0;
0178 30 30 2 300 -300 1 100 1 270 10 0 0 0 0 0 0 0 0 0 0 0;
0179 30 30 -3 300 -300 1 100 1 270 10 0 0 0 0 0 0 0 0 0 0 0;
0180 ];
0181 mpc.bus(3, BUS_TYPE) = PQ;
0182 r = runpf(mpc, mpopt);
0183 t_ok(r.success, [t 'success']);
0184 t_is(r.gen(2:4, QG), [-5; -5; 10] + [1; 2; 3]*1.989129794, 7, [t 'PV bus']);
0185 t_is(r.gen(5:7, QG), [1; 2; -3], 8, [t 'PQ bus']);
0186
0187
0188 t = sprintf('%s - network w/islands : AC PF : ', AC_alg{k});
0189 mpc = mpc1;
0190 r = runpf(mpc, mpopt);
0191 t_ok(r.success, [t 'success']);
0192 t_is(r.bus( 1:9, VA), bus_soln(:, VA), 7, [t 'voltage angles 1']);
0193 t_is(r.bus(10:18, VA), bus_soln(:, VA), 7, [t 'voltage angles 2']);
0194 Pg = [gen_soln(1, PG)-30; 30; gen_soln(2:3, PG)];
0195 t_is(r.gen(1:4, PG), Pg, 6, [t 'active power generation 1']);
0196 t_is(r.gen(5:8, PG), Pg, 6, [t 'active power generation 2']);
0197
0198 t = sprintf('%s - all buses isolated : ', AC_alg{k});
0199 mpc.bus(:, BUS_TYPE) = NONE;
0200 try
0201 r = runpf(mpc, mpopt);
0202 t_is(r.success, 0, 12, [t 'success = 0']);
0203 catch
0204 t_ok(0, [t 'unexpected fatal error']);
0205 end
0206
0207
0208 t = sprintf('%s - pf.enforce_q_lims == 0 : ', AC_alg{k});
0209 mpc = loadcase('case14');
0210 mpc.gen(1, QMIN) = -10;
0211 mpc.gen(:, QMAX) = [10; 30; 29; 15; 15];
0212 bt0 = mpc.bus(:, BUS_TYPE);
0213 bt = bt0;
0214 mpopt = mpoption(mpopt, 'pf.enforce_q_lims', 0);
0215 r = runpf(mpc, mpopt);
0216 t_ok(r.success, [t 'success']);
0217 t_is(r.bus(:, BUS_TYPE), bt, 12, [t 'bus type']);
0218 t_is(r.gen(:, QG), [-16.549300542; 43.557100134; 25.075348495; 12.730944405; 17.623451366], 6, [t 'Qg']);
0219
0220 t = sprintf('%s - pf.enforce_q_lims == 1 : ', AC_alg{k});
0221 mpopt = mpoption(mpopt, 'pf.enforce_q_lims', 1);
0222 r = runpf(mpc, mpopt);
0223 bt = bt0;
0224 bt([1 2 3 8]) = [PQ PQ REF PQ];
0225 t_is(r.success, 0, 12, [t 'success = 0']);
0226 t_is(r.bus(:, BUS_TYPE), bt, 12, [t 'bus type']);
0227 t_is(r.gen(:, QG), [-10; 30; 31.608422873; 16.420423190; 15], 4, [t 'Qg']);
0228
0229 t = sprintf('%s - pf.enforce_q_lims == 2 : ', AC_alg{k});
0230 mpopt = mpoption(mpopt, 'pf.enforce_q_lims', 2);
0231 r = runpf(mpc, mpopt);
0232 bt = bt0;
0233 bt([1 2 3 6 8]) = [REF PQ PQ PQ PQ];
0234 t_ok(r.success, [t 'success']);
0235 t_is(r.bus(:, BUS_TYPE), bt, 12, [t 'bus type']);
0236 t_is(r.gen(:, QG), [-6.30936644; 30; 29; 15; 15], 6, [t 'Qg']);
0237 end
0238
0239 t_end;
0240
0241 if have_feature('octave')
0242 warning(s1.state, file_in_path_warn_id);
0243 end