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-2015 by Power System Engineering Research Center (PSERC)
0022 %   by Ray Zimmerman, PSERC Cornell
0023 %
0024 %   $Id: fdpf.m 2644 2015-03-11 19:34:22Z ray $
0025 %
0026 %   This file is part of MATPOWER.
0027 %   Covered by the 3-clause BSD License (see LICENSE file for details).
0028 %   See http://www.pserc.cornell.edu/matpower/ for more info.
0029 
0030 %% default arguments
0031 if nargin < 7
0032     mpopt = mpoption;
0033 end
0034 
0035 %% options
0036 tol     = mpopt.pf.tol;
0037 max_it  = mpopt.pf.fd.max_it;
0038 if have_fcn('matlab') && have_fcn('matlab', 'vnum') < 7.3
0039     lu_vec = 0;     %% lu(..., 'vector') syntax not supported
0040 else
0041     lu_vec = 1;
0042 end
0043 
0044 %% initialize
0045 converged = 0;
0046 i = 0;
0047 V = V0;
0048 Va = angle(V);
0049 Vm = abs(V);
0050 
0051 %% set up indexing for updating V
0052 npv = length(pv);
0053 npq = length(pq);
0054 
0055 %% evaluate initial mismatch
0056 mis = (V .* conj(Ybus * V) - Sbus) ./ Vm;
0057 P = real(mis([pv; pq]));
0058 Q = imag(mis(pq));
0059 
0060 %% check tolerance
0061 normP = norm(P, inf);
0062 normQ = norm(Q, inf);
0063 if mpopt.verbose > 1
0064     fprintf('\niteration     max mismatch (p.u.)  ');
0065     fprintf('\ntype   #        P            Q     ');
0066     fprintf('\n---- ----  -----------  -----------');
0067     fprintf('\n  -  %3d   %10.3e   %10.3e', i, normP, normQ);
0068 end
0069 if normP < tol && normQ < tol
0070     converged = 1;
0071     if mpopt.verbose > 1
0072         fprintf('\nConverged!\n');
0073     end
0074 end
0075 
0076 %% reduce B matrices
0077 Bp = Bp([pv; pq], [pv; pq]);
0078 Bpp = Bpp(pq, pq);
0079 
0080 %% factor B matrices
0081 if lu_vec
0082     [Lp,  Up,  pp,  qp ] = lu(Bp,  'vector');
0083     [Lpp, Upp, ppp, qpp] = lu(Bpp, 'vector');
0084     [junk, iqp ] = sort(qp);
0085     [junk, iqpp] = sort(qpp);
0086     % [~, iqp ] = sort(qp);
0087     % [~, iqpp] = sort(qpp);
0088 else
0089     [Lp, Up, Pp] = lu(Bp);
0090     [Lpp, Upp, Ppp] = lu(Bpp);
0091 end
0092 
0093 %% do P and Q iterations
0094 while (~converged && i < max_it)
0095     %% update iteration counter
0096     i = i + 1;
0097 
0098     %%-----  do P iteration, update Va  -----
0099     if lu_vec
0100         dVa = -( Up \  (Lp \ P(pp)) );
0101         dVa = dVa(iqp);
0102     else
0103         dVa = -( Up \  (Lp \ (Pp * P)));
0104     end
0105 
0106     %% update voltage
0107     Va([pv; pq]) = Va([pv; pq]) + dVa;
0108     V = Vm .* exp(1j * Va);
0109 
0110     %% evalute mismatch
0111     mis = (V .* conj(Ybus * V) - Sbus) ./ Vm;
0112     P = real(mis([pv; pq]));
0113     Q = imag(mis(pq));
0114     
0115     %% check tolerance
0116     normP = norm(P, inf);
0117     normQ = norm(Q, inf);
0118     if mpopt.verbose > 1
0119         fprintf('\n  P  %3d   %10.3e   %10.3e', i, normP, normQ);
0120     end
0121     if normP < tol && normQ < tol
0122         converged = 1;
0123         if mpopt.verbose
0124             fprintf('\nFast-decoupled power flow converged in %d P-iterations and %d Q-iterations.\n', i, i-1);
0125         end
0126         break;
0127     end
0128 
0129     %%-----  do Q iteration, update Vm  -----
0130     if lu_vec
0131         dVm = -( Upp \ (Lpp \ Q(ppp)) );
0132         dVm = dVm(iqpp);
0133     else
0134         dVm = -( Upp \ (Lpp \ (Ppp * Q)) );
0135     end
0136 
0137     %% update voltage
0138     Vm(pq) = Vm(pq) + dVm;
0139     V = Vm .* exp(1j * Va);
0140 
0141     %% evalute mismatch
0142     mis = (V .* conj(Ybus * V) - Sbus) ./ Vm;
0143     P = real(mis([pv; pq]));
0144     Q = imag(mis(pq));
0145     
0146     %% check tolerance
0147     normP = norm(P, inf);
0148     normQ = norm(Q, inf);
0149     if mpopt.verbose > 1
0150         fprintf('\n  Q  %3d   %10.3e   %10.3e', i, normP, normQ);
0151     end
0152     if normP < tol && normQ < tol
0153         converged = 1;
0154         if mpopt.verbose
0155             fprintf('\nFast-decoupled power flow converged in %d P-iterations and %d Q-iterations.\n', i, i);
0156         end
0157         break;
0158     end
0159 end
0160 
0161 if mpopt.verbose
0162     if ~converged
0163         fprintf('\nFast-decoupled power flow did not converge in %d iterations.\n', i);
0164     end
0165 end