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 % $Id: fdpf.m 2453 2014-12-05 21:05:50Z 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.fd.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 0069 %% evaluate initial mismatch 0070 mis = (V .* conj(Ybus * V) - Sbus) ./ Vm; 0071 P = real(mis([pv; pq])); 0072 Q = imag(mis(pq)); 0073 0074 %% check tolerance 0075 normP = norm(P, inf); 0076 normQ = norm(Q, inf); 0077 if mpopt.verbose > 1 0078 fprintf('\niteration max mismatch (p.u.) '); 0079 fprintf('\ntype # P Q '); 0080 fprintf('\n---- ---- ----------- -----------'); 0081 fprintf('\n - %3d %10.3e %10.3e', i, normP, normQ); 0082 end 0083 if normP < tol && normQ < tol 0084 converged = 1; 0085 if mpopt.verbose > 1 0086 fprintf('\nConverged!\n'); 0087 end 0088 end 0089 0090 %% reduce B matrices 0091 Bp = Bp([pv; pq], [pv; pq]); 0092 Bpp = Bpp(pq, pq); 0093 0094 %% factor B matrices 0095 [Lp, Up, pp, qp ] = lu(Bp, 'vector'); 0096 [Lpp, Upp, ppp, qpp] = lu(Bpp, 'vector'); 0097 [junk, iqp ] = sort(qp); 0098 [junk, iqpp] = sort(qpp); 0099 % [~, iqp ] = sort(qp); 0100 % [~, iqpp] = sort(qpp); 0101 0102 %% do P and Q iterations 0103 while (~converged && i < max_it) 0104 %% update iteration counter 0105 i = i + 1; 0106 0107 %%----- do P iteration, update Va ----- 0108 dVa = -( Up \ (Lp \ P(pp)) ); 0109 dVa = dVa(iqp); 0110 0111 %% update voltage 0112 Va([pv; pq]) = Va([pv; pq]) + dVa; 0113 V = Vm .* exp(1j * Va); 0114 0115 %% evalute mismatch 0116 mis = (V .* conj(Ybus * V) - Sbus) ./ Vm; 0117 P = real(mis([pv; pq])); 0118 Q = imag(mis(pq)); 0119 0120 %% check tolerance 0121 normP = norm(P, inf); 0122 normQ = norm(Q, inf); 0123 if mpopt.verbose > 1 0124 fprintf('\n P %3d %10.3e %10.3e', i, normP, normQ); 0125 end 0126 if normP < tol && normQ < tol 0127 converged = 1; 0128 if mpopt.verbose 0129 fprintf('\nFast-decoupled power flow converged in %d P-iterations and %d Q-iterations.\n', i, i-1); 0130 end 0131 break; 0132 end 0133 0134 %%----- do Q iteration, update Vm ----- 0135 dVm = -( Upp \ (Lpp \ Q(ppp)) ); 0136 dVm = dVm(iqpp); 0137 0138 %% update voltage 0139 Vm(pq) = Vm(pq) + dVm; 0140 V = Vm .* exp(1j * Va); 0141 0142 %% evalute mismatch 0143 mis = (V .* conj(Ybus * V) - Sbus) ./ Vm; 0144 P = real(mis([pv; pq])); 0145 Q = imag(mis(pq)); 0146 0147 %% check tolerance 0148 normP = norm(P, inf); 0149 normQ = norm(Q, inf); 0150 if mpopt.verbose > 1 0151 fprintf('\n Q %3d %10.3e %10.3e', i, normP, normQ); 0152 end 0153 if normP < tol && normQ < tol 0154 converged = 1; 0155 if mpopt.verbose 0156 fprintf('\nFast-decoupled power flow converged in %d P-iterations and %d Q-iterations.\n', i, i); 0157 end 0158 break; 0159 end 0160 end 0161 0162 if mpopt.verbose 0163 if ~converged 0164 fprintf('\nFast-decoupled power flow did not converge in %d iterations.\n', i); 0165 end 0166 end