fxcall.m

这是一个很适合研究和学习用的电力系统仿真软件

function Fxcall(a)

global DAE Bus Pod

if ~a.n, return, end

V1 = DAE.V(a.bus1);
V2 = DAE.V(a.bus2);
t1 = DAE.a(a.bus1);
t2 = DAE.a(a.bus2);
ss = sin(t1-t2);
cc = cos(t1-t2);
Tr = a.con(:,9);

x1 = DAE.x(a.x1);
x1_max = a.con(:,10);
x1_min = a.con(:,11);

DB = dbtcsc(a);

DAE.Fx = DAE.Fx - sparse(a.x1,a.x1,1./Tr,DAE.n,DAE.n);
u = (DAE.x(a.x1) <= x1_max & DAE.x(a.x1) >= x1_min);

if a.pod
  DAE.Fx(a.x1(a.pod),Pod.Vs(Pod.tcsc)) = Pod.u(Pod.tcsc).*u(a.pod).*a.con(a.pod,16)./Tr(a.pod);
end

P1vs = DB.*V1.*V2.*ss;  
Q1vs = DB.*V1.*(V1-V2.*cc);
Q2vs = DB.*V2.*(V2-V1.*cc);

DAE.Gx = DAE.Gx + sparse(a.bus1,a.x1,u.*P1vs,2*Bus.n,DAE.n);  
DAE.Gx = DAE.Gx - sparse(a.bus2,a.x1,u.*P1vs,2*Bus.n,DAE.n);
DAE.Gx = DAE.Gx + sparse(a.bus1+Bus.n,a.x1,u.*Q1vs,2*Bus.n,DAE.n);
DAE.Gx = DAE.Gx + sparse(a.bus2+Bus.n,a.x1,u.*Q2vs,2*Bus.n,DAE.n);    

ty2 = a.con(:,3) == 2;
  
k = find(u.*ty2);
if k

  Kp = a.con(k,12);
  Ki = a.con(k,13);

  V1 = DAE.V(a.bus1);
  V2 = DAE.V(a.bus2);
  t1 = DAE.a(a.bus1);
  t2 = DAE.a(a.bus2);
  ss = sin(t1-t2);
  cc = cos(t1-t2);
  
  a1 = ss(k).*(a.B(k)+a.y(k));
  a3 = V1(k).*V2(k);
  a5 = a3.*cc(k).*(a.B(k)+a.y(k));
  
  DAE.Fy = DAE.Fy + sparse(a.x1(k),a.bus1(k)+Bus.n,Kp.*V2(k).*a1./Tr(k),DAE.n,2*Bus.n);
  DAE.Fy = DAE.Fy + sparse(a.x1(k),a.bus2(k)+Bus.n,Kp.*V1(k).*a1./Tr(k),DAE.n,2*Bus.n);
  DAE.Fy = DAE.Fy + sparse(a.x1(k),a.bus1(k),Kp.*a5./Tr(k),DAE.n,2*Bus.n);
  DAE.Fy = DAE.Fy - sparse(a.x1(k),a.bus2(k),Kp.*a5./Tr(k),DAE.n,2*Bus.n);

  DAE.Fy = DAE.Fy + sparse(a.x2(k),a.bus1(k)+Bus.n,Ki.*V2(k).*a1,DAE.n,2*Bus.n);
  DAE.Fy = DAE.Fy + sparse(a.x2(k),a.bus2(k)+Bus.n,Ki.*V1(k).*a1,DAE.n,2*Bus.n);
  DAE.Fy = DAE.Fy + sparse(a.x2(k),a.bus1(k),Ki.*a5,DAE.n,2*Bus.n);
  DAE.Fy = DAE.Fy - sparse(a.x2(k),a.bus2(k),Ki.*a5,DAE.n,2*Bus.n);
   
  DAE.Fx = DAE.Fx + sparse(a.x1(k),a.x1(k),Kp.*P1vs(k)./Tr(k),DAE.n,DAE.n);
  DAE.Fx = DAE.Fx + sparse(a.x1(k),a.x2(k),1./Tr(k),DAE.n,DAE.n);
  DAE.Fx = DAE.Fx + sparse(a.x2(k),a.x1(k),Ki.*P1vs(k),DAE.n,DAE.n);
  
  c = find(a.con(k,4) == 2);      
  if c
    
    k = k(c);

    DAE.Fy = DAE.Fy - sparse(a.x1(k),a.bus1(k)+Bus.n,Kp(c).*V2(k).*a1(c)./Tr(k),DAE.n,2*Bus.n);
    DAE.Fy = DAE.Fy - sparse(a.x1(k),a.bus2(k)+Bus.n,Kp(c).*V1(k).*a1(c)./Tr(k),DAE.n,2*Bus.n);
    DAE.Fy = DAE.Fy - sparse(a.x1(k),a.bus1(k),Kp(c).*a5(c)./Tr(k),DAE.n,2*Bus.n);
    DAE.Fy = DAE.Fy + sparse(a.x1(k),a.bus2(k),Kp(c).*a5(c)./Tr(k),DAE.n,2*Bus.n);

    DAE.Fy = DAE.Fy - sparse(a.x2(k),a.bus1(k)+Bus.n,Ki(c).*V2(k).*a1(c),DAE.n,2*Bus.n);
    DAE.Fy = DAE.Fy - sparse(a.x2(k),a.bus2(k)+Bus.n,Ki(c).*V1(k).*a1(c),DAE.n,2*Bus.n);
    DAE.Fy = DAE.Fy - sparse(a.x2(k),a.bus1(k),Ki(c).*a5(c),DAE.n,2*Bus.n);
    DAE.Fy = DAE.Fy + sparse(a.x2(k),a.bus2(k),Ki(c).*a5(c),DAE.n,2*Bus.n);

    DAE.Fx = DAE.Fx - sparse(a.x1(k),a.x1(k),Kp(c).*P1vs(k)./Tr(k),DAE.n,DAE.n);
    DAE.Fx = DAE.Fx - sparse(a.x2(k),a.x1(k),Ki(c).*P1vs(k),DAE.n,DAE.n);
    
  end
  
end