电力系统不对称故障计算的Matlab算法程序

NF=1;fori1=1:NFclear;n1=9;n2=9;n0=9;nl1=15;nl2=15;nl0=15;Lf=input('请输入短路类型（单相接地短路输入‘1’，两相短路输入‘2’，两相短路接地输入‘3’）:Lf=');f=7;If=0;zf=0;zg=0;B1=[010.0553i010;020.1228i010;030.1785i010;140.0476i01.0500;050.6514-0.2413i010;060.9174-0.2752i010;270.0474i01.0250;081.0492-0.4590i010;390.0484i01.0250;450.0854i011;460.0907i011;570.1248i011;690.1255i011;780.0696i011;890.0922i011];B2=[010.0553i010;020.1228i010;030.1785i010;140.0476i01.0500;050.6514-0.2413i010;060.9174-0.2752i010;270.0474i01.0250;081.0492-0.4590i010;390.0484i01.0250;450.0854i011;460.0907i011;570.1248i011;690.1255i011;780.0696i011;890.0922i011];B0=[010.0001010;020.0001010;030.0001010;140.0476i01.0500;050.6514-0.2413i010;060.9174-0.2752i010;270.0474i01.0250;081.0492-0.4590i010;390.0484i01.0250;450.0854i*3011;460.0907i*3011;570.1248i*3011;690.1255i*3011;780.0696i*3011;890.0922i*3011];A1=[111;1-1./2-i*sqrt(3)/2-1./2+i*sqrt(3)/2;1-1./2+i*sqrt(3)/2-1./2-i*sqrt(3)/2];Vcs=[1.0000;1.0000;1.0000;1.0346;1.0090;1.0139;1.0284;1.0043;1.0212];Z1=zeros(n1);Z2=zeros(n2);Z0=zeros(n0);Y1=zeros(n1);Y2=zeros(n2);Y0=zeros(n0);V1=zeros(n1,1);V2=zeros(n2,1);V0=zeros(n0,1);I1=zeros(nl1,1);I2=zeros(nl2,1);I0=zeros(nl0,1);form1=1:3m=0;ifm1==1nl=nl1;B=B1;Z=Z1;n=n1;elseifm1==2nl=nl2;B=B2;Z=Z2;n=n2;elseifm1==3nl=nl0;B=B0;Z=Z0;n=n0;endm=0;Z=zeros(n);V=zeros(n);I=zeros(nl);fork1=1:nlp=B(k1,1);q=B(k1,2);ifB(k1,6)==0k=1./B(k1,5);elsek=B(k1,5);endifp==0ifqmZ(q,q)=B(k1,3);m=m+1;elsefori=1:mZ(i,m+1)=-Z(i,q);Z(m+1,i)=-Z(q,i);endZ(m+1,m+1)=Z(q,q)+B(k1,3);fori=1:mforj=1:mZ(i,j)=Z(i,j)-Z(i,m+1)*Z(m+1,j)./Z(m+1,m+1);endZ(i,m+1)=0;endfori=1:m+1Z(m+1,i)=0;endendelseifqmfori=1:mZ(i,q)=Z(i,p)*k;Z(q,i)=Z(p,i)*k;endZ(q,q)=k^2*Z(p,p)+k^2*B(k1,3);m=m+1;elsefori=1:mZ(i,m+1)=k*Z(i,p)-Z(i,q);Z(m+1,i)=k*Z(p,i)-Z(q,i);endZ(m+1,m+1)=k^2*Z(p,p)+Z(q,q)-2*k*Z(p,q)+k^2*B(k1,3);fori=1:mforj=1:mZ(i,j)=Z(i,j)-Z(i,m+1)*Z(m+1,j)./Z(m+1,m+1);endZ(i,m+1)=0;endfori=1:m+1Z(m+1,i)=0;endendendendifm1==1Z1=Z;elseifm1==2Z2=Z;elseifm1==3Z0=Z;endendifIf==0Z1(f,f)=Z1(f,f);Z2(f,f)=Z2(f,f);Z0(f,f)=Z0(f,f);elseZ1(f,f)=Z1(f,f)+Z1(1,1)-2*Z1(f,1);Z2(f,f)=Z2(f,f)+Z2(1,1)-2*Z2(f,1);Z0(f,f)=Z0(f,f)+Z0(1,1)-2*Z0(f,1);endifLf==1Z4=Z2(f,f)+Z0(f,f)+3*zf;k2=1;k0=1;I1(f)=Vcs(f)./(Z1(f,f)+Z4);I2(f)=k2*I1(f);I0(f)=k0*I1(f);a1=I1(f);a2=I2(f);a0=I0(f);elseifLf==2Z4=Z2(f,f)+2*zf;k2=-1;k0=0;I1(f)=Vcs(f)./(Z1(f,f)+Z4);I2(f)=k2*I1(f);I0(f)=k0*I1(f);a1=I1(f);a2=I2(f);a0=I0(f);elseifLf==3Z4=zf+(Z2(f,f)+zf)*(Z0(f,f)+zf+3*zg)./(Z2(f,f)+Z0(f,f)+2*zf+3*zg);k2=-(Z0(f,f)+zf+3*zg)./(Z2(f,f)+Z0(f,f)+2*zf+3*zg);k0=-(Z2(f,f)+zf)./(Z2(f,f)+Z0(f,f)+2*zf+3*zg);I1(f)=Vcs(f)./(Z1(f,f)+Z4);I2(f)=k2*I1(f);I0(f)=k0*I1(f);a1=I1(f);a2=I2(f);a0=I0(f);elseifLf==4Z4=zf+(Z2(f,f)+zf)*(Z0(f,f)+zf)./(Z2(f,f)+Z0(f,f)+2*zf);k2=-(Z0(f,f)+zf)./(Z2(f,f)+Z0(f,f)+2*zf);k0=-(Z2(f,f)+zf)./(Z2(f,f)+Z0(f,f)+2*zf);I1(f)=Vcs(f)./(Z1(f,f)+Z4);I2(f)=k2*I1(f);I0(f)=k0*I1(f);a1=I1(f);a2=I2(f);a0=I0(f);elseifLf==5Z4=Z2(f,f)+Z0(f,f)+3*zf;k2=1;k0=1;I1(f)=Vcs(f)./(Z1(f,f)+Z4);I2(f)=k2*I1(f);I0(f)=k0*I1(f);a1=I1(f);a2=I2(f);a0=I0(f);enddisp('故障点处的正序电流标幺值为：');disp(a1);disp('故障点处的负序电流标幺值为：');disp(a2);disp('故障点处的零序电流标幺值为：');disp(a0);fori=1:n1V1(i)=Vcs(i)-I1(f)*Z1(i,f);endfori=1:n2V2(i)=-I2(f)*Z2(i,f);endfori=1:n0V0(i)=-I0(f)*Z0(i,f);enddisp('各节点的正序电压标幺值为（节点号从小到大排）：');disp(V1);disp('各节点的负序电压标幺值为（节点号从小到大排）：');disp(V2);disp('各节点的零序电压标幺值为（节点号从小到大排）：');disp(V0);fori=1:3ifi==1B=B1;V=V1;nl=nl1;elseifi==2B=B2;V=V2;nl=nl2;elseifi==3B=B0;V=V0;nl=nl0;endforj=1:nlifB(j,6)==0k=B(j,5);elsek=1./B(j,5);endp=B(j,1);q=B(j,2);ifp==0e=0;b=B(j,3);I(j)=(e-V(q)./k)./b;elseI(j)=(V(p)-V(q)./k)./B(j,3);endendifi==1I1=I;elseifi==2I2=I;elseifi==3I0=I;endenddisp('各支路的正序电流标幺值为（顺序同您输入B时一样）;');fori=1:nldisp(I1(i));enddisp('各支路的负序电流标幺值为（顺序同您输入B时一样）;');fori=1:nldisp(I2(i));enddisp('各支路的零序电流标幺值为（顺序同您输入B时一样）;');fori=1:nldisp(I0(i));endfori=1:3ifi==1s0=a0;s1=a1;s2=a2;elseifi==2s0=V0;s1=V1;s2=V2;elses0=I0;s1=I1;s2=I2;endforj=1:3CA=s0*A1(j,1)+s1*A1(j,2)+s2*A1(j,3);ifj==1A=CA;elseifj==2B=CA;elseifj==3C=CA;endendifi==1Iadl=A;Ibdl=B;Icdl=C;elseifi==2Vajy=A;Vbjy=B;Vcjy=C;elseIazl=A;Ibzl=B;Iczl=C;endenddisp('故障点处A相电流Iadl为：');disp(Iadl);disp('故障点处B相电流Ibdl为：');disp(Ibdl);disp('故障点处C相电流Icdl为：');disp(Icdl);disp('各节点处A相电压Vajy分别为(节点号从小到大排）：');disp(Vajy);disp('各节点处B相电压Vbjy分别为(节点号从小到大排）：');disp(Vbjy);disp('各节点处C相电压Vcjy分别为(节点号从小到大排）：');disp(Vcjy);disp('各支路A相电流Iazl分别为（顺序同您输入B时一样');fori=1:nldisp(Iazl(i));enddisp('各支路B相电流Ibzl分别为（顺序同您输入B时一样');fori=1:nldisp(Ibzl(i));enddisp('各支路C相电流Iczl分别为（顺序同您输入B时一样');fori=1:nldisp(Iczl(i));endend