Home > matpower7.1 > lib > calc_v_i_sum.m

calc_v_i_sum

PURPOSE ^

CALC_V_I_SUM Solves the power flow using the current summation method.

SYNOPSIS ^

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

DESCRIPTION ^

CALC_V_I_SUM  Solves the power flow using the current summation method.

   [V, Qpv, Sf, St, Sslack, iter, success] = calc_v_i_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 current summation taken from:
   D. Shirmohammadi, H. W. Hong, A. Semlyen and G. X. Luo,
   "A compensation-based power flow method for weakly meshed distribution
   and transmission networks," IEEE Transactions on Power Systems,
   vol. 3, no. 2, pp. 753-762, May 1988. https://doi.org/10.1109/59.192932
   and
   G. X. Luo and A. Semlyen, "Efficient load flow for large weakly meshed
   networks," IEEE Transactions on Power Systems, vol. 5, no. 4,
   pp. 1309-1316, Nov 1990. https://doi.org/10.1109/59.99382

   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_i_sum(Vslack,nb,nl,f,Zb,Ybf,Ybt,Yd,Sd,pv,Pg,Vg,mpopt)
0002 %CALC_V_I_SUM  Solves the power flow using the current summation method.
0003 %
0004 %   [V, Qpv, Sf, St, Sslack, iter, success] = calc_v_i_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 current summation taken from:
0017 %   D. Shirmohammadi, H. W. Hong, A. Semlyen and G. X. Luo,
0018 %   "A compensation-based power flow method for weakly meshed distribution
0019 %   and transmission networks," IEEE Transactions on Power Systems,
0020 %   vol. 3, no. 2, pp. 753-762, May 1988. https://doi.org/10.1109/59.192932
0021 %   and
0022 %   G. X. Luo and A. Semlyen, "Efficient load flow for large weakly meshed
0023 %   networks," IEEE Transactions on Power Systems, vol. 5, no. 4,
0024 %   pp. 1309-1316, Nov 1990. https://doi.org/10.1109/59.99382
0025 %
0026 %   See also RADIAL_PF.
0027 
0028 %% initialize
0029 tol      = mpopt.pf.tol;
0030 iter_max = mpopt.pf.radial.max_it;
0031 vcorr    = mpopt.pf.radial.vcorr == 1;
0032 Sd(pv) = Sd(pv) - Pg;
0033 V = Vslack * ones(nb,1);
0034 Vold = V;
0035 iter = 0;
0036 success = 0;
0037 % ZIP load model
0038 pw = mpopt.exp.sys_wide_zip_loads.pw;
0039 qw = mpopt.exp.sys_wide_zip_loads.qw;
0040 if isempty(pw)
0041     pw = [1 0 0];
0042 end
0043 if isempty(qw)
0044     qw = pw;
0045 end
0046 Sdz = real(Sd) * pw(3) + 1j * imag(Sd) * qw(3); % constant impedance
0047 Sdi = real(Sd) * pw(2) + 1j * imag(Sd) * qw(2); % constant current
0048 Sdp = real(Sd) * pw(1) + 1j * imag(Sd) * qw(1); % constant power
0049 % Add artificial branch at the top of the branch list, so that the branch
0050 % index for other branches is equal to the index of their receiving node.
0051  f = [0; f];
0052 Zb = [0; Zb];
0053 nl = nl + 1;
0054 %% make Zpv matrix, for calculation of the PV generators reactive powers
0055 if ~isempty(pv)
0056     Zpv = make_zpv(pv,nb,nl,f,Zb,Yd);
0057     Bpv = (imag(Zpv))^-1;
0058 end
0059 npv = length(pv);
0060 Qpv = zeros(npv,1);
0061 %% do backward-forward iterations
0062 if mpopt.verbose > 1
0063     fprintf('\n it    max V mismatch (p.u.)');
0064     fprintf('\n----  ----------------------');
0065 end
0066 while success == 0 && iter < iter_max
0067     iter = iter + 1;
0068     % calculate load demand using actual voltages
0069     Vm = abs(V);
0070     S = Sdp + Sdi.*Vm + Sdz.*Vm.^2 + conj(Yd).*Vm.^2;
0071     % backward sweep
0072     I = conj(S./V);
0073     for k = nl:-1:2
0074         i = f(k);
0075         I(i) = I(i) + I(k);
0076     end
0077     % forward sweep
0078     for k = 2:nl
0079         i = f(k);
0080         V(k) = V(i) - Zb(k) * I(k);
0081     end
0082     % check for convergence
0083     DU = abs(V - Vold);
0084     DU(isnan(DU)) = inf;
0085     if mpopt.verbose > 1
0086         fprintf('\n%3d        %10.3e', iter, max(DU));
0087     end
0088     if max(DU) > tol
0089         Vold = V;
0090         % update PV generators reactive powers
0091         if ~isempty(pv)
0092             DE = (Vg./abs(V(pv))-1).*real(V(pv)); % Rajicic (VCPF)
0093             DD = Bpv * DE;
0094             if vcorr
0095                 DC = DD .* imag(V(pv))./real(V(pv));
0096                 V_corr = make_vcorr(DC+1j*DD,pv,nb,nl,f,Zb);
0097                 V = V + V_corr;
0098             end
0099             DQ = DD .* abs(V(pv)).^2 ./ real(V(pv));
0100             Qpv = Qpv + DQ;
0101             Sdp(pv) = Sdp(pv) - 1j*DQ;
0102         end
0103     else
0104         success = 1;
0105     end
0106 end
0107 if mpopt.verbose
0108     if success
0109         fprintf('\nCurrent summation converged in %d iterations.\n', iter);
0110     else
0111         fprintf('\nCurrent summation did not converge in %d iterations.\n', iter);
0112     end
0113 end
0114 %% calculate branch flows
0115 % take out the first artificial branch
0116 Sslack = V(1)*conj(I(1));
0117 I = I(2:end);
0118 f = f(2:end);
0119 % calculate branch flows
0120 Sf = V(f).*conj(I) + conj(Ybf) .* abs(V(f)).^2;
0121 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