Home > matpower7.1 > lib > calc_v_y_sum.m

calc_v_y_sum

PURPOSE ^

CALC_V_Y_SUM Solves the power flow using the admittance summation method.

SYNOPSIS ^

function [V, Qpv, Sf, St, Sslack, iter, success] = calc_v_y_sum(Vslack,nb,nl,f,Zb,Ybf,Ybt,Yd,Sd,pv,Pg,Vg,mpopt)

DESCRIPTION ^

CALC_V_Y_SUM  Solves the power flow using the admittance summation method.

   [V, Qpv, Sf, St, Sslack, iter, success] = calc_v_y_sum(Vslack,nb,nl,f,Zb,Ybf,Ybt,Yd,Sd,pv,Pg,Vg,tol,iter_max)

   Solves for bus voltages, generator reactive power, branch active and
   reactive power flows and slack bus active and reactive power. The input
   data consist of slack bus voltage, vector "from bus" indices, branch
   impedance and shunt admittance, vector of bus shunt admittances and
   load demand, as well as vectors with indicies of PV buses with their
   specified voltages and active powers. It is assumed that the branches
   are ordered using the principle of oriented ordering: indicies of
   sending nodes are smaller then the indicies of the receiving nodes. The
   branch index is equal to the index of their receiving node. Branch
   addmittances are added in Yd and treated as constant admittance bus
   loads. The applied method is admittance summation taken from:
   Dragoslav Rajičić, Rubin Taleski, Two novel methods for radial and
   weakly meshed network analysis, Electric Power Systems Research,
   Volume 48, Issue 2, 15 December 1998, Pages 79-87
   https://doi.org/10.1016/S0378-7796(98)00067-4

   See also RADIAL_PF.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [V, Qpv, Sf, St, Sslack, iter, success] = calc_v_y_sum(Vslack,nb,nl,f,Zb,Ybf,Ybt,Yd,Sd,pv,Pg,Vg,mpopt)
0002 %CALC_V_Y_SUM  Solves the power flow using the admittance summation method.
0003 %
0004 %   [V, Qpv, Sf, St, Sslack, iter, success] = calc_v_y_sum(Vslack,nb,nl,f,Zb,Ybf,Ybt,Yd,Sd,pv,Pg,Vg,tol,iter_max)
0005 %
0006 %   Solves for bus voltages, generator reactive power, branch active and
0007 %   reactive power flows and slack bus active and reactive power. The input
0008 %   data consist of slack bus voltage, vector "from bus" indices, branch
0009 %   impedance and shunt admittance, vector of bus shunt admittances and
0010 %   load demand, as well as vectors with indicies of PV buses with their
0011 %   specified voltages and active powers. It is assumed that the branches
0012 %   are ordered using the principle of oriented ordering: indicies of
0013 %   sending nodes are smaller then the indicies of the receiving nodes. The
0014 %   branch index is equal to the index of their receiving node. Branch
0015 %   addmittances are added in Yd and treated as constant admittance bus
0016 %   loads. The applied method is admittance summation taken from:
0017 %   Dragoslav Rajičić, Rubin Taleski, Two novel methods for radial and
0018 %   weakly meshed network analysis, Electric Power Systems Research,
0019 %   Volume 48, Issue 2, 15 December 1998, Pages 79-87
0020 %   https://doi.org/10.1016/S0378-7796(98)00067-4
0021 %
0022 %   See also RADIAL_PF.
0023 
0024 %% initialize
0025 tol      = mpopt.pf.tol;
0026 iter_max = mpopt.pf.radial.max_it;
0027 vcorr    = mpopt.pf.radial.vcorr == 1;
0028 Sd(pv) = Sd(pv) - Pg;
0029 V = Vslack * ones(nb,1);
0030 Vold = V;
0031 iter = 0;
0032 success = 0;
0033 % ZIP load model
0034 pw = mpopt.exp.sys_wide_zip_loads.pw;
0035 qw = mpopt.exp.sys_wide_zip_loads.qw;
0036 if isempty(pw)
0037     pw = [1 0 0];
0038 end
0039 if isempty(qw)
0040     qw = pw;
0041 end
0042 Sdz = real(Sd) * pw(3) + 1j * imag(Sd) * qw(3); % constant impedance
0043 Sdi = real(Sd) * pw(2) + 1j * imag(Sd) * qw(2); % constant current
0044 Sdp = real(Sd) * pw(1) + 1j * imag(Sd) * qw(1); % constant power
0045 % Add artificial branch at the top of the branch list, so that the branch
0046 % index for other branches is equal to the index of their receiving node.
0047  f = [0; f];
0048 Zb = [0; Zb];
0049 nl = nl + 1;
0050 %% make Zpv matrix, for calculation of the PV generators reactive powers
0051 if ~isempty(pv)
0052     Zpv = make_zpv(pv,nb,nl,f,Zb,Yd);
0053     Bpv = (imag(Zpv))^-1;
0054 end
0055 npv = length(pv);
0056 Qpv = zeros(npv,1);
0057 %% do backward-forward iterations
0058 if mpopt.verbose > 1
0059     fprintf('\n it    max V mismatch (p.u.)');
0060     fprintf('\n----  ----------------------');
0061 end
0062 Ye = conj(Sdz) + Yd;
0063 D = zeros(nl,1);
0064 for k = nl:-1:2
0065     D(k) = 1 / (1 + Zb(k)*Ye(k));
0066     i = f(k);
0067     Ye(i) = Ye(i) + D(k)*Ye(k);
0068 end
0069 while success == 0 && iter < iter_max
0070     iter = iter + 1;
0071     % calculate load demand using actual voltages (Sdz and Yd are already
0072     % included in Ye)
0073     Vm = abs(V);
0074     S = Sdp + Sdi.*Vm;
0075     % backward sweep
0076     Je = conj(S./V);
0077     for k = nl:-1:2
0078         i = f(k);
0079         Je(i) = Je(i) + D(k)*Je(k);
0080     end
0081     % forward sweep
0082     for k = 2:nl
0083         i = f(k);
0084         V(k) = D(k) * (V(i)-Zb(k)*Je(k));
0085     end
0086     % check for convergence
0087     DU = abs(V - Vold);
0088     DU(isnan(DU)) = inf;
0089     if mpopt.verbose > 1
0090         fprintf('\n%3d        %10.3e', iter, max(DU));
0091     end
0092     if max(DU) > tol
0093         Vold = V;
0094         % update PV generators reactive powers
0095         if ~isempty(pv)
0096             DE = (Vg./abs(V(pv))-1).*real(V(pv)); % Rajicic (VCPF)
0097             DD = Bpv * DE;
0098             if vcorr
0099                 DC = DD .* imag(V(pv))./real(V(pv));
0100                 V_corr = make_vcorr(DC+1j*DD,pv,nb,nl,f,Zb);
0101                 V = V + V_corr;
0102             end
0103             DQ = DD .* abs(V(pv)).^2 ./ real(V(pv));
0104             Qpv = Qpv + DQ;
0105             Sdp(pv) = Sdp(pv) - 1j*DQ;
0106         end
0107     else
0108         success = 1;
0109     end
0110 end
0111 if mpopt.verbose
0112     if success
0113         fprintf('\nAdmittance summation converged in %d iterations.\n', iter);
0114     else
0115         fprintf('\nAdmittance summation did not converge in %d iterations.\n', iter);
0116     end
0117 end
0118 %% calculate branch flows
0119 % take out the first artificial branch
0120 Sslack = V(1)*conj(Je(1)) + conj(Ye(1))*abs(V(1))^2;
0121 f = f(2:end);
0122 Zb = Zb(2:end);
0123 % calculate branch flows
0124 I = (V(f) - V(2:end))./Zb;
0125 Sf = V(f).*conj(I) + conj(Ybf) .* abs(V(f)).^2;
0126 St = V(2:end).*conj(I) - conj(Ybt) .* abs(V(2:end)).^2;

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