Home > matpower5.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 %   $Id: newtonpf.m 2229 2013-12-11 01:28:09Z ray $
0022 %   by Ray Zimmerman, PSERC Cornell
0023 %   Copyright (c) 1996-2011 by Power System Engineering Research Center (PSERC)
0024 %
0025 %   This file is part of MATPOWER.
0026 %   See http://www.pserc.cornell.edu/matpower/ for more info.
0027 %
0028 %   MATPOWER is free software: you can redistribute it and/or modify
0029 %   it under the terms of the GNU General Public License as published
0030 %   by the Free Software Foundation, either version 3 of the License,
0031 %   or (at your option) any later version.
0032 %
0033 %   MATPOWER is distributed in the hope that it will be useful,
0034 %   but WITHOUT ANY WARRANTY; without even the implied warranty of
0035 %   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
0036 %   GNU General Public License for more details.
0037 %
0038 %   You should have received a copy of the GNU General Public License
0039 %   along with MATPOWER. If not, see <http://www.gnu.org/licenses/>.
0040 %
0041 %   Additional permission under GNU GPL version 3 section 7
0042 %
0043 %   If you modify MATPOWER, or any covered work, to interface with
0044 %   other modules (such as MATLAB code and MEX-files) available in a
0045 %   MATLAB(R) or comparable environment containing parts covered
0046 %   under other licensing terms, the licensors of MATPOWER grant
0047 %   you additional permission to convey the resulting work.
0048 
0049 %% default arguments
0050 if nargin < 7
0051     mpopt = mpoption;
0052 end
0053 
0054 %% options
0055 tol     = mpopt.pf.tol;
0056 max_it  = mpopt.pf.nr.max_it;
0057 
0058 %% initialize
0059 converged = 0;
0060 i = 0;
0061 V = V0;
0062 Va = angle(V);
0063 Vm = abs(V);
0064 
0065 %% set up indexing for updating V
0066 npv = length(pv);
0067 npq = length(pq);
0068 j1 = 1;         j2 = npv;           %% j1:j2 - V angle of pv buses
0069 j3 = j2 + 1;    j4 = j2 + npq;      %% j3:j4 - V angle of pq buses
0070 j5 = j4 + 1;    j6 = j4 + npq;      %% j5:j6 - V mag of pq buses
0071 
0072 %% evaluate F(x0)
0073 mis = V .* conj(Ybus * V) - Sbus;
0074 F = [   real(mis([pv; pq]));
0075         imag(mis(pq))   ];
0076 
0077 %% check tolerance
0078 normF = norm(F, inf);
0079 if mpopt.verbose > 1
0080     fprintf('\n it    max P & Q mismatch (p.u.)');
0081     fprintf('\n----  ---------------------------');
0082     fprintf('\n%3d        %10.3e', i, normF);
0083 end
0084 if normF < tol
0085     converged = 1;
0086     if mpopt.verbose > 1
0087         fprintf('\nConverged!\n');
0088     end
0089 end
0090 
0091 %% do Newton iterations
0092 while (~converged && i < max_it)
0093     %% update iteration counter
0094     i = i + 1;
0095     
0096     %% evaluate Jacobian
0097     [dSbus_dVm, dSbus_dVa] = dSbus_dV(Ybus, V);
0098     
0099     j11 = real(dSbus_dVa([pv; pq], [pv; pq]));
0100     j12 = real(dSbus_dVm([pv; pq], pq));
0101     j21 = imag(dSbus_dVa(pq, [pv; pq]));
0102     j22 = imag(dSbus_dVm(pq, pq));
0103     
0104     J = [   j11 j12;
0105             j21 j22;    ];
0106 
0107     %% compute update step
0108     dx = -(J \ F);
0109 
0110     %% update voltage
0111     if npv
0112         Va(pv) = Va(pv) + dx(j1:j2);
0113     end
0114     if npq
0115         Va(pq) = Va(pq) + dx(j3:j4);
0116         Vm(pq) = Vm(pq) + dx(j5:j6);
0117     end
0118     V = Vm .* exp(1j * Va);
0119     Vm = abs(V);            %% update Vm and Va again in case
0120     Va = angle(V);          %% we wrapped around with a negative Vm
0121 
0122     %% evalute F(x)
0123     mis = V .* conj(Ybus * V) - Sbus;
0124     F = [   real(mis(pv));
0125             real(mis(pq));
0126             imag(mis(pq))   ];
0127 
0128     %% check for convergence
0129     normF = norm(F, inf);
0130     if mpopt.verbose > 1
0131         fprintf('\n%3d        %10.3e', i, normF);
0132     end
0133     if normF < tol
0134         converged = 1;
0135         if mpopt.verbose
0136             fprintf('\nNewton''s method power flow converged in %d iterations.\n', i);
0137         end
0138     end
0139 end
0140 
0141 if mpopt.verbose
0142     if ~converged
0143         fprintf('\nNewton''s method power flow did not converge in %d iterations.\n', i);
0144     end
0145 end

Generated on Mon 26-Jan-2015 15:21:31 by m2html © 2005