DSBR_DV Computes partial derivatives of power flows w.r.t. voltage. [DSF_DVA, DSF_DVM, DST_DVA, DST_DVM, SF, ST] = DSBR_DV(BRANCH, YF, YT, V) returns four matrices containing partial derivatives of the complex branch power flows at "from" and "to" ends of each branch w.r.t voltage magnitude and voltage angle respectively (for all buses). If YF is a sparse matrix, the partial derivative matrices will be as well. Optionally returns vectors containing the power flows themselves. The following explains the expressions used to form the matrices: If = Yf * V; Sf = diag(Vf) * conj(If) = diag(conj(If)) * Vf Partials of V, Vf & If w.r.t. voltage angles dV/dVa = j * diag(V) dVf/dVa = sparse(1:nl, f, j * V(f)) = j * sparse(1:nl, f, V(f)) dIf/dVa = Yf * dV/dVa = Yf * j * diag(V) Partials of V, Vf & If w.r.t. voltage magnitudes dV/dVm = diag(V./abs(V)) dVf/dVm = sparse(1:nl, f, V(f)./abs(V(f)) dIf/dVm = Yf * dV/dVm = Yf * diag(V./abs(V)) Partials of Sf w.r.t. voltage angles dSf/dVa = diag(Vf) * conj(dIf/dVa) + diag(conj(If)) * dVf/dVa = diag(Vf) * conj(Yf * j * diag(V)) + conj(diag(If)) * j * sparse(1:nl, f, V(f)) = -j * diag(Vf) * conj(Yf * diag(V)) + j * conj(diag(If)) * sparse(1:nl, f, V(f)) = j * (conj(diag(If)) * sparse(1:nl, f, V(f)) - diag(Vf) * conj(Yf * diag(V))) Partials of Sf w.r.t. voltage magnitudes dSf/dVm = diag(Vf) * conj(dIf/dVm) + diag(conj(If)) * dVf/dVm = diag(Vf) * conj(Yf * diag(V./abs(V))) + conj(diag(If)) * sparse(1:nl, f, V(f)./abs(V(f))) Derivations for "to" bus are similar. Example: [Ybus, Yf, Yt] = makeYbus(baseMVA, bus, branch); [dSf_dVa, dSf_dVm, dSt_dVa, dSt_dVm, Sf, St] = ... dSbr_dV(branch, Yf, Yt, V); For more details on the derivations behind the derivative code used in MATPOWER information, see: [TN2] R. D. Zimmerman, "AC Power Flows, Generalized OPF Costs and their Derivatives using Complex Matrix Notation", MATPOWER Technical Note 2, February 2010. http://www.pserc.cornell.edu/matpower/TN2-OPF-Derivatives.pdf
0001 function [dSf_dVa, dSf_dVm, dSt_dVa, dSt_dVm, Sf, St] = dSbr_dV(branch, Yf, Yt, V) 0002 %DSBR_DV Computes partial derivatives of power flows w.r.t. voltage. 0003 % [DSF_DVA, DSF_DVM, DST_DVA, DST_DVM, SF, ST] = DSBR_DV(BRANCH, YF, YT, V) 0004 % returns four matrices containing partial derivatives of the complex 0005 % branch power flows at "from" and "to" ends of each branch w.r.t voltage 0006 % magnitude and voltage angle respectively (for all buses). If YF is a 0007 % sparse matrix, the partial derivative matrices will be as well. Optionally 0008 % returns vectors containing the power flows themselves. The following 0009 % explains the expressions used to form the matrices: 0010 % 0011 % If = Yf * V; 0012 % Sf = diag(Vf) * conj(If) = diag(conj(If)) * Vf 0013 % 0014 % Partials of V, Vf & If w.r.t. voltage angles 0015 % dV/dVa = j * diag(V) 0016 % dVf/dVa = sparse(1:nl, f, j * V(f)) = j * sparse(1:nl, f, V(f)) 0017 % dIf/dVa = Yf * dV/dVa = Yf * j * diag(V) 0018 % 0019 % Partials of V, Vf & If w.r.t. voltage magnitudes 0020 % dV/dVm = diag(V./abs(V)) 0021 % dVf/dVm = sparse(1:nl, f, V(f)./abs(V(f)) 0022 % dIf/dVm = Yf * dV/dVm = Yf * diag(V./abs(V)) 0023 % 0024 % Partials of Sf w.r.t. voltage angles 0025 % dSf/dVa = diag(Vf) * conj(dIf/dVa) 0026 % + diag(conj(If)) * dVf/dVa 0027 % = diag(Vf) * conj(Yf * j * diag(V)) 0028 % + conj(diag(If)) * j * sparse(1:nl, f, V(f)) 0029 % = -j * diag(Vf) * conj(Yf * diag(V)) 0030 % + j * conj(diag(If)) * sparse(1:nl, f, V(f)) 0031 % = j * (conj(diag(If)) * sparse(1:nl, f, V(f)) 0032 % - diag(Vf) * conj(Yf * diag(V))) 0033 % 0034 % Partials of Sf w.r.t. voltage magnitudes 0035 % dSf/dVm = diag(Vf) * conj(dIf/dVm) 0036 % + diag(conj(If)) * dVf/dVm 0037 % = diag(Vf) * conj(Yf * diag(V./abs(V))) 0038 % + conj(diag(If)) * sparse(1:nl, f, V(f)./abs(V(f))) 0039 % 0040 % Derivations for "to" bus are similar. 0041 % 0042 % Example: 0043 % [Ybus, Yf, Yt] = makeYbus(baseMVA, bus, branch); 0044 % [dSf_dVa, dSf_dVm, dSt_dVa, dSt_dVm, Sf, St] = ... 0045 % dSbr_dV(branch, Yf, Yt, V); 0046 % 0047 % For more details on the derivations behind the derivative code used 0048 % in MATPOWER information, see: 0049 % 0050 % [TN2] R. D. Zimmerman, "AC Power Flows, Generalized OPF Costs and 0051 % their Derivatives using Complex Matrix Notation", MATPOWER 0052 % Technical Note 2, February 2010. 0053 % http://www.pserc.cornell.edu/matpower/TN2-OPF-Derivatives.pdf 0054 0055 % MATPOWER 0056 % Copyright (c) 1996-2016, Power Systems Engineering Research Center (PSERC) 0057 % by Ray Zimmerman, PSERC Cornell 0058 % 0059 % This file is part of MATPOWER. 0060 % Covered by the 3-clause BSD License (see LICENSE file for details). 0061 % See http://www.pserc.cornell.edu/matpower/ for more info. 0062 0063 %% define named indices into bus, gen, branch matrices 0064 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ... 0065 TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ... 0066 ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch; 0067 0068 %% define 0069 f = branch(:, F_BUS); %% list of "from" buses 0070 t = branch(:, T_BUS); %% list of "to" buses 0071 nl = length(f); 0072 nb = length(V); 0073 0074 %% compute currents 0075 If = Yf * V; 0076 It = Yt * V; 0077 0078 Vnorm = V ./ abs(V); 0079 if issparse(Yf) %% sparse version (if Yf is sparse) 0080 diagVf = sparse(1:nl, 1:nl, V(f), nl, nl); 0081 diagIf = sparse(1:nl, 1:nl, If, nl, nl); 0082 diagVt = sparse(1:nl, 1:nl, V(t), nl, nl); 0083 diagIt = sparse(1:nl, 1:nl, It, nl, nl); 0084 diagV = sparse(1:nb, 1:nb, V, nb, nb); 0085 diagVnorm = sparse(1:nb, 1:nb, Vnorm, nb, nb); 0086 0087 dSf_dVa = 1j * (conj(diagIf) * sparse(1:nl, f, V(f), nl, nb) - diagVf * conj(Yf * diagV)); 0088 dSf_dVm = diagVf * conj(Yf * diagVnorm) + conj(diagIf) * sparse(1:nl, f, Vnorm(f), nl, nb); 0089 dSt_dVa = 1j * (conj(diagIt) * sparse(1:nl, t, V(t), nl, nb) - diagVt * conj(Yt * diagV)); 0090 dSt_dVm = diagVt * conj(Yt * diagVnorm) + conj(diagIt) * sparse(1:nl, t, Vnorm(t), nl, nb); 0091 else %% dense version 0092 diagVf = diag(V(f)); 0093 diagIf = diag(If); 0094 diagVt = diag(V(t)); 0095 diagIt = diag(It); 0096 diagV = diag(V); 0097 diagVnorm = diag(Vnorm); 0098 temp1 = zeros(nl, nb); temp1(sub2ind([nl,nb], (1:nl)', f)) = V(f); 0099 temp2 = zeros(nl, nb); temp2(sub2ind([nl,nb], (1:nl)', f)) = Vnorm(f); 0100 temp3 = zeros(nl, nb); temp3(sub2ind([nl,nb], (1:nl)', t)) = V(t); 0101 temp4 = zeros(nl, nb); temp4(sub2ind([nl,nb], (1:nl)', t)) = Vnorm(t); 0102 0103 dSf_dVa = 1j * (conj(diagIf) * temp1 - diagVf * conj(Yf * diagV)); 0104 dSf_dVm = diagVf * conj(Yf * diagVnorm) + conj(diagIf) * temp2; 0105 dSt_dVa = 1j * (conj(diagIt) * temp3 - diagVt * conj(Yt * diagV)); 0106 dSt_dVm = diagVt * conj(Yt * diagVnorm) + conj(diagIt) * temp4; 0107 end 0108 0109 if nargout > 4 0110 Sf = V(f) .* conj(If); 0111 St = V(t) .* conj(It); 0112 end