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