电力系统潮流及短路电流计算程序

班级：姓名：学号：作业要求编写程序计算图1所示算例系统的潮流及三相短路电流。

潮流计算：方法不限，计算系统的节点电压和相角。

短路电流：4号母线发生金属性三相短路时（z f=0），分别按照精确算法和近似算法计算短路电流、系统中各节点电压以及网络中各支路的电流分布，并对两种情况下的计算结果进行比较。

电路图及参数图1 3机9节点系统表1 9节点系统支路参数计算步骤进行系统正常运行状态的潮流计算，求得(0)i U &形成不含发电机和负荷的节点导纳矩阵Y N ；将发电机表示为电流源iI &(/i di E jx ''=&)和导纳i y (1/di jx '=)的并联组合；节点负荷用恒阻抗的接地支路表示；形成包括所有发电机支路和负荷支路的节点导纳矩阵Y ，即在Y N 中的发电机节点和负荷节点的自导纳上分别增加发电机导纳i y 和负荷导纳,LD iy （*,,22LD i LDi LDiLD ii i S P jQ y V V -==）； 利用1Z Y-=，计算节点阻抗矩阵，从而得到阻抗矩阵中的第f 列；利用公式(6-7)或(6-10)计算短路电流；利用公式(6-8)或(6-11)计算系统中各节点电压；利用公式(6-9)计算变压器支路的电流；对输电线路利用П型等值电路计算支路电流。

计算结果节点导纳矩阵Yn ： Columns 1 through 50 0 0 0 + 00 0 0 0 00 0 0 0 00 + 0 0 +0 0 0 +0 0 0 + 00 0 + 0 0 +0 0 0 0 00 0 0 + 0 0 Columns 6 through 90 0 0 00 0 + 0 00 0 0 0 ++ 0 0 00 + 0 00 0 +0 + 00 + ++ 0 +电压幅值：电压相角：节点有功：节点无功：修正后的节点导纳矩阵Y：Columns 1 through 50 0 0 0 + 00 0 0 0 00 0 0 0 00 + 0 0 +0 0 0 +0 0 0 + 00 0 + 0 0 +0 0 0 0 00 0 0 + 0 0 Columns 6 through 90 0 0 00 0 + 0 00 0 0 0 ++ 0 0 00 + 0 00 0 +0 + 00 + ++ 0 +节点阻抗矩阵Z的第4列：+++++++++精确计算结果：短路电流：模值：相角：节点电压模值：支路电流：i j Iij1 4 7 9 5 +4 6 +5 7 +6 9 +7 8 +8 9 +近似计算结果：短路电流：模值：相角：节点电压模值：五、程序流程图六、程序及输入文件文件：文件：l=9;%支路数n=9;%节点数m=6;%PQ节点数Yn=zeros(n);%初始化节点导纳矩阵Y DATA1=xlsread('',1);%计算节点导纳矩阵Yfor k=1:li=DATA1(k,1);j=DATA1(k,2);R=DATA1(k,3);X=DATA1(k,4);B2=DATA1(k,5);Yn(i,i)=Yn(i,i)+1i*B2+1/(R+1i*X);Yn(j,j)=Yn(j,j)+1i*B2+1/(R+1i*X);Yn(i,j)=Yn(i,j)-1/(R+1i*X);Yn(j,i)=Yn(j,i)-1/(R+1i*X);enddisp('节点导纳矩阵Yn：');disp(Yn);G=real(Yn);B=imag(Yn);DATA2=xlsread('',2);P=zeros(1,n);Q=zeros(1,n);U=ones(1,n);P(2:n)=DATA2(2:n,3);Q(4:n)=DATA2(4:n,4);U(1:3)=DATA2(1:3,5);%设置节点电压初值e(1)=DATA2(1,5);e(2:n)=;f(1:n)=;%设置迭代次数t=0;tmax=10;while t<=tmax%计算f(x)a(1:n)=;c(1:n)=;for i=2:nfor j=1:na(i)=a(i)+G(i,j)*e(j)-B(i,j)*f(j);c(i)=c(i)+G(i,j)*f(j)+B(i,j)*e(j);endendfor i=2:ndeltaP(i)=P(i)-e(i)*a(i)-f(i)*c(i);endfor j=4:ndeltaQ(j)=Q(j)-f(j)*a(j)+e(j)*c(j);endfor k=2:3deltaU2(k)=U(k)*U(k)-e(k)*e(k)-f(k)*f(k);endfx=[deltaP(2:n) deltaQ(4:n) deltaU2(2:3)]';%计算雅克比矩阵Jfor i=2:nfor j=2:nif i~=jH(i,j)=-(G(i,j)*e(i)+B(i,j)*f(i));N(i,j)=B(i,j)*e(i)-G(i,j)*f(i);elseH(i,j)=-a(i)-(G(i,i)*e(i)+B(i,i)*f(i));N(i,j)=-c(i)+(B(i,i)*e(i)-G(i,i)*f(i));endendendfor i=4:nfor j=2:nif i~=jM(i,j)=B(i,j)*e(i)-G(i,j)*f(i);L(i,j)=G(i,j)*e(i)+B(i,j)*f(i);elseM(i,j)=c(i)+(B(i,i)*e(i)-G(i,i)*f(i));L(i,j)=-a(i)+(G(i,i)*e(i)+B(i,i)*f(i));endendendfor i=2:3for j=2:nif i~=jR(i,j)=0;S(i,j)=0;elseR(i,j)=-2*e(i);S(i,j)=-2*f(i);endendendJ=[H(2:n,2:n) N(2:n,2:n);M(4:n,2:n) L(4:n,2:n);R(2:3,2:n) S(2:3,2:n)];if max(abs(fx))<%输出结果break;else%求解修正方程获得dxdx=-J^(-1)*fx;dx=dx';e(2:n)=e(2:n)+dx(1:n-1);f(2:n)=f(2:n)+dx(n:2*(n-1));t=t+1;endendif t>tmaxstr='潮流计算不收敛';disp(str);elsea(1:n)=;c(1:n)=;for i=1:nfor j=1:na(i)=a(i)+G(i,j)*e(j)-B(i,j)*f(j);c(i)=c(i)+G(i,j)*f(j)+B(i,j)*e(j);endendfor i=1:nU(i)=e(i)+1i*f(i);amp(i)=abs(U(i));arg(i)=angle(U(i));P(i)=e(i)*a(i)+f(i)*c(i);Q(i)=f(i)*a(i)-e(i)*c(i);enddisp('电压幅值：');disp(amp);disp('电压相角：');disp(arg);disp('节点有功：');disp(P);disp('节点无功：');disp(Q);end%计算短路电流f=4;zf=;%修正节点导纳矩阵Xd=DATA2(1:3,6);E=DATA2(1:3,7);for i=1:3Ii(i)=E(i)/(1i*Xd(i));endY=Yn;for i=1:3Y(i,i)=Y(i,i)+1/(1i*Xd(i));endfor j=4:nY(j,j)=Y(j,j)+(-P(j)+1i*Q(j))/(U(j)*U(j)); enddisp('修正后的节点导纳矩阵Y：');disp(Y);disp('节点阻抗矩阵Z的第4列：');disp(Z(:,4));%精确计算disp('精确计算结果：');U0=U;If=U0(f)/(Z(f,f)+zf);amp=abs(If);arg=atand(imag(If)/real(If));disp('短路电流：');disp('模值：');disp(amp);disp('相角：');disp(arg);for i=1:nU(i)=U0(i)-Z(i,f)*If;amp=abs(U);enddisp('节点电压模值：');disp(amp);disp('支路电流：');str=['i ''j '' Iij'];disp(str);for k=1:li=DATA1(k,1);j=DATA1(k,2);r=DATA1(k,3);x=DATA1(k,4);z=r+1i*x;I=(U(i)-U(j))/z;str=[num2str(i) ' ' num2str(j) ' ' num2str(I)];disp(str);end%近似计算disp('近似计算结果：');U0(1:n)=;If=U0(f)/(Z(f,f)+zf);amp=abs(If);arg=atand(imag(If)/real(If));disp('短路电流：');disp('模值：');disp(amp);disp('相角：');disp(arg);U(i)=U0(i)-Z(i,f)*If;amp=abs(U);enddisp('节点电压模值：'); disp(amp);。