(完整版)短路計(jì)算matlab程序

1、% % Power System Analysis' Project: Matlab Program% clear all close allclca=exp(j*2*pi/3);A=1 1 1 ;1 a2 a ;1 a a2;Sref=900e6;Uref=300e3;Iref=Sref/(sqrt(3)*Uref);% data from the different documents:Positive;Negative;Zero;% run power flow:opt = mpoption( 'ENFORCE_Q_LIMS' , 1); results=runp

2、f(case9gs,opt);Vp=results.bus(:,8).*(cos(pi/180*results.bus(:,9)+j*sin(pi/180*result s.bus(:,9);clc;% Menu:Bucle=1;while (Bucle=1)% Menu:disp( '' )disp( '% Program to calculate the symetrical and unsymetrical short circuit faults %');disp( ' 1) Symetrical fault (three phase short

3、 circuit fault)' ); disp( ' 2) Unsymetrical fault (A single phase to ground fault)');disp( ' 0) For exit the program' );obtion=input( 'Choose the option you want to be calculated by typing its number:' , 's' );disp( '' ) % Calculation of the fault switch (

4、obtion)%case '1'disp('You Have choose: SYMETRICAL FAULT');disp(' ' );disp('Please follow the menu');% Ask for the fault's bus number and the Rf value:disp( ' ' );disp(disp(DATA*' ) disp(VoltageTokke= sprintf(Tokke %g*' ) bus=input ( 'Please, nu

5、mber the bus where the three fault occur. Bus number:' );%bus=bus+1; % Since node zero is also analyzed, the bus number X will have a X+1 indexphaseshortcircuitRf=input( 'Please write the value of the arc resistance Rf=' );Rf perphase.% Calculation:If=zeros(8,1);If(bus)=Vp(bus)/(Zpos(bus

6、,bus)+Rf);Vs=Vp-Zpos*If;% Showing the Results: disp( ' ' );disp(disp(*' )RESULTSdisp(FaultCurrent=sprintf('The short circuit currents is: %gA. withangle %g?, Iref*abs(If(bus),180/pi*angle(If(bus);disp(FaultCurrent);figure(); compass(If(bus), 'blue' ); hold oncompass(If(bus)*a

7、, 'red' ); compass(If(bus)*a2, 'green' ); text(-abs(If(bus),-abs(If(bus),);disp( disp ( disp' ' )'The Voltages at all the buses in the system are:' ( '' )disp(' Bus VAngleDegree')disp('Short Circuit current')V %g?,Uref*abs(Vs(1),180/pi*angle(Vs

8、(1);disp(VoltageTokke) VoltageVinje= sprintf('Vinje %gV %g?,Uref*abs(Vs(2),180/pi*angle(Vs(2);disp(VoltageVinje) VoltageSonga= sprintf('Songa %gV %g?,Uref*abs(Vs(3),180/pi*angle(Vs(3);disp(VoltageSonga) VoltageVemork= sprintf('Vemork %gV %g?,Uref*abs(Vs(4),180/pi*angle(Vs(4);disp(Voltage

9、Vemork) VoltageRjukan= sprintf('Rjunkan %gV %g?,Uref*abs(Vs(5),180/pi*angle(Vs(5);disp(VoltageRjukan) VoltageFlesaker=sprintf('Flesaker %gV %g?,Uref*abs(Vs(6),180/pi*angle(Vs(6);disp(VoltageFlesaker) VoltageTveiten= sprintf('Tveiten %gV %g?,Uref*abs(Vs(7),180/pi*angle(Vs(7);disp(VoltageT

10、veiten) VoltageRod= sprintf( ' Rod %gV %g?,Uref*abs(Vs(8),180/pi*angle(Vs(8);disp(VoltageRod) disp ( '' )'Tokkes Voltage vectors');figure() subplot(2,4,1) compass(Vs(1), 'b' ); hold on compass(Vs(1)*a, 'r' ); compass(Vs(1)*a2, 'g' ); text(-abs(Vs(1),-abs(V

11、s(1),subplot(2,4,2) compass(Vs(2); hold on compass(Vs(2)*a, 'r' ); compass(Vs(2)*a2, 'g' );text(-abs(Vs(2),-abs(Vs(2),subplot(2,4,3)compass(Vs(3);hold oncompass(Vs(3)*a, 'r' ); compass(Vs(3)*a2, 'g' ); text(-abs(Vs(3),-abs(Vs(3),'Vinjes Voltage vectors''So

12、ngas Voltage vectors'););subplot(2,4,4) compass(Vs(4); hold on compass(Vs(4)*a, 'r' ); compass(Vs(4)*a2, 'g' );text(-abs(Vs(4),-abs(Vs(4),subplot(2,4,5) compass(Vs(5);hold oncompass(Vs(5)*a, 'r' ); compass(Vs(5)*a2, 'g' ); text(-abs(Vs(5),-abs(Vs(5),'Vemorks V

13、oltage vectors''Rjukans Voltage vectors'););subplot(2,4,6) compass(Vs(6);hold oncompass(Vs(6)*a, 'r' ); compass(Vs(6)*a2, 'g' ); text(-abs(Vs(6),-abs(Vs(6),subplot(2,4,7)compass(Vs(7);hold oncompass(Vs(7)*a, 'r' ); compass(Vs(7)*a2, 'g' ); text(-abs(Vs(7),

14、-abs(Vs(7),subplot(2,4,8) compass(Vs(8); hold on compass(Vs(8)*a, 'r' ); compass(Vs(8)*a2, 'g' );text(-abs(Vs(8),-abs(Vs(8),'Flesakers Voltage vectors''Tveitens Voltage vectors''Rods Voltage vectors');disp(););disp( '*disp()CURRENTSdisp(Q=inputpoints? Y/N:

15、'k=0;whilek=k+1;B1(k)=input(B2(k)=input ()( 'Do, 's' );Q= 'Y'youwant to know any of the currents'Write the first bus number:''Write the second bus number:'if Ypos(B1(k),B2(k)=0' ' ) 'Warning: *There is no connection )' ' )disp(disp( between

16、*'disp(k=k-1;elseIbranch(k)=(Vs(B1(k)-Vs(B2(k)*(-Ypos(B1(k),B2(k);endQ=inputend( 'Do you want to know any more ofthedisp(through););currents?two connectedY/N:' , 's' );disp( '*' )disp(CURRENTS RESULTSdisp(disp( '' ) if k=0disp( 'There are no currents results&#

17、39; end);for K=1:kCurrent=sprintf ( 'From bus no %g to bus no %g the current is: %gA %g?,B1(K),B2(K),Iref*abs(Ibranch(K),180/pi*angle(Ibranch( K);disp (Current)flowingfigure();compass(Ibranch(K), 'blue' );hold on compass(Ibranch(K)*a, 'red' ); compass(Ibranch(K)*a2, 'green

18、9; ); text(-abs(Ibranch(K),-abs(Ibranch(K),selected nodes' );'Current betweenend%case '2'disp('You Have choose: UNSYMETRICAL FAULT' );disp('The program only calcutes single line to ground faults, it isconsider to happen in phase "a"' )disp(' ' )disp(

19、'Please follow the menu');% Ask for the fault's bus number and the Rf value:disp( ' ' );disp(disp(DATA*' )disp( '*) bus=input ('Please, number the bus where the single phase to groundfault occur. Bus number:' );%bus=bus+1; % Since node zero is also analyzed, the b

20、us number X will have a X+1 indexRf=input( 'Please write the value of the arc resistance Rf=' );Rf perphase.% Calculation of the currents:Ifphaseground=3*Vp(bus)/(Zpos(bus,bus)+Zneg(bus,bus)+Zzero(bus,bus)+3*R f);%Fase aIfa=zeros(8,1);Ifa(bus)=Ifphaseground;Ifa1=Ifa/3;Ifa2=Ifa/3;Ifa0=Ifa/3;f

21、or k=1:8Ifault=A*Ifa0(k);Ifa1(k);Ifa2(k); IF(k,:)=Ifault'end% Voltages;%Fase a: V1=Vp-Zpos*Ifa1; V2=-Zneg*Ifa2; V0=-Zzero*Ifa0;for k=1:8U=A*V0(k);V1(k);V2(k) V(k,:)=transpose(U) end% Showing the Results: disp( ' ' );disp(RESULTSdisp(*' ) disp(*' ) disp( ' ' )FaultCurrent=

22、sprintf('The short circuit currents is: %g A. withangle %g ?, Iref*abs(Ifa(bus),180/pi*angle(Ifa(bus);disp(FaultCurrent);figure();compass(IF(bus,1), 'red' ); hold oncompass(IF(bus,2),'blue' );'Short Circuit current' );compass(IF(bus,3),'green' );text(-abs(IF(bus,1

23、),-abs(IF(bus,1),disp( ' ' )disp ( 'The Voltages at all the buses in the system are:'disp( '' )disp( ' Bus Va Vb Vc' );disp(' )VoltageTokke= sprintf( 'Tokke %gV and %g? %gV and %g? %gV and %g?' ,Uref/sqrt(3)*abs(V(1,1),180/pi*angle(V(1,1),Uref/sqrt(3)*abs(

24、V(1,2) ),180/pi*angle(V(1,2),Uref/sqrt(3)*abs(V(1,3),180/pi*angle(V(1,3);d isp(VoltageTokke)VoltageVinje= sprintf( 'Vinje %gV and %g? %gV and %g? %gV and %g?' ,Uref/sqrt(3)*abs(V(1,1),180/pi*angle(V(2,1),Uref/sqrt(3)*abs(V(2,2) ),180/pi*angle(V(2,2),Uref/sqrt(3)*abs(V(2,3),180/pi*angle(V(2,3

25、);d isp(VoltageVinje)VoltageSonga= sprintf( 'Songa %gV and %g? %gV and %g? %gV and %g?' ,Uref/sqrt(3)*abs(V(3,1),180/pi*angle(V(2,1),Uref/sqrt(3)*abs(V(3,2) ),180/pi*angle(V(3,2),Uref/sqrt(3)*abs(V(3,3),180/pi*angle(V(3,3);d isp(VoltageSonga)VoltageVemork= sprintf( 'Vemork %gV and %g? %g

26、V and %g? %gV and %g?' ,Uref/sqrt(3)*abs(V(4,1),180/pi*angle(V(2,1),Uref/sqrt(3)*abs(V(4,2) ),180/pi*angle(V(4,2),Uref/sqrt(3)*abs(V(4,3),180/pi*angle(V(4,3);d isp(VoltageVemork)VoltageRjukan= sprintf( 'Rjukan %gV and %g? %gV and %g? %gV and %g?' ,Uref/sqrt(3)*abs(V(5,1),180/pi*angle(V(2

27、,1),Uref/sqrt(3)*abs(V(5,2) ),180/pi*angle(V(5,2),Uref/sqrt(3)*abs(V(5,3),180/pi*angle(V(5,3);d isp(VoltageRjukan)VoltageFlesaker=sprintf( 'Flesaker %gV and %g? %gV and %g? %gV and %g?' ,Uref/sqrt(3)*abs(V(6,1),180/pi*angle(V(2,1),Uref/sqrt(3)*abs(V(6,2) ),180/pi*angle(V(6,2),Uref/sqrt(3)*ab

28、s(V(6,3),180/pi*angle(V(6,3);d isp(VoltageFlesaker)VoltageTveiten= sprintf( 'Tveiten %gV and %g? %gV and %g? %gV and %g?' ,Uref/sqrt(3)*abs(V(7,1),180/pi*angle(V(2,1),Uref/sqrt(3)*abs(V(7,2) ),180/pi*angle(V(7,2),Uref/sqrt(3)*abs(V(7,3),180/pi*angle(V(7,3);d isp(VoltageTveiten)VoltageRod= sp

29、rintf( 'Rod%gV and %g? %gVand %g? %gV and %g?' ,Uref/sqrt(3)*abs(V(8,1),180/pi*angle(V(2,1),Uref/sqrt(3)*abs(V(8,2) ),180/pi*angle(V(8,2),Uref/sqrt(3)*abs(V(8,3),180/pi*angle(V(8,3);d isp(VoltageRod)disp( '' )figure()subplot(2,4,1)compass(V(1,2),'b');hold oncompass(V(1,1),

30、9;r');'Tokkes Voltage vectors');compass(V(1,3),'g');text(-abs(V(1,2),-abs(V(1,2), subplot(2,4,2) compass(V(2,2); hold oncompass(V(2,1),'r' );compass(V(2,3),'g' );text(-abs(V(2,2),-abs(V(2,2),'Vinjes Voltage vectors');subplot(2,4,3) compass(V(3,2); hold onc

31、ompass(V(3,1),'r' );compass(V(4,3),'g' );text(-abs(V(3,2),-abs(V(3,2),'Songas Voltage vectors');subplot(2,4,4) compass(V(4,2); hold oncompass(V(4,1),'r' );compass(V(4,3),'g' );text(-abs(V(4,2),-abs(V(4,2),'Vemorks Voltage vectors');subplot(2,4,5) compa

32、ss(V(5,2); hold oncompass(V(5,1),'r' );compass(V(5,3),'g' );text(-abs(V(5,2),-abs(V(5,2),'Rjukans Voltage vectors');subplot(2,4,6) compass(V(6,2); hold oncompass(V(6,1),'r' );compass(V(6,3),'g' );text(-abs(V(6,2),-abs(V(6,2),'Flesakers Voltage vectors'

33、);subplot(2,4,7) compass(V(7,2); hold oncompass(V(7,1),'r' );compass(V(7,3),'g' );text(-abs(V(7,2),-abs(V(7,3),'Tveitens Voltage vectors');subplot(2,4,8) compass(V(8,2); hold oncompass(V(8,1),'r' );compass(V(8,3),'g' );text(-abs(V(8,2),-abs(V(8,2),'Rods Vo

34、ltage vectors');disp(disp(disp()CURRENTSdisp( '*Q=input( 'Do you want to know any of the currentspoints? Y/N:', 's' );k=0;whileQ= 'Y'through two connectedk=k+1;B1(k)=input('Write the first bus number:'B2(k)=input ( 'Write the second bus number:' if Ypo

35、s(B1(k),B2(k)=0disp( ' ' ););disp(between*'disp(k=k-1;else'Warning: *There is no connection )' ' )Ibranch0=(V0(B1(k)-V0(B2(k)*(-Yzero(B1(k),(B2(k)Ibranch1=(V1(B1(k)-V1(B2(k)*(-Ypos(B1(k),(B2(k)Ibranch2=(V2(B1(k)-V2(B2(k)*(-Yneg(B1(k),(B2(k)IbranchK=A*Ibranch0;Ibranch1;Ibranch2Ibranch(k,:)=transpose(Ib