Symfault.m
上传用户:eighthdate
上传日期:2014-05-24
资源大小:270k
文件大小:4k
源码类别:

其他行业

开发平台:

Matlab

  1. % The program symfault is designed for the balanced three-phase
  2. % fault analysis of a power system network. The program requires
  3. % the bus impedance matrix Zbus. Zbus may be defined by the
  4. % user, obtained by the inversion of Ybus or it may be
  5. % determined either from the function Zbus = zbuild(zdata)
  6. % or the function Zbus = zbuildpi(linedata, gendata, yload).
  7. % The program prompts the user to enter the faulted bus number
  8. % and the fault impedance Zf. The prefault bus voltages are
  9. % defined by the reserved Vector V. The array V may be defined or
  10. % it is returned from the power flow programs lfgauss, lfnewton,
  11. % decouple or perturb. If V does not exist the prefault bus voltages
  12. % are automatically set to 1.0 per unit. The program obtains the
  13. % total fault current, the postfault bus voltages and line currents.
  14. %
  15. % Copyright (C) 1998 H. Saadat
  16. function symfaul(zdata, Zbus, V)
  17. nl = zdata(:,1); nr = zdata(:,2); R = zdata(:,3);
  18. X = zdata(:,4);
  19. nc = length(zdata(1,:));
  20.   if nc > 4
  21.   BC = zdata(:,5);
  22.   elseif nc ==4, BC = zeros(length(zdata(:,1)), 1);
  23.   end
  24. ZB = R + j*X;
  25. nbr=length(zdata(:,1)); nbus = max(max(nl), max(nr));
  26. if exist('V') == 1
  27.     if length(V) == nbus
  28.     V0 = V;
  29.     else, end
  30. else,  V0 = ones(nbus, 1) + j*zeros(nbus, 1);
  31. end
  32. fprintf('Three-phase balanced fault analysis n')
  33. ff = 999;
  34. while ff > 0
  35. nf = input('Enter Faulted Bus No. -> ');
  36.      while nf <= 0 | nf > nbus
  37.      fprintf('Faulted bus No. must be between 1 & %g n', nbus)
  38.      nf = input('Enter Faulted Bus No. -> ');
  39.      end
  40. fprintf('nEnter Fault Impedance Zf = R + j*X in ')
  41. Zf = input('complex form (for bolted fault enter 0). Zf = ');
  42. fprintf('  n')
  43. fprintf('Balanced three-phase fault at bus No. %gn', nf)
  44. If = V0(nf)/(Zf + Zbus(nf, nf));
  45. Ifm = abs(If); Ifmang=angle(If)*180/pi;
  46. fprintf('Total fault current = %8.4f per unit nn', Ifm)
  47. %fprintf(' p.u. nn', Ifm)
  48. fprintf('Bus Voltages during fault in per unit nn')
  49. fprintf('     Bus     Voltage       Anglen')
  50. fprintf('     No.     Magnitude     degreesn')
  51. for n = 1:nbus
  52.     if n==nf
  53.     Vf(nf) = V0(nf)*Zf/(Zf + Zbus(nf,nf)); Vfm = abs(Vf(nf)); angv=angle(Vf(nf))*180/pi;
  54.     else, Vf(n) = V0(n) - V0(n)*Zbus(n,nf)/(Zf + Zbus(nf,nf));
  55.     Vfm = abs(Vf(n)); angv=angle(Vf(n))*180/pi;
  56.     end
  57.     fprintf('   %4g',  n), fprintf('%13.4f', Vfm),fprintf('%13.4fn', angv)
  58. end
  59. fprintf('  n')
  60. fprintf('Line currents for fault at bus No.  %gnn', nf)
  61. fprintf('     From      To     Current     Anglen')
  62. fprintf('     Bus       Bus    Magnitude   degreesn')
  63. for n= 1:nbus
  64.    %Ign=0;
  65.     for I = 1:nbr
  66.       if nl(I) == n | nr(I) == n
  67.          if nl(I) ==n       k = nr(I);
  68.          elseif nr(I) == n  k = nl(I);
  69.          end
  70.            if k==0
  71.             Ink = (V0(n) - Vf(n))/ZB(I);
  72.             Inkm = abs(Ink); th=angle(Ink);
  73.                %if th <= 0
  74.                 if real(Ink) > 0
  75.                 fprintf('      G   '), fprintf('%7g',n), fprintf('%12.4f', Inkm)
  76.                 fprintf('%12.4fn', th*180/pi)
  77.                 elseif real(Ink) ==0 & imag(Ink) < 0
  78.                 fprintf('      G   '), fprintf('%7g',n), fprintf('%12.4f', Inkm)
  79.                 fprintf('%12.4fn', th*180/pi)
  80.                else, end
  81.            Ign=Ink;
  82.            elseif k ~= 0
  83.             Ink = (Vf(n) - Vf(k))/ZB(I)+BC(I)*Vf(n);
  84.             %Ink = (Vf(n) - Vf(k))/ZB(I);
  85.             Inkm = abs(Ink); th=angle(Ink);
  86.             %Ign=Ign+Ink;
  87.             %if th <= 0
  88.              if real(Ink) > 0
  89.                   fprintf('%7g', n), fprintf('%10g', k),
  90.                   fprintf('%12.4f', Inkm), fprintf('%12.4fn', th*180/pi)
  91.                   elseif real(Ink) ==0 & imag(Ink) < 0
  92.                   fprintf('%7g', n), fprintf('%10g', k),
  93.                   fprintf('%12.4f', Inkm), fprintf('%12.4fn', th*180/pi)
  94.                   else, end
  95.               else, end
  96.       else, end
  97.    end
  98.  if n==nf
  99.  fprintf('%7g',n), fprintf('         F'),  fprintf('%12.4f', Ifm)
  100.  fprintf('%12.4fn', Ifmang)
  101.  else, end
  102. end
  103.    resp=0;
  104.    while strcmp(resp, 'n')~=1 & strcmp(resp, 'N')~=1 & strcmp(resp, 'y')~=1 & strcmp(resp, 'Y')~=1
  105.    resp = input('Another fault location? Enter ''y'' or ''n'' within single quote -> ');
  106.    if strcmp(resp, 'n')~=1 & strcmp(resp, 'N')~=1 & strcmp(resp, 'y')~=1 & strcmp(resp, 'Y')~=1
  107.    fprintf('n Incorrect reply, try again nn'), end
  108.    end
  109.    if resp == 'y' | resp == 'Y'
  110.    nf = 999;
  111.    else ff = 0; end
  112. end   % end for while