fdpf

FDPF  Solves the power flow using a fast decoupled method.
   [V, CONVERGED, I] = FDPF(YBUS, SBUS, V0, BP, BPP, REF, PV, PQ, MPOPT)
   solves for bus voltages given the full system admittance matrix (for
   all buses), the complex bus power injection vector (for all buses),
   the initial vector of complex bus voltages, the FDPF matrices B prime
   and B double prime, and column vectors with the lists of bus indices
   for the swing bus, PV buses, and PQ buses, respectively. The bus voltage
   vector contains the set point for generator (including ref bus)
   buses, and the reference angle of the swing bus, as well as an initial
   guess for remaining magnitudes and angles. MPOPT is a MATPOWER options
   vector which can be used to set the termination tolerance, maximum
   number of iterations, and output options (see MPOPTION for details).
   Uses default options if this parameter is not given. Returns the
   final complex voltages, a flag which indicates whether it converged
   or not, and the number of iterations performed.

   See also RUNPF.

0001 function [V, converged, i] = fdpf(Ybus, Sbus, V0, Bp, Bpp, ref, pv, pq, mpopt)
0002 %FDPF  Solves the power flow using a fast decoupled method.
0003 %   [V, CONVERGED, I] = FDPF(YBUS, SBUS, V0, BP, BPP, REF, PV, PQ, MPOPT)
0004 %   solves for bus voltages given the full system admittance matrix (for
0005 %   all buses), the complex bus power injection vector (for all buses),
0006 %   the initial vector of complex bus voltages, the FDPF matrices B prime
0007 %   and B double prime, and column vectors with the lists of bus indices
0008 %   for the swing bus, PV buses, and PQ buses, respectively. The bus voltage
0009 %   vector contains the set point for generator (including ref bus)
0010 %   buses, and the reference angle of the swing bus, as well as an initial
0011 %   guess for remaining magnitudes and angles. MPOPT is a MATPOWER options
0012 %   vector which can be used to set the termination tolerance, maximum
0013 %   number of iterations, and output options (see MPOPTION for details).
0014 %   Uses default options if this parameter is not given. Returns the
0015 %   final complex voltages, a flag which indicates whether it converged
0016 %   or not, and the number of iterations performed.
0017 %
0018 %   See also RUNPF.
0019 
0020 %   MATPOWER
0021 %   Copyright (c) 1996-2016, Power Systems Engineering Research Center (PSERC)
0022 %   by Ray Zimmerman, PSERC Cornell
0023 %
0024 %   This file is part of MATPOWER.
0025 %   Covered by the 3-clause BSD License (see LICENSE file for details).
0026 %   See http://www.pserc.cornell.edu/matpower/ for more info.
0027 
0028 %% default arguments
0029 if nargin < 7
0030     mpopt = mpoption;
0031 end
0032 
0033 %% options
0034 tol     = mpopt.pf.tol;
0035 max_it  = mpopt.pf.fd.max_it;
0036 if have_fcn('matlab') && have_fcn('matlab', 'vnum') < 7.3
0037     lu_vec = 0;     %% lu(..., 'vector') syntax not supported
0038 else
0039     lu_vec = 1;
0040 end
0041 
0042 %% initialize
0043 converged = 0;
0044 i = 0;
0045 V = V0;
0046 Va = angle(V);
0047 Vm = abs(V);
0048 
0049 %% set up indexing for updating V
0050 npv = length(pv);
0051 npq = length(pq);
0052 
0053 %% evaluate initial mismatch
0054 mis = (V .* conj(Ybus * V) - Sbus(Vm)) ./ Vm;
0055 P = real(mis([pv; pq]));
0056 Q = imag(mis(pq));
0057 
0058 %% check tolerance
0059 normP = norm(P, inf);
0060 normQ = norm(Q, inf);
0061 if mpopt.verbose > 1
0062     fprintf('\niteration     max mismatch (p.u.)  ');
0063     fprintf('\ntype   #        P            Q     ');
0064     fprintf('\n---- ----  -----------  -----------');
0065     fprintf('\n  -  %3d   %10.3e   %10.3e', i, normP, normQ);
0066 end
0067 if normP < tol && normQ < tol
0068     converged = 1;
0069     if mpopt.verbose > 1
0070         fprintf('\nConverged!\n');
0071     end
0072 end
0073 
0074 %% reduce B matrices
0075 Bp = Bp([pv; pq], [pv; pq]);
0076 Bpp = Bpp(pq, pq);
0077 
0078 %% factor B matrices
0079 if lu_vec
0080     [Lp,  Up,  pp,  qp ] = lu(Bp,  'vector');
0081     [Lpp, Upp, ppp, qpp] = lu(Bpp, 'vector');
0082     [junk, iqp ] = sort(qp);
0083     [junk, iqpp] = sort(qpp);
0084     % [~, iqp ] = sort(qp);
0085     % [~, iqpp] = sort(qpp);
0086 else
0087     [Lp, Up, Pp] = lu(Bp);
0088     [Lpp, Upp, Ppp] = lu(Bpp);
0089 end
0090 
0091 %% do P and Q iterations
0092 while (~converged && i < max_it)
0093     %% update iteration counter
0094     i = i + 1;
0095 
0096     %%-----  do P iteration, update Va  -----
0097     if lu_vec
0098         dVa = -( Up \  (Lp \ P(pp)) );
0099         dVa = dVa(iqp);
0100     else
0101         dVa = -( Up \  (Lp \ (Pp * P)));
0102     end
0103 
0104     %% update voltage
0105     Va([pv; pq]) = Va([pv; pq]) + dVa;
0106     V = Vm .* exp(1j * Va);
0107 
0108     %% evalute mismatch
0109     mis = (V .* conj(Ybus * V) - Sbus(Vm)) ./ Vm;
0110     P = real(mis([pv; pq]));
0111     Q = imag(mis(pq));
0112     
0113     %% check tolerance
0114     normP = norm(P, inf);
0115     normQ = norm(Q, inf);
0116     if mpopt.verbose > 1
0117         fprintf('\n  P  %3d   %10.3e   %10.3e', i, normP, normQ);
0118     end
0119     if normP < tol && normQ < tol
0120         converged = 1;
0121         if mpopt.verbose
0122             fprintf('\nFast-decoupled power flow converged in %d P-iterations and %d Q-iterations.\n', i, i-1);
0123         end
0124         break;
0125     end
0126 
0127     %%-----  do Q iteration, update Vm  -----
0128     if lu_vec
0129         dVm = -( Upp \ (Lpp \ Q(ppp)) );
0130         dVm = dVm(iqpp);
0131     else
0132         dVm = -( Upp \ (Lpp \ (Ppp * Q)) );
0133     end
0134 
0135     %% update voltage
0136     Vm(pq) = Vm(pq) + dVm;
0137     V = Vm .* exp(1j * Va);
0138 
0139     %% evalute mismatch
0140     mis = (V .* conj(Ybus * V) - Sbus(Vm)) ./ Vm;
0141     P = real(mis([pv; pq]));
0142     Q = imag(mis(pq));
0143     
0144     %% check tolerance
0145     normP = norm(P, inf);
0146     normQ = norm(Q, inf);
0147     if mpopt.verbose > 1
0148         fprintf('\n  Q  %3d   %10.3e   %10.3e', i, normP, normQ);
0149     end
0150     if normP < tol && normQ < tol
0151         converged = 1;
0152         if mpopt.verbose
0153             fprintf('\nFast-decoupled power flow converged in %d P-iterations and %d Q-iterations.\n', i, i);
0154         end
0155         break;
0156     end
0157 end
0158 
0159 if mpopt.verbose
0160     if ~converged
0161         fprintf('\nFast-decoupled power flow did not converge in %d iterations.\n', i);
0162     end
0163 end