牛顿-拉夫逊迭代法电力网潮流计算方法与程序

1、牛顿 -拉夫逊迭代法电力网潮流计算方法与程序编写 佘名寰牛顿-拉夫逊迭代法在电力网潮流计算中因其收敛性较好获得广泛运用，该算法的难点是需反复计算功率方程中雅可比矩阵各个元素表达式。本文简叙了牛顿-拉夫逊迭代法潮流计算的基本公式，通过例题介绍了用牛顿-拉夫逊法计算电力网潮流电压的MATLAB 程序。 程序采用MATLAB 语言的符号矩阵简化了雅可比矩阵系数的计算。本文可供电力系统电气技术人员和大专院校电力类专业师生参考。2 牛顿 - 拉夫逊迭代法潮流计算2.1 极坐标型式节点功率方程式由电源注入或从负载流出节点的电流统称节点电流，由节点电压和节点电流所求得的节点注入功率为：Sti = Vi I?

2、iIi=Z>jVj(I,j=1,2,.n)Pti=Vi ZVMGjcos 由+Bjsin Qij)Qti=V i EnVj (Gj sin 仇j - Bij cos 仇j)(i,j=1,2 .n)(2-01 )式中 Pti,QtI节点注入的有功和无功功率Vi ， Vj 节点 i, j 电压幅值6j= 0i- 0j节点i, j电压的相角差Gij, Bij 节点导纳矩阵的元素，YIJ=GIJ+jB IJ节点功率平衡关系为:Pgi-PlI= P tiQgi -Q li= Q tiAPi=P gi-Pii-Pti=0AQi=Q gi-Q li-Q ti=0(2-02 )Pgi, Qgi 节点i发

3、电机输入有功和无功功率Pli , Qli节点i负荷有功和无功功率APi, AQi-节点i不平衡功率不平衡功率的微分d( APi), d( AQi)d(?Pi)=-(翁?vi+ ?V2t 羡?1+ 缓?2 )(i=1,2,n)d(?qi)=-(察？vi+ U?v2 .+ 等？1+ 学？2 )(2-03 )i? v？ v2? 1? 2对于n个节点系统可得如下矩阵形式修正方程式；?P1 ?1? P1 ?2?P1 ?n? P1 ?V1? P1 ?V2?P1 ?Vn?P1?P2? P2? P2? P2? P2?1?2c ?n?V1?V2? ?Vn八,1?P2?:?2?Pn? Pn?Pn? Pn? Pn?

4、Pn:?Pn?1?2?n?V1?V2? ?Vn?n?Q1 =?Q1?Q1?Q1?Q1?Q1? ?Q19w(2-04 )1,V11Q C?1?2?n?V1?V2? ?Vn?V, Q2?q2? q2?q2? q2? q2?q2, V2?Qn?1?2?n?V1?V2? ?Vn?：?Vn?Qn?Qn?Qn?Qn?Qn?Qn?1?2?n?V1?V2?Vn n-1-个,式中偏微分矩阵为雅可比矩阵。若系统总共有n个节点，其中r个为PV节点，去掉平衡节点则电压角度。变量有电压幅值v变量有n-1-r个，方程式共有2n-r-2个。2.2 牛顿 -拉夫逊迭代法解算过程；指定各节点电压的初值；v1,v2, -vn;

5、?1,?2,？n.一般设V=1, 9=0 ；将所给定的电压值代入式(2-01 ), (2-02 )计算各节点不平衡功率；将电压初值代入求雅可比矩阵系数；解方程式(2-04 ),求节点电压幅值和相角的修正值？Vi，?i；用所求得的电压修正值修正各节点电压：v(1 ) =v(0) + ?v(0)?(1) = ?(0) + ?(0)用新的节点电压重复以上五个步骤计算直到计算精度合符要求。牛顿 -拉夫逊迭代法计算的难点是求雅可比矩阵系数，本文采用MATLAB 语言符号矩阵运算jacobian(f,x)命令求解函数f(x)对x的雅可比矩阵。2.3 牛顿 -拉夫逊迭代法计算程序我们通过例题说明牛顿-拉夫逊

6、迭代法计算潮流程序。【例 2.3.1 】 网络接线如图2-1 所示。支路阻抗分别为Z12=j0.1,Z 23=j0.1,Z 31=j0.1; 三条支路两端的对地电纳皆是j0.01 。给定的注入功率分别为：S1=-1.7192-j0.7346;S2=0.6661 。节点 1 是负荷节点即PQ 节点； 节点 2 给定的电压大小为V2=1.05, 其属于 PV 节点； 节点 3是平衡节点，其电压保持定值，大小为V3=1.0 。Z131V3=1.0P1=1.7192Qi=0.7346Z21=j0.1Z23=j0.1Z13=j0.1YC1=J0.02YC2=J0.02YC3=J0.02-可编辑修改-图2

7、-1例2.3.1网络接线图2.3.1.1 程序参变量说明Y=G+j*B节点导纳矩阵，G为电导，B为电纳；P节点有功功率给定值，正为发电机注入，负为负荷功率；q节点无功功率给定值，正为发电机注入，负为负荷功率；u 节点电压幅值的符号矩阵；delt节点电压相角的符号矩阵；pt节点有功功率符号表达式；qt节点无功功率符号表达式；pp节点有功功率不平衡值符号表达式；qq节点无功功率不平衡值符号表达式；uu节点电压幅值数值矩阵；dd节点电压相角数值矩阵；PP节点不平衡功率的数值矩阵N1网络独立节点总数；N2网络PV节点总数;Sm 节点功率矩阵；Smn 支路功率矩阵；J1, J2, J 节点不平衡功率雅可

8、比符号矩阵JJ 节点不平衡功率雅可比数值矩阵2.3.1.2 源程序*源程序如下：%* NUE78.m example 2-3-1 fig 2-1%The following Program for load flow calculation is based on MATLAB7.0 clear u=sym('u1,u2,u3');delt=sym('d1,d2,d3');G=zeros(3);B=-19.98,10,10;10,-19.98,10;10,10,-19.98;Y=G+j*B;p(1)=-1.7192;q(1)=-0.7346;p(2)=0.666

9、1;k=0;precision=1;N1=2;%the N1 is the amount of the PQ and PV busN2=1;%the N2 is the amount of the PV busfor m=1:N1for n=1:N1+1pt(n)=u(m)*u(n)*(G(m,n)*cos(delt(m)-delt(n)+B(m,n)*sin(delt(m)-delt(n);endpp(m)=p(m)-sum(pt);end for m=1:N1-N2for n=1:N1+1qt(n)=u(m)*u(n)*(G(m,n)*sin(delt(m)-delt(n)-B(m,n)*

10、cos(delt(m)-delt(n);endqq(m)=q(m)-sum(qt);endJ1=jacobian(pp,'d1,d2,u1');J2=jacobian(qq,'d1,d2,u1');J=vertcat(J1,J2);uu=1.0,1.05,1.0;dd=0,0,0;while precision>0.00001u1=uu(1);u2=uu(2);u3=uu(3);d1=dd(1);d2=dd(2);d3=dd(3);for m=1:N1PP(m)=eval(pp(m);endfor m=1:N1-N2PP(N1+m)=eval(qq(m);

11、endJJ=eval(J);-可编辑修改-du=-inv(JJ)*PP'precision=max(abs(du);for n=1:N1dd(n)=dd(n)+du(n);endfor n=1:N1-N2uu(n)=uu(n)+du(N1+n);endk=k+1;endk-1,dd,uu%*% the following program is used to calculate the Sm and S mnfor n=1:N1+1U(n)=uu(n)*(cos(dd(n)+j*sin(dd(n);endUm=conj(U');I=Y*Um;Sm=diag(Um)*conj(I

12、)for m=1:N1+1for n=1:N1+1Smn(m,n)=U(m)*(conj(U(m)-conj(U(n)*conj(-Y(m,n);endend2.3.1.3 程序说明节点导纳矩阵： Y11 =Y22=Y 33 =0+j(-10-10+0.02)=0-j19.98Y12 =Y 23=Y 13 =0+j10.0节点电压符号矩阵：一共三个节点，独立节点数N1=2 , PV节点数N2=1设节点电压幅值符号矩阵为u=sym('u1,u2,u3');节点电压相角符号矩阵为delt=sym('d1,d2,d3');p(1) , q(1) , p(2)为节点功率

13、给定值，节点2为PV节点，q(2)待求；pt(n),pp(m)相关语句是按公式(2-01),(2-02)建立APi符号表达式；APi=pp(1)=-2149/1250-10*u1*u2*sin(d1-d2)-10*u1*u3*sin(d1-d3);AP2=pp(2)= 6661/10000+10*u1*u2*sin(d1-d2)-10*u2*u3*sin(d2-d3)qt(n),qq(m)相关语句是建立 AQi符号表达式，只针对 PQ节点；AQ1 = qq(1)= -3673/5000-999/50*u1A2+10*u1*u2*cos(d1-d2)+10*u1*u3*cos(d1-d3) J1

14、 ， J2， J 语句是建立公式(2-03 )中的雅可比符号矩阵；J =-10*u1*u2*cos(d1-d2)-10*u1*u3*cos(d1-d3),10*u1*u2*cos(d1-d2),-10*u2*sin(d1-d2)-10*u3*sin(d1-d3)10*u1*u2*cos(d1-d2),-10*u1*u2*cos(d1-d2)-10*u2*u3*cos(d2-d3),10*u2*sin(d1-d2)-10*u1*u2*sin(d1-d2)-10*u1*u3*sin(d1-d3), 10*u1*u2*sin(d1-d2), -999/25*u1+10*u2*cos(d1-d2)+1

15、0*u3*cos(d1-d3) while 循环语句是将所假设的节点电压幅值和相角代入符号矩阵，求出节点不平衡功率和雅可比数值矩阵PP, JJ ,解方程(2-04 )，求出节点电压幅值和相角的修正值？Vi，?i,反复迭代,直到精度合符要求。 Um , Smn 语句计算节点注入功率和支路功率。2.3.1.4 程序计算结果：各次迭代节点电压幅值(uu) 和相角 (dd) 的变化：Kuu 1uu 2uu3dd 1dd 2dd 310.98901.05001.0000-0.0909-0.0137020.98501.05001.0000-0.0923-0.0138030.98491.05001.0000-0.0923-0.01380迭代次数k-1=3节点电压相角(弧度)dd =节点 123-0.0923-0.01380节点电压幅值(标么值)uu =节点123节点功率矩阵（正号表示功率注入）0.98491.05001.0000Sm =节点1-1.7192 - 0.7346i20.6661 + 1.2189i31.0531 - 0.3265i支路功率（正号表示功率方向与支路方向一致）Smn =节点1-0.8111 - 0.6089i-0.9081 - 0.1063i20.8111 + 0.7149i0-0.1450 + 0.5260i30.9081