1、下載可編輯目錄：一、軟件需求說明書 .3二、概要設計說明書 .41、編寫潮流計算程序 .42、數據的輸入測試 .43、運行得出結果 .44、進行實驗結果驗證 .4三、詳細設計說明書 .5.專業 .整理 .下載可編輯1、數據導入模塊 .52、節點導納矩陣模塊 .53、編號判斷模塊 .54、收斂條件判定模塊 .55、雅可比矩陣模塊 .56、迭代計算模塊 .57、計算輸出參數模塊 .5四、程序代.專業 .整理 .下載可編輯碼 .6五、最測試例 .151、輸入結果 .152、輸出結果 .153、結果驗證 .15一、軟件需求說明書本次設計利用 MATLAB/C+/C（使用 MATLAB）編程工具編寫潮流

2、計算，實現對節點電壓和功率分布的求取。潮流方程的求解基本方法是迭代，包括牛頓- 拉夫遜法，以及P-Q 分解法，本次設計采用牛頓迭代法。牛頓迭代法（ Newton'smethod）又稱為牛頓 - 拉夫遜方法，它是牛頓在17.專業 .整理 .下載可編輯世紀提出的一種在實數域和復數域上近似求解方程的方法。多數方程不存在求根公式，因此求精確根非常困難， 甚至不可能， 從而尋找方程的近似根就顯得特別重要。方法使用函數f(x) 的泰勒級數的前面幾項來尋找方程f(x)=0 的根。牛頓迭代法是求方程根的重要方法之一，其最大優點是在方程f(x)=0 的單根附近具有平方收斂，而且該法還可以用來求方程的重根

3、、復根。牛頓迭代法是取x0 之后，在這個基礎上，找到比x0 更接近的方程的跟，一步一步迭代， 從而找到更接近方程根的近似跟。牛頓迭代法是求方程根的重要方法之一，其最大優點是在方程f(x) = 0的單根附近具有平方收斂，而且該法還可以用來求方程的重根、復根。電力系統潮流計算，一般來說，各個母線所供負荷的功率是已知的， 各個節點電壓是未知的 （平衡節點外） 可以根據網絡結構形成節點導納矩陣， 然后由節點導納矩陣列寫功率方程， 由于功率方程里功率是已知的，電壓的幅值和相角是未知的， 這樣潮流計算的問題就轉化為求解非線性方程組的問題了。 為了便于用迭代法解方程組， 需要將上述功率方程改寫成功率平衡方程

4、，并對功率平衡方程求偏導， 得出對應的雅可比矩陣， 給未知節點賦電壓初值，一般為額定電壓，將初值帶入功率平衡方程，得到功率不平衡量，這樣由功率不平衡量、雅可比矩陣、節點電壓不平衡量（未知的）構成了誤差方程，解誤差方程，得到節點電壓不平衡量， 節點電壓加上節點電壓不平衡量構成新的節點電壓初值，將新的初值帶入原來的功率平衡方程，并重新形成雅可比矩陣，然后計算新的電壓不平衡量， 這樣不斷迭代， 不斷修正， 一般迭代三到五次就能收斂。二、概要設計說明書1、編寫潮流計算程序本程序主要分為七個模塊： 數據導入模塊、 節點導納矩陣模塊、 編號判斷模塊、收斂條件判定模塊、雅可比矩陣模塊、迭代計算模塊、計算輸出

5、參數模塊。下圖為潮流迭代框圖。.專業 .整理 .下載可編輯節點編號輸入原始數據數據交互形成節點導納矩陣YB設定初值， k=0計算失配功率Pi、QiMax|Pi,Qi|< ?否k=k+1形成雅可比矩陣 J是解修正方程，得到U(k)、 Q(k)進行收斂后的有關計算輸出結果2、數據的輸入測試本次設計是將電力系統分析課本P88例題進行潮流計算。3、運行得出結果得出電壓、有功功率、無功功率、角度的數據。4、進行實驗結果驗證三、詳細設計說明書1、數據導入模塊：利用 Excel 輸入已知節點、 支路數據，通過“讀取” 將數據導入 MATLAB中。2、節點導納矩陣模塊：.專業 .整理 .下載可編輯利用已

6、知的電阻電抗及導納的值，根據導納的計算公式， 計算出節點的自導納及節點間的互導納的值， 按照節點編號組成導納矩陣， 利用 MATLAB“real ”和“ imag”調用導納矩陣中的實部和虛部，分別形成電導和電納的矩陣。3 、編號判斷模塊：當首節點在變壓器左側，設為1，位于變壓器右側，設為2，既非 1 也非 0 為不含變壓器；節點類型為 PQ時，為 1，節點類型為 PV時，為 2，節點類型為 V? 時，為 3。4 、收斂條件判定模塊：根據節點的類型賦初值， 并進行失配功率的初步計算， 判斷是否符合收- 4斂條件 maxPi ，Qi< 10 。如不符合，則進行后續的計算。5 、雅可比矩陣模塊

7、：根據節點類型確定雅克比矩陣的階數，然后根據 n 維非線性方程組的修正方程求出雅克比矩陣。6、迭代計算模塊：解修正方程，并進行收斂條件判斷， 如不符合條件則進行下一次迭代，以一直到符合條件為止。7、計算輸出參數模塊：- 4當滿足收斂條件 maxPi ，Qi< 10 時，結束迭代計算， 計算電壓、有功功率、無功功率、電壓角度等參數。四、程序代碼clear%清除變量；filename='E:MATLABzhilu.xlsx'a=xlsread(filename);%讀取支路信息；filename='E:MATLABjiedian.xlsx'c=xlsread(

8、filename);%讀取節點信息.專業 .整理 .下載可編輯b=zeros(a(1,7);%定義節點導納矩陣G=zeros(a(1,7);B=zeros(a(1,7);for i1=1:(a(1,7)if (a(i1,5)=0)%首節點在變壓器左側b(a(i1,1),a(i1,1)=b(a(i1,1),a(i1,1)+1/(a(i1,3)+1j*a(i1,4); b(a(i1,2),a(i1,2)=b(a(i1,2),a(i1,2)+1/(a(i1,3)+1j*a(i1,4)*a(i1,6)2);b(a(i1,1),a(i1,2)=b(a(i1,1),a(i1,2)+(-1/(a(i1,3)

9、+1j*a(i1,4)*a(i1,6);b(a(i1,2),a(i1,1)=b(a(i1,2),a(i1,1)+(-1/(a(i1,3)+1j*a(i1,4)*a(i1,6);%進行導納計算G(a(i1,1),a(i1,1)=real(b(a(i1,1),a(i1,1);G(a(i1,2),a(i1,2)=real(b(a(i1,2),a(i1,2);G(a(i1,1),a(i1,2)=real(b(a(i1,1),a(i1,2);G(a(i1,2),a(i1,1)=real(b(a(i1,2),a(i1,1);B(a(i1,1),a(i1,1)=imag(b(a(i1,1),a(i1,1);

10、B(a(i1,2),a(i1,2)=imag(b(a(i1,2),a(i1,2);B(a(i1,1),a(i1,2)=imag(b(a(i1,1),a(i1,2);B(a(i1,2),a(i1,1)=imag(b(a(i1,2),a(i1,1);elseif (a(i1,5)=1)%首節點在變壓器右側b(a(i1,1),a(i1,1)=b(a(i1,1),a(i1,1)+1/(a(i1,3)+1j*a(i1,4)*a(i1,6)2);b(a(i1,2),a(i1,2)=b(a(i1,2),a(i1,2)+1/(a(i1,3)+1j*a(i1,4); b(a(i1,1),a(i1,2)=b(a(

11、i1,1),a(i1,2)+(-1/(a(i1,3)+1j*a(i1,4)*a(i1,6);b(a(i1,2),a(i1,1)=b(a(i1,2),a(i1,1)+(-1/(a(i1,3)+1j*a(i1,4)*a(i1,6); %進行導納計算G(a(i1,1),a(i1,1)=real(b(a(i1,1),a(i1,1);G(a(i1,2),a(i1,2)=real(b(a(i1,2),a(i1,2);G(a(i1,1),a(i1,2)=real(b(a(i1,1),a(i1,2);G(a(i1,2),a(i1,1)=real(b(a(i1,2),a(i1,1);B(a(i1,1),a(i1

12、,1)=imag(b(a(i1,1),a(i1,1);B(a(i1,2),a(i1,2)=imag(b(a(i1,2),a(i1,2);B(a(i1,1),a(i1,2)=imag(b(a(i1,1),a(i1,2);B(a(i1,2),a(i1,1)=imag(b(a(i1,2),a(i1,1);else % 不含變壓器支路 b(a(i1,1),a(i1,1)=b(a(i1,1),a(i1,1)+1/(a(i1,3)+1j*a(i1,4)+1j*a(i1,5); b(a(i1,2),a(i1,2)=b(a(i1,2),a(i1,2)+1/(a(i1,3)+1j*a(i1,4)+1j*a(i1

13、,5); b(a(i1,1),a(i1,2)=b(a(i1,1),a(i1,2)+(-1/(a(i1,3)+1j*a(i1,4);b(a(i1,2),a(i1,1)=b(a(i1,2),a(i1,1)+(-1/(a(i1,3)+1j*a(i1,4);%進行導納計算.專業 .整理 .下載可編輯G(a(i1,1),a(i1,1)=real(b(a(i1,1),a(i1,1);G(a(i1,2),a(i1,2)=real(b(a(i1,2),a(i1,2);G(a(i1,1),a(i1,2)=real(b(a(i1,1),a(i1,2);G(a(i1,2),a(i1,1)=real(b(a(i1,2

14、),a(i1,1);B(a(i1,1),a(i1,1)=imag(b(a(i1,1),a(i1,1);B(a(i1,2),a(i1,2)=imag(b(a(i1,2),a(i1,2);B(a(i1,1),a(i1,2)=imag(b(a(i1,1),a(i1,2);B(a(i1,2),a(i1,1)=imag(b(a(i1,2),a(i1,1);endend%計數各個節點個數PQjd=0;PVjd=0;VOjd=0;for i2=1:a(1,7)if c(i2,6)=1PQjd=PQjd+1;endif c(i2,6)=2PVjd=PVjd+1;endif c(i2,6)=3VOjd=VOjd

15、+1;endend%注入功率的計算zhuru_PQ=zeros(PQjd*2+PVjd,1);k=1;for i3=1:a(1,7)if c(i3,6)=1zhuru_PQ(k,1)=c(i3,2)-c(i3,4);zhuru_PQ(k+a(1,7)-1,1)=0;for i10=1:a(1,7)zhuru_PQ(k+a(1,7)-1,1)=zhuru_PQ(k+a(1,7)-1,1)+B(i3,i10);endzhuru_PQ(k+a(1,7)-1,1)=zhuru_PQ(k+a(1,7)-1,1)-c(i3,5)k=k+1;endif c(i3,6)=2zhuru_PQ(k,1)=c(i3

16、,2)-c(i3,4);k=k+1;.專業 .整理 .下載可編輯endendn=a(1,7);%后續矩陣變量定義angle_u=zeros(PQjd*2+PVjd,1);%迭代角度、電壓矩陣angle_u(n:PQjd*2+PVjd,1)=1;%電壓初值設為1，角度初值設為0delta_au=zeros(PQjd*2+PVjd,1);delta_au(n:PQjd*2+PVjd,1)=0;%角度、電壓修正量向量%失配功率的計算k=1;k_=PQjd+1;delta_PQ= zeros(PQjd*2+PVjd,1);for i4=1:n%PQ節點的計算if c(i4,6)=1delta_PQ(k

17、,1)=zhuru_PQ(k,1);delta_PQ(k+n-1,1)=zhuru_PQ(k+n-1,1);i6=1;for i5=1:nif c(i5,6)=1delta_PQ(k,1)=delta_PQ(k,1)-angle_u(n-1+k)*angle_u(n-1+i6)*(G(i4,i5)*cos(angle_u(k, 1)-angle_u(i6,1)+B(i4,i5)*sin(angle_u(k,1)-angle_u(i6,1);delta_PQ(k+n-1,1)=delta_PQ(k+n-1,1)-angle_u(n-1+k)*angle_u(n-1+i6)*(G(i4,i5)*s

18、in(an gle_u(k,1)-angle_u(i6,1)-B(i4,i5)*cos(angle_u(k,1)-angle_u(i6,1);i6=i6+1;endif c(i5,6)=2delta_PQ(k,1)=delta_PQ(k,1)-angle_u(n-1+k)*c(i5,1)*(G(i4,i5)*cos(angle_u(k,1)-angle _u(i6,1)+B(i4,i5)*sin(angle_u(k,1)-angle_u(i6,1);delta_PQ(k+n-1,1)=delta_PQ(k+n-1,1)-angle_u(n-1+k)*c(i5,1)*(G(i4,i5)*sin(

19、angle_u(k, 1)-angle_u(i6,1)-B(i4,i5)*cos(angle_u(k,1)-angle_u(i6,1);i6=i6+1;endif c(i5,6)=3delta_PQ(k,1)=delta_PQ(k,1)-angle_u(n-1+k)*c(i5,1)*(G(i4,i5)*cos(angle_u(k,1)+B(i4 ,i5)*sin(angle_u(k,1);.專業 .整理 .下載可編輯delta_PQ(k+n-1,1)=delta_PQ(k+n-1,1)-angle_u(n-1+k)*c(i5,1)*(G(i4,i5)*sin(angle_u(k, 1)-B(i

20、4,i5)*cos(angle_u(k,1);endendk=k+1;end%PV節點的計算if c(i4,6)=2delta_PQ(k_,1)=zhuru_PQ(k_,1);i6=1;for i5=1:nif c(i5,6)=1delta_PQ(k_,1)=delta_PQ(k_,1)-c(i4,1)*angle_u(n-1+i6)*(G(i4,i5)*cos(angle_u(k_,1)-a ngle_u(i6,1)+B(i4,i5)*sin(angle_u(k_,1)-angle_u(i6,1);i6=i6+1;endif c(i5,6)=2delta_PQ(k_,1)=delta_PQ(

21、k_,1)-c(i4,1)*c(i5,1)*(G(i4,i5)*cos(angle_u(k_,1)-angle_u(i 6,1)+B(i4,i5)*sin(angle_u(k_,1)-angle_u(i6,1);i6=i6+1;endif c(i5,6)=3delta_PQ(k_,1)=delta_PQ(k_,1)-c(i4,1)*c(i5,1)*(G(i4,i5)*cos(angle_u(k_,1)+B(i4,i5)*sin(angle_u(k_,1);endendk_=k_+1;endendzhongjian_delta_PQ=abs(delta_PQ); %取絕對值jingdu=max(

22、zhongjian_delta_PQ); % 為迭代精度的判斷做賦值準備 while (jingdu>0.0001)%雅各比矩陣的形成J=zeros(2*PQjd+PVjd);%矩陣定義U=ones(n,1);%為方便運算將電壓取出放入一個新的矩陣for i=1:PQjdU(i,1)=angle_u(i+n-1,1);.專業 .整理 .下載可編輯endangle=zeros(n,1);%為方便運算將角度取出放入一個新的矩陣for i=1:n-1angle(i,1)=angle_u(i,1);end%H%H部分生成for i1=1:n-1for i2=1:n-1if i1=i2J(i1,i

23、2)=-U(i1,1)*U(i2,1)*(G(i1,i2)*sin(angle(i1,1)-angle(i2,1)-B(i1,i2)*cos(angle (i1,1)-angle(i2,1);endif i1=i2for i3=1:nif i3=i1J(i1,i2)=J(i1,i2)+U(i1,1)*U(i3,1)*(G(i1,i3)*sin(angle(i1,1)-angle(i3,1)-B(i1,i3)*c os(angle(i1,1)-angle(i3,1);endendendendend%N%N部分生成for i1=1:n-1for i2=1:PQjdif i1=i2J(i1,i2+n

24、-1)=-U(i1,1)*(G(i1,i2)*cos(angle(i1,1)-angle(i2,1)-B(i1,i2)*sin(angle(i1, 1)-angle(i2,1);endif i1=i2J(i1,i2+n-1)=(-2)*U(i1,1)*U(i1,1)*G(i1,i1);for i3=1:nif i3=i1J(i1,i2+n-1)=J(i1,i2+n-1)-U(i3,1)*(G(i1,i3)*cos(angle(i1,1)-angle(i3,1)+B(i1,i3)*s.專業 .整理 .下載可編輯in(angle(i1,1)-angle(i3,1);endendendendend%

25、K%K部分生成for i1=1:PQjdfor i2=1:n-1if i1=i2J(i1+n-1,i2)=U(i1,1)*U(i2,1)*(G(i1,i2)*cos(angle(i1,1)-angle(i2,1)-B(i1,i2)*sin(an gle(i1,1)-angle(i2,1);endif i1=i2for i3=1:nif i3=i1J(i1+n-1,i2)=J(i1+n-1,i2)-U(i1,1)*U(i3,1)*(G(i1,i3)*cos(angle(i1,1)-angle(i3,1)+B( i1,i3)*sin(angle(i1,1)-angle(i3,1);endenden

26、dendend%L%L部分生成for i1=1:PQjdfor i2=1:PQjdif i1=i2J(i1+n-1,i2+n-1)=-U(i1,1)*(G(i1,i2)*sin(angle(i1,1)-angle(i2,1)-B(i1,i2)*cos(angle (i1,1)-angle(i2,1);endif i1=i2J(i1+n-1,i2+n-1)=2*U(i1,1)*U(i1,1)*B(i1,i1);for i3=1:nif i3=i1J(i1+n-1,i2+n-1)=J(i1+n-1,i2+n-1)-U(i3,1)*(G(i1,i3)*sin(angle(i1,1)-angle(i3

27、,1)-B( i1,i3)*cos(angle(i1,1)-angle(i3,1);end.專業 .整理 .下載可編輯endendendenddelta_au=-inv(J)*delta_PQ;%解修正方程angle_u=angle_u+delta_au;%形成新的初始值%再次計算失配功率k=1;k_=PQjd+1;delta_PQ= zeros(PQjd*2+PVjd,1);for i4=1:n%if c(i4,6)=1delta_PQ(k,1)=zhuru_PQ(k,1);delta_PQ(k+n-1,1)=zhuru_PQ(k+n-1,1);i6=1;for i5=1:nif c(i5,

28、6)=1delta_PQ(k,1)=delta_PQ(k,1)-angle_u(n-1+k)*angle_u(n-1+i6)*(G(i4,i5)*cos(angle_u(k, 1)-angle_u(i6,1)+B(i4,i5)*sin(angle_u(k,1)-angle_u(i6,1);delta_PQ(k+n-1,1)=delta_PQ(k+n-1,1)-angle_u(n-1+k)*angle_u(n-1+i6)*(G(i4,i5)*sin(an gle_u(k,1)-angle_u(i6,1)-B(i4,i5)*cos(angle_u(k,1)-angle_u(i6,1);i6=i6+

29、1;endif c(i5,6)=2delta_PQ(k,1)=delta_PQ(k,1)-angle_u(n-1+k)*c(i5,1)*(G(i4,i5)*cos(angle_u(k,1)-angle _u(i6,1)+B(i4,i5)*sin(angle_u(k,1)-angle_u(i6,1);delta_PQ(k+n-1,1)=delta_PQ(k+n-1,1)-angle_u(n-1+k)*c(i5,1)*(G(i4,i5)*sin(angle_u(k, 1)-angle_u(i6,1)-B(i4,i5)*cos(angle_u(k,1)-angle_u(i6,1);i6=i6+1;e

30、ndif c(i5,6)=3delta_PQ(k,1)=delta_PQ(k,1)-angle_u(n-1+k)*c(i5,1)*(G(i4,i5)*cos(angle_u(k,1)+B(i4 ,i5)*sin(angle_u(k,1);.專業 .整理 .下載可編輯delta_PQ(k+n-1,1)=delta_PQ(k+n-1,1)-angle_u(n-1+k)*c(i5,1)*(G(i4,i5)*sin(angle_u(k, 1)-B(i4,i5)*cos(angle_u(k,1);endendk=k+1;end%if c(i4,6)=2delta_PQ(k_,1)=zhuru_PQ(k_,1);i6=1;for i5=1:nif c(i5,6)=1delta_PQ(k_,1)=delta_PQ(k_,1)-c(i4,1)*angle_u(n-1+i6)*(G(i4,i5)*cos(angle_u(k_,1)-a ngle_u(i6,1)+B(i4,i5)*sin(angle_u(