fcall.m

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

function p = fcall(p)

global DAE Pod Cluster

if ~p.n, return, end

% Updating reference voltages
Vref = p.Vref;
Vref(p.pod) = Vref(p.pod) + DAE.x(Pod.Vs(Pod.svc));

if p.cluster
  Vref(p.cluster) = Vref(p.cluster) + DAE.x(Cluster.Vs(Cluster.svc));
end

V = DAE.V(p.bus);

if ~isempty(p.ty1)
  bcv = DAE.x(p.bcv);
  Tr = p.con(p.ty1,6);
  Kr = p.con(p.ty1,7);
  bcv_max = p.con(p.ty1,9);
  bcv_min = p.con(p.ty1,10);
  DAE.f(p.bcv) = (Kr.*(Vref(p.ty1)-V(p.ty1))-bcv)./Tr;
  a = find(bcv >= bcv_max & DAE.f(p.bcv) > 0);
  if ~isempty(a), DAE.f(p.bcv(a)) = 0; end
  a = find(bcv <= bcv_min & DAE.f(p.bcv) < 0);
  if ~isempty(a), DAE.f(p.bcv(a)) = 0; end
  DAE.x(p.bcv) = min(bcv_max,DAE.x(p.bcv));
  DAE.x(p.bcv) = max(bcv_min,DAE.x(p.bcv));
end

if ~isempty(p.ty2)
  alpha = DAE.x(p.alpha);
  vm = DAE.x(p.vm);
  a_max = p.con(p.ty2,9);
  a_min = p.con(p.ty2,10);
  T2 = p.con(p.ty2,6);
  K = p.con(p.ty2,7);
  Kd = p.con(p.ty2,11);
  T1 = p.con(p.ty2,12);
  Km = p.con(p.ty2,13);
  Tm = p.con(p.ty2,14);
  xl = p.con(p.ty2,15);
  xc = p.con(p.ty2,16);
  DAE.f(p.vm) = (Km.*V(p.ty2)-vm)./Tm;
  DAE.f(p.alpha) = -Kd./T2.*alpha + K.*T1./T2./Tm.*(vm-Km.*V(p.ty2)) ...
      + K./T2.*(Vref(p.ty2)-vm);
  a = find(alpha >= a_max & DAE.f(p.alpha) > 0);
  if ~isempty(a), DAE.f(p.alpha(a)) = 0; end
  a = find(alpha <= a_min & DAE.f(p.alpha) < 0);
  if ~isempty(a), DAE.f(p.alpha(a)) = 0; end
  DAE.x(p.alpha) = min(a_max,DAE.x(p.alpha));
  DAE.x(p.alpha) = max(a_min,DAE.x(p.alpha));
end

p.Be = bsvc(p);