gycall.m

电力系统分析计算程序

function Gycall(a)

global DAE

if ~a.n, return, end

V1 = a.u.*DAE.y(a.v1);
V2 = a.u.*DAE.y(a.v2);  
V12 = V1.*V2;
y = admittance(a);
g = real(y);
b = imag(y);
m = DAE.x(a.m);
m2 = m.*m;

[s12,c12] = angles(a);

k1 = m.*(c12.*g+s12.*b);
k2 = m.*(c12.*g-s12.*b);
k3 = m.*(s12.*g-c12.*b);
k4 = m.*(s12.*g+c12.*b);

DAE.Gy = DAE.Gy + ...
          sparse(a.bus1,a.v1,2.*V1.*g-V2.*k1,DAE.m,DAE.m) + ...
          sparse(a.bus1,a.v2,-V1.*k1,DAE.m,DAE.m) + ...
          sparse(a.bus2,a.v1,-V2.*k2,DAE.m,DAE.m) + ...
          sparse(a.bus2,a.v2,2*m2.*V2.*g-V1.*k2,DAE.m,DAE.m);

DAE.Gy = DAE.Gy + ...
          sparse(a.bus1,a.bus1, V12.*k3,DAE.m,DAE.m) - ...
          sparse(a.bus1,a.bus2, V12.*k3,DAE.m,DAE.m) + ...
          sparse(a.bus2,a.bus1, V12.*k4,DAE.m,DAE.m) - ...
          sparse(a.bus2,a.bus2, V12.*k4,DAE.m,DAE.m);

DAE.Gy = DAE.Gy + ...
          sparse(a.v1,a.v1,-2.*V1.*b-V2.*k3,DAE.m,DAE.m) + ...
          sparse(a.v1,a.v2,-V1.*k3,DAE.m,DAE.m) + ...
          sparse(a.v2,a.v1, V2.*k4,DAE.m,DAE.m) + ...
          sparse(a.v2,a.v2,-2*m2.*V2.*b+V1.*k4,DAE.m,DAE.m);

DAE.Gy = DAE.Gy - ...
          sparse(a.v1,a.bus1,V12.*k1,DAE.m,DAE.m) + ...
          sparse(a.v1,a.bus2,V12.*k1,DAE.m,DAE.m) + ...
          sparse(a.v2,a.bus1,V12.*k2,DAE.m,DAE.m) - ...
          sparse(a.v2,a.bus2,V12.*k2,DAE.m,DAE.m);