电力系统分析大作业-教材

1、电力系统分析程序设计专业班级电力系统分析大作业要求计算所给系统的潮流，设发电机 G1 的端电压为 1 p.u. ，发出的有功、无功可 调；发电机 G2的端电压为 1 p.u. ，按指定的有功 P=0.5 p.u. 发电，取=10-4。SB=100MV，A UB=Uav。G1G2j0.1905j0.19051:1.0522L10.02+j0.061:1.05220.025+j0.08j0.05j0.05j0.03+j0.10.18+j0.120.8055+j0.5320电力系统分析大作业一牛顿拉夫逊法(1) 程序代码clear;clc;%节点编号， 5 号为平衡节点， 4 号为 PV 节点， 1

2、、2、3 号为 PQ 节点 %输入原始数据，求节点导纳矩阵y=0;y(1,1) = 1/(0.03+0.1i)+1/(0.025+0.08i)+0.09i+0.07i;y(1,2) = -1/(0.025+0.08i);y(1,3) = -1/(0.03+0.1i);y(2,1) = -1/(0.025+0.08i);y(2,2) = 1/(0.02+0.06i)+1/(0.025+0.08i)+0.05i+0.07i+1/1.05222*1/(0.1905i);y(2,3) = -1/(0.02+0.06i);y(2,4) = 4.9889i;y(3,1) = -1/(0.03+0.1i);

3、y(3,2) = -1/(0.02+0.06i);y(3,3) = 7.7523-28.7757i;y(3,5) = 4.9889i;y(4,2) = 4.9889i;y(4,4) = -5.2493i;y(5,3) = 4.9889i;电力系统分析大作业y(5,5) = -5.2493i;Y=y; display(导纳矩阵为 : );%Y 为导纳矩阵YG=real(Y);B=imag(Y);%原始节点功率Pg=zeros(1,5);Qg=zeros(1,5);Pd=zeros(1,5);Qd=zeros(1,5);S(1)=-0.8055-0.5320i;S(2)=-0.18-0.12i;S

4、(3)=0;S(4)=0.5;S(5)=0;P=real(S);Q=imag(S);%赋初值U = ones(1,5);a = zeros(1,5);电力系统分析大作业dx=ones(7,1);fx=ones(7,1);count=0;% 计算迭代次数while (max(fx)1e-4)%收敛判据for i=1:4for j=1:4H(i,j)=0;dP(i)=0;endendfor i=1:3for j=1:3L(i,j)=0;dQ(i)=0;endendfor i=1:4for j=1:3N(i,j)=0;endend for i=1:3电力系统分析大作业for j=1:4K(i,j)=

5、0;endend%P%失配功率表达式for i=1:4x=0;for j=1:5%dP(i)=dP(i)-U(i)*U(j)*(G(i,j)*cos(a(i)-a(j) + B(i,j)*sin(a(i)-a(j) ); x=x+U(i)*U(j)*(G(i,j)*cos(a(i)-a(j)+B(i,j)*sin(a(i)-a(j);enddP(i)=P(i)-x;endfor i=1:3y=0;for j=1:5y=y+(U(i)*U(j)*(G(i,j)*sin(a(i)-a(j)-B(i,j)*cos(a(i)-a(j);enddQ(i)=Q(i)-y;end fx=dP,dQ; %求雅

6、克比矩阵电力系统分析大作业% 当 i=j 时候求 H,N,K,L 如下： for i=1:4 for j=1:5 if i=jH(i,i)=H(i,i)+U(i)*U(j)*(G(i,j)*sin(a(i)-a(j)-B(i,j)*cos(a(i)-a(j) );end endendfor i=1:4 for j=1:4if i=jH(i,j)=-U(i)*U(j)*(G(i,j)*sin(a(i)-a(j)-B(i,j)*cos(a(i)-a(j); end end end for i=1:3 x=0;for j=1:5if i=jx=x + (- (U(j)*( G(i,j)*cos(a(

7、i)+a(j) + B(i,j)*sin(a(i)-a(j) ) ); end end电力系统分析大作业N(i,i)=-2*U(i)*G(i,i)+x; end for i=1:4 for j=1:3 if i=jN(i,j)=-U(i)*( G(i, j)*cos(a(i)-a(j) + B(i,j)*sin(a(i)-a(j) ); end end end for i=1:3 for j=1:5 if i=jK(i,i)=K(i,i)-U(i)*U(j)*(G(i,j)*cos(a(i)-a(j) + B(i,j)*sin(a(i)-a(j) ); end end end for i=1:

8、3 for j=1:4 if i=jK(i,j)=U(i)*U(j)*( G(i,j)*cos(a(i)-a(j) + B(i,j)*sin(a(i)-a(j) ); end end电力系统分析大作业endfor i=1:3x=0;for j=1:5 if i=jx=x+(-( U(j)*(G(i,j)*sin(a(i)-a(j) - B(i,j)*cos(a(i)-a(j); end endL(i,i)=2*U(i)*B(i,i)+x;endfor i=1:3for j=1:3 if i=jL(i,j)=-U(i)*( G(i,j)*sin(a(i)-a(j) - B(i,j)*cos(a(

9、i)-a(j) ); end endendJ=H,N;K,L; %J 为雅克比矩阵 dx=-(inv(J)*fx);for i=1:4da(i)=dx(i);end电力系统分析大作业for i=1:3 dU(i)=dx(i+4);endfor i=1:4a(i)=a(i)+da(i); end for i=1:3U(i)=U(i)+dU(i); end count=count+1; end%求节点注入的净功率i=5;for j=1:5P(i)=P(i)+U(i)*U(j)*(G(i,j)*cos(a(i)-a(j)+B(i,j)*sin(a(i)-a(j);Q(i)=Q(i)+U(i)*U(j

10、)*(G(i,j)*sin(a(i)-a(j)-B(i,j)*cos(a(i)-a(j);endS(5)=P(5)+Q(5)*sqrt(-1);if count=3display(最终结果为 : );电力系统分析大作业H,N,K,L display(雅可比矩阵如下 : );Jdisp(/*/)a./pi*180disp(/*/)fx,dx,U,a display（迭代次数为 : ）;count display（节点注入净功率如下 : ）;S%求节点注入电流 display（节电注入电流如下 : ）;I=Y display（迭代次数为 : ）;countend10电力系统分析大作业运算结果如下：

11、导纳矩阵为 :最终结果为 :H =-20.602111.38949.2127011.5923-32.239815.59475.05289.395815.6419-30.1181005.05280-5.0528N =-5.61273.84733.02673.2930-8.72645.15802.50655.0456-8.09920-0.49140K =7.0121-3.9156-3.09650-3.26619.0432-5.2770-0.5001-2.4861-5.13528.1121011电力系统分析大作业L =-19.699511.19079.005011.6878-31.442415.24

12、319.473215.3690-29.4399雅可比矩阵如下 :ans =-7.4732 -5.8532 -5.5850 -0.2006 0fx =-0.0017-0.0021-0.0020-0.0001-0.000312电力系统分析大作业-0.0004-0.0004 dx =-0.0013-0.0013-0.0012-0.0013-0.0003-0.0002-0.0002U =0.9916 1.01751.02291.00001.0000a =-0.1304 -0.1022-0.0975-0.00350迭代次数为 :count =3节点注入净功率如下 :S =-0.8055-0.5320i

13、-0.1800-0.1200i00.5000节电注入电流如下0.4966 + 0.1706i13电力系统分析大作业迭代次数为 : count =3二 P Q 解耦迭代 程序代码 n=5;m=3;% 各参数初始化Y=zeros(n,n); % 导纳矩阵 U=ones(n,1); % 电压矢量 a=zeros(n,1); % 相角矢量 Ps=zeros(n,1);Qs=zeros(n,1);P=zeros(n,1);Q=zeros(n,1); p=zeros(n-1,1);q=zeros(m,1); aa=zeros(n-1,1);14电力系统分析大作业u=zeros(m,1);k=0;%节点导纳

14、矩阵的生成I=sqrt(-1);%由于后面用到 ij 作为循环变量，故重定义虚数单位Y=6.311-I*20.3722 -2.7523+I*9.1743 -3.5587+I*11.3879 0 0;-3.5587+I*11.3879 8.5587-I*30.9993 -5+I*15 I*4.9889 0;-2.7523+I*9.1743 -5+I*15 7.7523-I*28.7557 0 I*4.9889;0 I*4.9889 0 -I*5.2493 0;0 0 I*4.9889 0 -I*5.2493;%初值的设定U=1;1;1;1;1;a=0;0;0;0;0;Ps=-0.8055;-0.

15、18;0;0.5;0;Qs=-0.532;-0.12;0;0;0;%求失配功率 P(k) 和 Q(k)for i=1:n-1s=0;for j=1:ns=s+U(j)*(real(Y(i,j)*cos(a(i)-a(j)+imag(Y(i,j)*sin(a(i)-a(j); endP(i)=U(i)*s;end15电力系统分析大作业for i=1:n-1s=0;for j=1:n s=s+U(j)*(real(Y(i,j)*sin(a(i)-a(j)-imag(Y(i,j)*cos(a(i)-a(j);endQ(i)=U(i)*s;endfor i=1:n-1detp(i)=Ps(i)-P(i

16、);endfor i=1:mdetq(i)=Qs(i)-Q(i);end%生成 B和 B 矩阵Bp=zeros(n-1,n-1);Bpp=zeros(m,m);for i=1:n-1;for j=1:n-1;Bp(i,j)=imag(Y(i,j);end;end;for i=1:m;16电力系统分析大作业for j=1:m;Bpp(i,j)=imag(Y(i,j);end;end;k=0;while (max(detp)=0.0001 | max(detq)=0.0001) % 设定精确度 aa=(-1)*inv(Bp)*detp(1)/U(1);detp(2)/U(2);detp(3)/U(

17、3);detp(4)/U(4); u=(-1)*inv(Bpp)*detq(1)/U(1);detq(2)/U(2);detq(3)/U(3);% 求解修正方程，得到修正量 aa,ufor i=1:n-1 %计算 a(节点电压相角 ),U( 节点电压大小 ) a(i)=a(i)+aa(i);endfor i=1:mU(i)=U(i)+u(i);end% detp;detq % 迭代过程中失配功率变化情况k=k+1;%( 记录迭代次数 )for i=1:n-1 %计算功率误差 p,qs=0;for j=1:ns=s+U(j)*(real(Y(i,j)*cos(a(i)-a(j)+imag(Y(i

18、,j)*sin(a(i)-a(j);end17电力系统分析大作业P(i)=U(i)*s;endfor i=1:n-1s=0;for j=1:n s=s+U(j)*(real(Y(i,j)*sin(a(i)-a(j)-imag(Y(i,j)*cos(a(i)-a(j); endQ(i)=U(i)*s;endfor i=1:n-1detp(i)=Ps(i)-P(i);endfor i=1:mdetq(i)=Qs(i)-Q(i);endendA=a*180/pi;%对相角进行转换(由弧度制装换为角度)display(最终的计算结果如下 : );AUdisplay( 生成的节点导纳矩阵为： );18电力系统分析大作业display（算得矩阵 B为： ）;Bp display（算得矩阵 B”为： ）;Bpp display（ 各节点的电压为 :）;for i=1:nfprintf（ 节点 %d: %d （%d） nn,i,U（i）,A（i）; end d