fm_fault.m

一个较好的MATLAB潮流程序

function fm_fault(flag,t)
% FM_FAULT define three phase faults at buses
%
% FM_FAULT(FLAG,T)
%     FLAG   0  fault initialization
%            1  fault interventions
%     T      integration time
%
%Author:    Federico Milano
%Date:      11-Nov-2002
%Update:    13-Jan-2003
%Update:    11-Sep-2003
%Version:   1.1.0
%
%E-mail:    [e-mail]
%Web-site:  http://thunderbox.uwaterloo.ca/~fmilano
%
% Copyright (C) 2002-2005 Federico Milano
%
% This toolbox is free software; you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation; either version 2.0 of the License, or
% (at your option) any later version.
%
% This toolbox is distributed in the hope that it will be useful, but
% WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANDABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
% General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with this toolbox; if not, write to the Free Software
% Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307,
% USA.

global Fault DAE Bus Shunt Line jay Settings Syn

% Fault.dat:
%          1.  fault conductance
%          2.  fault susceptance
%          3.  old bus conductance
%          4.  old bus susceptance
%          5.  fault check:
%                          0  if fault is not yet met
%                          1  if fault is met
%                         -1  if fault is ended

switch flag
 case 0 % fault data initializations

  Fault.dat = zeros(Fault.n,5);
  Fault.V = DAE.V;
  Fault.ang = DAE.a;
  if Syn.n
    Fault.delta = mean(DAE.x(Syn.delta));
  else
    Fault.delta = 0;
  end

  x = Fault.con(:,7)+jay*Fault.con(:,8);
  x(find(abs(x) == 0)) = jay*1e-6;

  y = 1./x;
  Fault.dat(:,1) = real(y);
  Fault.dat(:,2) = imag(y);
  Fault.dat(:,3) = Shunt.g(Fault.bus);
  Fault.dat(:,4) = Shunt.b(Fault.bus);

 case 1 % faults computations

  for i = 1:Fault.n

    h = Fault.bus(i);

    if t == Fault.con(i,5) % fault intervention

      fm_disp(['Applying fault(s) at t = ',num2str(t),' s'])

      % update admittance matrix
      Shunt.g(h) = Fault.dat(i,3) + Fault.dat(i,1);
      Shunt.b(h) = Fault.dat(i,4) + Fault.dat(i,2);
      fm_y

      % update algebraic variables
      conv = fm_nrlf(40, 1e-4);

    elseif t == Fault.con(i,6) % fault clearance

      fm_disp(['Clearing fault(s) at t = ',num2str(t),' s'])

      % update admittance matrix
      Shunt.g(h) = Fault.dat(i,3);
      Shunt.b(h) = Fault.dat(i,4);
      fm_y

      if Syn.n
        DAE.a = mean(DAE.x(Syn.delta))-Fault.delta+Fault.ang;
      else
        DAE.a = Fault.ang;
      end

      % update algebraic variables
      conv = fm_nrlf(40, 1e-4);

    end

  end

end