Home > matpower6.0 > newtonpf.m

newtonpf

PURPOSE ^

NEWTONPF Solves the power flow using a full Newton's method.

SYNOPSIS ^

function [V, converged, i] = newtonpf(Ybus, Sbus, V0, ref, pv, pq, mpopt)

DESCRIPTION ^

NEWTONPF  Solves the power flow using a full Newton's method.
   [V, CONVERGED, I] = NEWTONPF(YBUS, SBUS, V0, 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, 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 struct 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.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [V, converged, i] = newtonpf(Ybus, Sbus, V0, ref, pv, pq, mpopt)
0002 %NEWTONPF  Solves the power flow using a full Newton's method.
0003 %   [V, CONVERGED, I] = NEWTONPF(YBUS, SBUS, V0, 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, and column vectors with
0007 %   the lists of bus indices for the swing bus, PV buses, and PQ buses,
0008 %   respectively. The bus voltage vector contains the set point for
0009 %   generator (including ref bus) buses, and the reference angle of the
0010 %   swing bus, as well as an initial guess for remaining magnitudes and
0011 %   angles. MPOPT is a MATPOWER options struct which can be used to
0012 %   set the termination tolerance, maximum number of iterations, and
0013 %   output options (see MPOPTION for details). Uses default options if
0014 %   this parameter is not given. Returns the final complex voltages, a
0015 %   flag which indicates whether it converged or not, and the number of
0016 %   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.nr.max_it;
0036 
0037 %% initialize
0038 converged = 0;
0039 i = 0;
0040 V = V0;
0041 Va = angle(V);
0042 Vm = abs(V);
0043 
0044 %% set up indexing for updating V
0045 npv = length(pv);
0046 npq = length(pq);
0047 j1 = 1;         j2 = npv;           %% j1:j2 - V angle of pv buses
0048 j3 = j2 + 1;    j4 = j2 + npq;      %% j3:j4 - V angle of pq buses
0049 j5 = j4 + 1;    j6 = j4 + npq;      %% j5:j6 - V mag of pq buses
0050 
0051 %% evaluate F(x0)
0052 mis = V .* conj(Ybus * V) - Sbus(Vm);
0053 F = [   real(mis([pv; pq]));
0054         imag(mis(pq))   ];
0055 
0056 %% check tolerance
0057 normF = norm(F, inf);
0058 if mpopt.verbose > 1
0059     fprintf('\n it    max P & Q mismatch (p.u.)');
0060     fprintf('\n----  ---------------------------');
0061     fprintf('\n%3d        %10.3e', i, normF);
0062 end
0063 if normF < tol
0064     converged = 1;
0065     if mpopt.verbose > 1
0066         fprintf('\nConverged!\n');
0067     end
0068 end
0069 
0070 %% do Newton iterations
0071 while (~converged && i < max_it)
0072     %% update iteration counter
0073     i = i + 1;
0074     
0075     %% evaluate Jacobian
0076     [dSbus_dVm, dSbus_dVa] = dSbus_dV(Ybus, V);
0077     [dummy, neg_dSd_dVm] = Sbus(Vm);
0078     dSbus_dVm = dSbus_dVm - neg_dSd_dVm;
0079     
0080     j11 = real(dSbus_dVa([pv; pq], [pv; pq]));
0081     j12 = real(dSbus_dVm([pv; pq], pq));
0082     j21 = imag(dSbus_dVa(pq, [pv; pq]));
0083     j22 = imag(dSbus_dVm(pq, pq));
0084     
0085     J = [   j11 j12;
0086             j21 j22;    ];
0087 
0088     %% compute update step
0089     dx = -(J \ F);
0090 
0091     %% update voltage
0092     if npv
0093         Va(pv) = Va(pv) + dx(j1:j2);
0094     end
0095     if npq
0096         Va(pq) = Va(pq) + dx(j3:j4);
0097         Vm(pq) = Vm(pq) + dx(j5:j6);
0098     end
0099     V = Vm .* exp(1j * Va);
0100     Vm = abs(V);            %% update Vm and Va again in case
0101     Va = angle(V);          %% we wrapped around with a negative Vm
0102 
0103     %% evalute F(x)
0104     mis = V .* conj(Ybus * V) - Sbus(Vm);
0105     F = [   real(mis(pv));
0106             real(mis(pq));
0107             imag(mis(pq))   ];
0108 
0109     %% check for convergence
0110     normF = norm(F, inf);
0111     if mpopt.verbose > 1
0112         fprintf('\n%3d        %10.3e', i, normF);
0113     end
0114     if normF < tol
0115         converged = 1;
0116         if mpopt.verbose
0117             fprintf('\nNewton''s method power flow converged in %d iterations.\n', i);
0118         end
0119     end
0120 end
0121 
0122 if mpopt.verbose
0123     if ~converged
0124         fprintf('\nNewton''s method power flow did not converge in %d iterations.\n', i);
0125     end
0126 end

Generated on Fri 16-Dec-2016 12:45:37 by m2html © 2005