fcall.m

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

function p = fcall(p)

global DAE Pod

if ~p.n, return, end

ty1 = find(p.con(:,3) == 1);
ty2 = find(p.con(:,3) == 2);

x0 = p.x0;
x1 = DAE.x(p.x1);

if ty2    
  [Ps,Qs,Pr,Qr,Is,Ir] = fex_lineflows;
  Kp = p.con(:,12);
  Ki = p.con(:,13);
  tp = ty2(find(p.con(ty2,4) == 1));
  ta = ty2(find(p.con(ty2,4) == 2));   
  if tp
    x0(tp) = -Kp(tp).*(p.Pref(tp)-Ps(p.line(tp))) + DAE.x(p.x2(tp));
    DAE.f(p.x2(tp)) = -Ki(tp).*(p.Pref(tp)-Ps(p.line(tp)));
  end
  if ta
    x0(ta) = -Kp(ta).*(p.Pref(ta)-Ps(p.line(ta))-Pr(p.line(ta))) ...
             + DAE.x(p.x2(ta));
    DAE.f(p.x2(ta)) = -Ki(ta).*(p.Pref(ta)-Ps(p.line(ta))-Pr(p.line(ta)));
  end
end

% Updating signal x0
if p.pod
  x0(p.pod) = x0(p.pod) + p.con(p.pod,16).*DAE.x(Pod.Vs(Pod.tcsc));
end

x1_max = p.con(:,10);
x1_min = p.con(:,11);

u = ~(x1 >= x1_max & DAE.f(p.x1) > 0) & ...
    ~(x1 <= x1_min & DAE.f(p.x1) < 0);
DAE.f(p.x1) = u.*(x0-x1)./p.con(:,9);
DAE.x(p.x1) = min(x1_max,x1);
DAE.x(p.x1) = max(x1_min,DAE.x(p.x1));

% update B
p.B = btcsc(p);