CPF_CORRECTOR Solves the corrector step of a continuation power flow using a
full Newton method with selected parameterization scheme.
[V, CONVERGED, I, LAM] = CPF_CORRECTOR(YBUS, SBUS, V0, REF, PV, PQ, ...
LAM0, SXFR, VPRV, LPRV, Z, STEP, PARAMETERIZATION, MPOPT)
solves for bus voltages and lambda 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, the number
of iterations performed, and the final lambda.
The extra continuation inputs are LAM0 (initial predicted lambda),
SXFR ([delP+j*delQ] transfer/loading vector for all buses), VPRV
(final complex V corrector solution from previous continuation step),
LAMPRV (final lambda corrector solution from previous continuation step),
Z (normalized predictor for all buses), and STEP (continuation step size).
The extra continuation output is LAM (final corrector lambda).
See also RUNCPF.0001 function [V, converged, i, lam] = cpf_corrector(Ybus, Sbus, V0, ref, pv, pq, ... 0002 lam0, Sxfr, Vprv, lamprv, z, step, parameterization, mpopt) 0003 %CPF_CORRECTOR Solves the corrector step of a continuation power flow using a 0004 % full Newton method with selected parameterization scheme. 0005 % [V, CONVERGED, I, LAM] = CPF_CORRECTOR(YBUS, SBUS, V0, REF, PV, PQ, ... 0006 % LAM0, SXFR, VPRV, LPRV, Z, STEP, PARAMETERIZATION, MPOPT) 0007 % solves for bus voltages and lambda given the full system admittance 0008 % matrix (for all buses), the complex bus power injection vector (for 0009 % all buses), the initial vector of complex bus voltages, and column 0010 % vectors with the lists of bus indices for the swing bus, PV buses, and 0011 % PQ buses, respectively. The bus voltage vector contains the set point 0012 % for generator (including ref bus) buses, and the reference angle of the 0013 % swing bus, as well as an initial guess for remaining magnitudes and 0014 % angles. MPOPT is a MATPOWER options struct which can be used to 0015 % set the termination tolerance, maximum number of iterations, and 0016 % output options (see MPOPTION for details). Uses default options if 0017 % this parameter is not given. Returns the final complex voltages, a 0018 % flag which indicates whether it converged or not, the number 0019 % of iterations performed, and the final lambda. 0020 % 0021 % The extra continuation inputs are LAM0 (initial predicted lambda), 0022 % SXFR ([delP+j*delQ] transfer/loading vector for all buses), VPRV 0023 % (final complex V corrector solution from previous continuation step), 0024 % LAMPRV (final lambda corrector solution from previous continuation step), 0025 % Z (normalized predictor for all buses), and STEP (continuation step size). 0026 % The extra continuation output is LAM (final corrector lambda). 0027 % 0028 % See also RUNCPF. 0029 0030 % MATPOWER 0031 % $Id: cpf_corrector.m 2238 2013-12-12 15:32:29Z ray $ 0032 % by Ray Zimmerman, PSERC Cornell, 0033 % Shrirang Abhyankar, Argonne National Laboratory, 0034 % and Alexander Flueck, IIT 0035 % Copyright (c) 1996-2013 by Power System Engineering Research Center (PSERC) 0036 % 0037 % Modified by Alexander J. Flueck, Illinois Institute of Technology 0038 % 2001.02.22 - corrector.m (ver 1.0) based on newtonpf.m (MATPOWER 2.0) 0039 % 0040 % Modified by Shrirang Abhyankar, Argonne National Laboratory 0041 % Updated to be compatible with MATPOWER version 4.1) 0042 % 0043 % This file is part of MATPOWER. 0044 % See http://www.pserc.cornell.edu/matpower/ for more info. 0045 % 0046 % MATPOWER is free software: you can redistribute it and/or modify 0047 % it under the terms of the GNU General Public License as published 0048 % by the Free Software Foundation, either version 3 of the License, 0049 % or (at your option) any later version. 0050 % 0051 % MATPOWER is distributed in the hope that it will be useful, 0052 % but WITHOUT ANY WARRANTY; without even the implied warranty of 0053 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 0054 % GNU General Public License for more details. 0055 % 0056 % You should have received a copy of the GNU General Public License 0057 % along with MATPOWER. If not, see <http://www.gnu.org/licenses/>. 0058 % 0059 % Additional permission under GNU GPL version 3 section 7 0060 % 0061 % If you modify MATPOWER, or any covered work, to interface with 0062 % other modules (such as MATLAB code and MEX-files) available in a 0063 % MATLAB(R) or comparable environment containing parts covered 0064 % under other licensing terms, the licensors of MATPOWER grant 0065 % you additional permission to convey the resulting work. 0066 0067 %% default arguments 0068 if nargin < 14 0069 mpopt = mpoption; 0070 end 0071 0072 %% options 0073 tol = mpopt.pf.tol; 0074 max_it = mpopt.pf.nr.max_it; 0075 0076 %% initialize 0077 converged = 0; 0078 i = 0; 0079 V = V0; 0080 Va = angle(V); 0081 Vm = abs(V); 0082 lam = lam0; %% set lam to initial lam0 0083 0084 %% set up indexing for updating V 0085 npv = length(pv); 0086 npq = length(pq); 0087 nb = length(V); %% number of buses 0088 j1 = 1; j2 = npv; %% j1:j2 - V angle of pv buses 0089 j3 = j2 + 1; j4 = j2 + npq; %% j3:j4 - V angle of pq buses 0090 j5 = j4 + 1; j6 = j4 + npq; %% j5:j6 - V mag of pq buses 0091 j7 = j6 + 1; j8 = j6 + 1; %% j7:j8 - lambda 0092 0093 %% evaluate F(x0, lam0), including Sxfr transfer/loading 0094 mis = V .* conj(Ybus * V) - Sbus - lam*Sxfr; 0095 F = [ real(mis([pv; pq])); 0096 imag(mis(pq)) ]; 0097 0098 %% evaluate P(x0, lambda0) 0099 P = cpf_p(parameterization, step, z, V, lam, Vprv, lamprv, pv, pq); 0100 0101 %% augment F(x,lambda) with P(x,lambda) 0102 F = [ F; 0103 P ]; 0104 0105 %% check tolerance 0106 normF = norm(F, inf); 0107 if mpopt.verbose > 1 0108 fprintf('\n it max P & Q mismatch (p.u.)'); 0109 fprintf('\n---- ---------------------------'); 0110 fprintf('\n%3d %10.3e', i, normF); 0111 end 0112 if normF < tol 0113 converged = 1; 0114 if mpopt.verbose > 1 0115 fprintf('\nConverged!\n'); 0116 end 0117 end 0118 0119 %% do Newton iterations 0120 while (~converged && i < max_it) 0121 %% update iteration counter 0122 i = i + 1; 0123 0124 %% evaluate Jacobian 0125 [dSbus_dVm, dSbus_dVa] = dSbus_dV(Ybus, V); 0126 0127 j11 = real(dSbus_dVa([pv; pq], [pv; pq])); 0128 j12 = real(dSbus_dVm([pv; pq], pq)); 0129 j21 = imag(dSbus_dVa(pq, [pv; pq])); 0130 j22 = imag(dSbus_dVm(pq, pq)); 0131 0132 J = [ j11 j12; 0133 j21 j22; ]; 0134 0135 dF_dlam = -[real(Sxfr([pv; pq])); imag(Sxfr(pq))]; 0136 [dP_dV, dP_dlam] = cpf_p_jac(parameterization, z, V, lam, Vprv, lamprv, pv, pq); 0137 0138 %% augment J with real/imag -Sxfr and z^T 0139 J = [ J dF_dlam; 0140 dP_dV dP_dlam ]; 0141 0142 %% compute update step 0143 dx = -(J \ F); 0144 0145 %% update voltage 0146 if npv 0147 Va(pv) = Va(pv) + dx(j1:j2); 0148 end 0149 if npq 0150 Va(pq) = Va(pq) + dx(j3:j4); 0151 Vm(pq) = Vm(pq) + dx(j5:j6); 0152 end 0153 V = Vm .* exp(1j * Va); 0154 Vm = abs(V); %% update Vm and Va again in case 0155 Va = angle(V); %% we wrapped around with a negative Vm 0156 0157 %% update lambda 0158 lam = lam + dx(j7:j8); 0159 0160 %% evalute F(x, lam) 0161 mis = V .* conj(Ybus * V) - Sbus - lam*Sxfr; 0162 F = [ real(mis(pv)); 0163 real(mis(pq)); 0164 imag(mis(pq)) ]; 0165 0166 %% evaluate P(x, lambda) 0167 P = cpf_p(parameterization, step, z, V, lam, Vprv, lamprv, pv, pq); 0168 0169 %% augment F(x,lambda) with P(x,lambda) 0170 F = [ F; 0171 P ]; 0172 0173 %% check for convergence 0174 normF = norm(F, inf); 0175 if mpopt.verbose > 1 0176 fprintf('\n%3d %10.3e', i, normF); 0177 end 0178 if normF < tol 0179 converged = 1; 0180 if mpopt.verbose 0181 fprintf('\nNewton''s method corrector converged in %d iterations.\n', i); 0182 end 0183 end 0184 end 0185 0186 if mpopt.verbose 0187 if ~converged 0188 fprintf('\nNewton''s method corrector did not converge in %d iterations.\n', i); 0189 end 0190 end