fxcall.m

电力系统分析计算程序

function Fxcall(p)

global DAE

if ~p.n, return, end

m = DAE.x(p.m);
h = p.con(:,7);
k = p.con(:,8);
mmax = p.con(:,9);
mmin = p.con(:,10);

u = m < mmax & m > mmin & p.u & ~p.con(:,11);

V1 = u.*DAE.y(p.v1);
V2 = u.*DAE.y(p.v2);  
V12 = V1.*V2;
V22 = 2.*m.*V2.*V2;
y = admittance(p);
g = real(y);
b = imag(y);
[s12,c12] = angles(p);

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

DAE.Gx = DAE.Gx + ...
         sparse(p.bus1,p.m, -V12.*k1,DAE.m,DAE.n) + ...
         sparse(p.v1,p.m,-V12.*k3,DAE.m,DAE.n) + ...
         sparse(p.bus2,p.m, V22.*g-V12.*k2,DAE.m,DAE.n) + ...
         sparse(p.v2,p.m,-V22.*b+V12.*k4,DAE.m,DAE.n);

DAE.Fx = DAE.Fx - sparse(p.m,p.m,h+(~u),DAE.n,DAE.n);

% voltage control
u1 = u.*(p.con(:,16) == 1 | p.con(:,16) == 3);
if sum(u1)
  DAE.Fy = DAE.Fy + sparse(p.m,p.vr,u1.*k,DAE.n,DAE.m);
end

% reactive power control
u2 = u.*(p.con(:,16) == 2);
if sum(u2)
  ku = u2.*k.*m;
  a1 = -ku.*(V2.*k4);
  a2 = -ku.*(-2*m.*V2.*b+V1.*k4);
  a3 = -ku.*(V12.*k2);
  
  DAE.Fx = DAE.Fx - sparse(p.m,p.m,ku.*k4,DAE.n,DAE.n);

  DAE.Fy = DAE.Fy + ... 
           sparse(p.m,p.v1,a1,DAE.n,DAE.m) + ...
           sparse(p.m,p.v2,a2,DAE.n,DAE.m) + ...
           sparse(p.m,p.bus1,a3,DAE.n,DAE.m) - ...
           sparse(p.m,p.bus2,a3,DAE.n,DAE.m);
end