課 程 設(shè) 計(jì) 任 務(wù) 書

1、 山東交通學(xué)院電力系統(tǒng)分析課程設(shè)計(jì)報(bào)告書 院(部)別 信息科學(xué)與電氣工程學(xué)院 班 級 電氣114 學(xué) 號 110812436 姓 名 趙文兵 指導(dǎo)教師 栗玉霞 時 間 2014.06.9-2013.06.13 課 程 設(shè) 計(jì) 任 務(wù) 書題 目 電力系統(tǒng)分析課程設(shè)計(jì) 題目C 學(xué) 院 信息科學(xué)與電氣工程學(xué)院 專 業(yè) 電氣工程及其自動化 班 級 電氣114（電力方向） 學(xué)生姓名 趙文兵 學(xué) 號 110812436 6 月 9 日至 6 月 13 日 共 1 周 6 月 16 日至 6 月 20 日 共 1 周指導(dǎo)教師(簽字) 負(fù)責(zé)人(簽字) 年 月 日一、 設(shè)計(jì)內(nèi)容及要求1.設(shè)計(jì)目的 (1). 掌握

2、電力系統(tǒng)潮流計(jì)算的基本原理； (2). 掌握并能熟練運(yùn)用一門計(jì)算機(jī)語言（MATLAB語言或FORTRAN或C語言或C+語言）； (3). 采用計(jì)算機(jī)語言對潮流計(jì)算進(jìn)行計(jì)算機(jī)編程計(jì)算。2 .設(shè)計(jì)要求 (1). 程序源代碼； (2). 給定題目的輸入，輸出文件； ( 3). 程序說明； (4). 給定系統(tǒng)的程序計(jì)算過程； 3 .設(shè)計(jì)內(nèi)容 (1).根據(jù)電力系統(tǒng)網(wǎng)絡(luò)推導(dǎo)電力網(wǎng)絡(luò)數(shù)學(xué)模型，寫出節(jié)點(diǎn)導(dǎo)納矩陣； (2).賦予各節(jié)點(diǎn)電壓變量（直角坐標(biāo)系形式）初值后，求解不平衡量； (3).形成雅可比矩陣； (4).求解修正量后，重新修改初值，從2開始重新循環(huán)計(jì)算； ( 5).求解的電壓變量達(dá)到所要求的精度時，

3、再計(jì)算各支路功率分布、功率損耗和平衡節(jié) 點(diǎn)功率; (6).上機(jī)編程調(diào)試；連調(diào)； ( 7).計(jì)算分析給定系統(tǒng)潮流分析并與手工計(jì)算結(jié)果作比較分析。 (8).準(zhǔn)備計(jì)算機(jī)演示答辯，書寫該課程設(shè)計(jì)說明書.二、 設(shè)計(jì)原始資料題目C：試用牛頓-拉夫遜潮流計(jì)算方法求解教材P88圖4-10網(wǎng)絡(luò)的潮流分布。求解精度為10e-6。求解思路 此電力系統(tǒng)是一個5節(jié)點(diǎn)，4支路的電力網(wǎng)絡(luò)。其中包含3個PQ節(jié)點(diǎn)，一個PV節(jié)點(diǎn)，和一個平衡節(jié)點(diǎn)。綜合比較牛頓拉夫遜法（直角坐標(biāo)、極坐標(biāo)）、PQ分解法等多種求解方法的特點(diǎn)，最后確定采用牛頓拉夫遜法（極坐標(biāo)）。因?yàn)榇朔椒ㄋ杞獾姆匠探M最少。三、設(shè)計(jì)完成后提交的文件和圖表1計(jì)算說明書部分

4、附錄有設(shè)計(jì)報(bào)告及手算結(jié)果2圖紙部分：（1）設(shè)計(jì)計(jì)算程序以及程序注解（2）計(jì)算詳細(xì)結(jié)果可見附錄四、進(jìn)程安排1.第一天上午：選題，查資料，制定設(shè)計(jì)方案；2.第一天下午第三天下午：手算完成潮流計(jì)算的要求；3.第四天上午第五天上午：編程完成潮流計(jì)算，并對照手算結(jié)果，分析誤差4.第五天下午：答辯，交設(shè)計(jì)報(bào)告五、主要參考資料1.電力系統(tǒng)分析（第三版） 于永源主編，中國電力出版社，2007年2.電力系統(tǒng)分析何仰贊 溫增銀編著，華中科技大學(xué)出版社，2002年版；3.電力系統(tǒng)分析韓楨祥主編，浙江大學(xué)出版社，2001年版；4.電力系統(tǒng)穩(wěn)態(tài)分析陳珩 編，水利電力出版社； 成績評定表指導(dǎo)教師成績答辯小組成績總評成績

5、摘要潮流計(jì)算是電力系統(tǒng)最基本最常用的計(jì)算。根據(jù)系統(tǒng)給定的運(yùn)行條件，網(wǎng)絡(luò)接線及元件參數(shù)，通過潮流計(jì)算可以確定各母線的電壓，包括電壓的幅值和相角，各元件流過的功率，整個系統(tǒng)的功率損耗等一系列數(shù)據(jù)。牛頓拉夫遜Newton-Raphson法是數(shù)學(xué)上解非線性方程組的有效方法，有較好的收斂性。將N-R法用于潮流計(jì)算是以導(dǎo)納矩陣為基礎(chǔ)的，由于利用了導(dǎo)納矩陣的對稱性，稀疏性及節(jié)點(diǎn)編號順序優(yōu)劃等技巧，使N-R法在收斂性，占用內(nèi)存，計(jì)算速度等方面的優(yōu)點(diǎn)都超過了阻抗法。本文運(yùn)用MATLAB軟件進(jìn)行潮流計(jì)算，對給定題目進(jìn)行分析計(jì)算，并應(yīng)用了電力系統(tǒng)潮流計(jì)算仿真軟件DDTR與利用程序計(jì)算的結(jié)果進(jìn)行比較，使計(jì)算的結(jié)果更加

6、準(zhǔn)確。利用成形的程序?qū)ο到y(tǒng)中出現(xiàn)的各種情況，例如負(fù)荷的變化以及線路上所發(fā)生的變化進(jìn)行計(jì)算，并對母線上不滿足范圍的數(shù)據(jù)進(jìn)行調(diào)控，最終進(jìn)行合理的系統(tǒng)潮流。目錄1潮流計(jì)算11.1潮流計(jì)算概述11.2潮流計(jì)算的要求21.3潮流計(jì)算的優(yōu)勢21.4潮流計(jì)算的用途32 MATLAB簡介42.1MATLAB概述43牛頓-拉夫遜法概述53.1牛頓-拉夫遜基本原理53.2直角牛頓-拉夫遜法潮流計(jì)算求解過程64設(shè)計(jì)過程124.1題目C124.2牛頓拉夫遜法程序流程圖244.3設(shè)計(jì)程序25心得體會28參考文獻(xiàn)29MATLAB概述目前電子計(jì)算機(jī)已廣泛應(yīng)用于電力系統(tǒng)的分析計(jì)算，潮流計(jì)算是其基本應(yīng)用之一?，F(xiàn)有很多潮流計(jì)算方

7、法。對潮流計(jì)算方法有五方面的要求：（1）計(jì)算速度快；（2）內(nèi)存需要少；（3）計(jì)算結(jié)果有良好的可靠性和可信性；（4）適應(yīng)性好，亦即能處理變壓器變比調(diào)整、系統(tǒng)元件的不同描述和與其它程序配合的能力強(qiáng)；（5）簡單。MATLAB是一種交互式、面向?qū)ο蟮某绦蛟O(shè)計(jì)語言，廣泛應(yīng)用于工業(yè)界與學(xué)術(shù)界，主要用于矩陣運(yùn)算，同時在數(shù)值分析、自動控制模擬、數(shù)字信號處理、動態(tài)分析、繪圖等方面也具有強(qiáng)大的功能。MATLAB程序設(shè)計(jì)語言結(jié)構(gòu)完整，且具有優(yōu)良的移植性，它的基本數(shù)據(jù)元素是不需要定義的數(shù)組。它可以高效率地解決工業(yè)計(jì)算問題，特別是關(guān)于矩陣和矢量的計(jì)算。MATLAB與C語言和FORTRAN語言相比更容易被掌握。通過M語言

8、，可以用類似數(shù)學(xué)公式的方式來編寫算法，大大降低了程序所需的難度并節(jié)省了時間，從而可把主要的精力集中在算法的構(gòu)思而不是編程上。另外，MATLAB提供了一種特殊的工具：工具箱（TOOLBOXES）.這些工具箱主要包括：信號處理（SIGNALPROCESSING）、控制系統(tǒng)（CONTROLSYSTEMS）、神經(jīng)網(wǎng)絡(luò)（NEURALNETWORKS）、模糊邏輯(FUZZYLOGIC)、小波(WAVELETS)和模擬（SIMULATION）等等。不同領(lǐng)域、不同層次的用戶通過相應(yīng)工具的學(xué)習(xí)和應(yīng)用，可以方便地進(jìn)行計(jì)算、分析及設(shè)計(jì)工作。MATLAB設(shè)計(jì)中，原始數(shù)據(jù)的填寫格式是很關(guān)鍵的一個環(huán)節(jié)，它與程序使用的方便

9、性和靈活性有著直接的關(guān)系。原始數(shù)據(jù)輸入格式的設(shè)計(jì)，主要應(yīng)從使用的角度出發(fā)，原則是簡單明了，便于修改。牛頓-拉夫遜基本原理潮流計(jì)算的目標(biāo)是求取電力系統(tǒng)在給定運(yùn)行狀態(tài)的計(jì)算。即節(jié)點(diǎn)電壓和功率分布，用以檢查系統(tǒng)各元件是否過負(fù)荷。各點(diǎn)電壓是否滿足要求，功率的分布和分配是否合理以及功率損耗等。對現(xiàn)有電力系統(tǒng)的運(yùn)行和擴(kuò)建，對新的電力系統(tǒng)進(jìn)行規(guī)劃設(shè)計(jì)以及對電力系統(tǒng)進(jìn)行靜態(tài)和暫態(tài)穩(wěn)定分析都是以潮流計(jì)算為基礎(chǔ)。潮流計(jì)算結(jié)果可用如電力系統(tǒng)穩(wěn)態(tài)研究，安全估計(jì)或最優(yōu)潮流等對潮流計(jì)算的模型和方法有直接影響。實(shí)際電力系統(tǒng)的潮流技術(shù)那主要采用牛頓-拉夫遜法。牛頓-拉夫遜法(簡稱牛頓法)在數(shù)學(xué)上是求解非線性代數(shù)方程式的有效方

10、法。其要點(diǎn)是把非線性方程式的求解過程變成反復(fù)地對相應(yīng)的線性方程式進(jìn)行求解的過程。即通常所稱的逐次線性化過程。對于非線性代數(shù)方程組： 即 (3-1-1)在待求量x的某一個初始估計(jì)值附近，將上式展開成泰勒級數(shù)并略去二階及以上的高階項(xiàng)，得到如下的經(jīng)線性化的方程組： (3-1-2)上式稱之為牛頓法的修正方程式。由此可以求得第一次迭代的修正量 (3-1-3)將和相加，得到變量的第一次改進(jìn)值。接著就從出發(fā)，重復(fù)上述計(jì)算過程。因此從一定的初值出發(fā)，應(yīng)用牛頓法求解的迭代格式為： (3-1-4) (3-1-5)上兩式中：是函數(shù)對于變量x的一階偏導(dǎo)數(shù)矩陣，即雅可比矩陣J;k為迭代次數(shù)。有上式可見，牛頓法的核心便是

11、反復(fù)形式并求解修正方程式。牛頓法當(dāng)初始估計(jì)值和方程的精確解足夠接近時，收斂速度非?？?，具有平方收斂特性。牛頓潮流算法突出的優(yōu)點(diǎn)是收斂速度快，若選擇到一個較好的初值，算法將具有平方收斂特性，一般迭代45次便可以收斂到一個非常精確的解。而且其迭代次數(shù)與所計(jì)算網(wǎng)絡(luò)的規(guī)模基本無關(guān)。牛頓法也具有良好的收斂可靠性，對于對以節(jié)點(diǎn)導(dǎo)納矩陣為基礎(chǔ)的高斯法呈病態(tài)的系統(tǒng)，牛頓法也能可靠收斂。牛頓法所需的內(nèi)存量及每次迭代所需時間均較高斯法多。牛頓法的可靠收斂取決于有一個良好的啟動初值。如果初值選擇不當(dāng)，算法有可能根本不收斂或收斂到一個無法運(yùn)行的節(jié)點(diǎn)上。對于正常運(yùn)行的系統(tǒng)，各節(jié)點(diǎn)電壓一般均在額定值附近，偏移不會太大，并

12、且各節(jié)點(diǎn)間的相位角差也不大，所以對各節(jié)點(diǎn)可以采用統(tǒng)一的電壓初值(也稱為平直電壓)，如假定： 或 (3-1-6) 這樣一般能得到滿意的結(jié)果。但若系統(tǒng)因無功緊張或其它原因?qū)е码妷嘿|(zhì)量很差或有重載線路而節(jié)點(diǎn)間角差很大時，仍用上述初始電壓就有可能出現(xiàn)問題。解決這個問題的辦法可以用高斯法迭代12次，以此迭代結(jié)果作為牛頓法的初值。也可以先用直流法潮流求解一次以求得一個較好的角度初值，然后轉(zhuǎn)入牛頓法迭代。直角牛頓-拉夫遜法潮流計(jì)算求解過程以下討論的是用直角坐標(biāo)形式的牛頓拉夫遜法潮流的求解過程。當(dāng)采用直角坐標(biāo)時，潮流問題的待求量為各節(jié)點(diǎn)電壓的實(shí)部和虛部兩個分量由于平衡節(jié)點(diǎn)的電壓向量是給定的，因此待求兩共需要2

13、(n-1)個方程式。事實(shí)上，除了平衡節(jié)點(diǎn)的功率方程式在迭代過程中沒有約束作用以外，其余每個節(jié)點(diǎn)都可以列出兩個方程式。對PQ節(jié)點(diǎn)來說，是給定的，因而可以寫出 （3-2-1）對PV節(jié)點(diǎn)來說，給定量是，因此可以列出 （3-2-2）求解過程大致可以分為以下步驟：（1）形成節(jié)點(diǎn)導(dǎo)納矩陣Y（2）將各節(jié)點(diǎn)電壓設(shè)初值U，（3）將節(jié)點(diǎn)初值代入式(2-2-1)或式(2-2-2)，求出修正方程式的常數(shù)項(xiàng)向量（4）將節(jié)點(diǎn)電壓初值代入求式，求出雅可比矩陣元素（5）求解修正方程，求修正向量（6）求取節(jié)點(diǎn)電壓的新值（7）檢查是否收斂，如不收斂，則以各節(jié)點(diǎn)電壓的新值作為初值自第3步重新開始進(jìn)行狹義次迭代，否則轉(zhuǎn)入下一步（8）

14、計(jì)算支路功率分布，PV節(jié)點(diǎn)無功功率和平衡節(jié)點(diǎn)柱入功率。以直角坐標(biāo)系形式表示. 迭代推算式 采用直角坐標(biāo)時,節(jié)點(diǎn)電壓相量及復(fù)數(shù)導(dǎo)納可表示為: (3-2-3)將以上二關(guān)系式代入上式中,展開并分開實(shí)部和虛部;假定系統(tǒng)中的第1,2,m號為PQ節(jié)點(diǎn),第m+1,m+2,n-1為PV節(jié)點(diǎn),根據(jù)節(jié)點(diǎn)性質(zhì)的不同,得到如下迭代推算式: 對于PQ節(jié)點(diǎn) (3-2-4)對于PV節(jié)點(diǎn) (3-2-5)對于平衡節(jié)點(diǎn) 平衡節(jié)點(diǎn)只設(shè)一個,電壓為已知,不參見迭代,其電壓為: (3-2-6). 修正方程式迭代式共包括2(n-1)個方程.選定電壓初值及變量修正量符號,代入方程并按泰勒級數(shù)展開,略去二次方程及以后各項(xiàng),得到一組線性方程組

15、或線性化了的方程組，常稱修正方程組: (3-2-7) (3-2-8).雅可比矩陣各元素的算式式(3-2-8)中, 雅可比矩陣中的各元素可通過對式(3-2-4)和(3-2-5)進(jìn)行偏導(dǎo)而求得.當(dāng)時, 雅可比矩陣中非對角元素為 (3-2-9)當(dāng)時,雅可比矩陣中對角元素為: (3-2-10)由式(2-2-9)和(2-2-10)看出,雅可比矩陣的特點(diǎn):矩陣中各元素是節(jié)點(diǎn)電壓的函數(shù),在迭代過程中,這些元素隨著節(jié)點(diǎn)電壓的變化而變化;導(dǎo)納矩陣中的某些非對角元素為零時,雅可比矩陣中對應(yīng)的元素也是為零.若,則必有;雅可比矩陣不是對稱矩陣;雅可比矩陣各元素的表示如下式(2-2-11): 牛頓拉夫遜法程序流程圖潮流

16、計(jì)算的計(jì)算機(jī)算法程序：%本程序 的功能是用牛頓拉夫遜法進(jìn)行潮流計(jì)算% B1矩陣：1、支路首端號；2、末端號；3、支路阻抗；4、支路首段對地電納% 5、支路的變比；6、支路首端處于K側(cè)為1，1側(cè)為0不含變壓器；% 7、支路末段對地電納% B2矩陣：1、該節(jié)點(diǎn)發(fā)電機(jī)功率；2、該節(jié)點(diǎn)負(fù)荷功率；3、節(jié)點(diǎn)電壓初始值% 4、PV節(jié)點(diǎn)電壓V的給定值；5、節(jié)點(diǎn)所接的無功補(bǔ)償設(shè)備的容量% 6、節(jié)點(diǎn)分類標(biāo)號：1為平衡節(jié)點(diǎn)（應(yīng)為1號節(jié)點(diǎn)）；2為PQ節(jié)點(diǎn)；% 3為PV節(jié)點(diǎn)；clear;n=4;%input(請輸入節(jié)點(diǎn)數(shù)：n=);nl=4;%input(請輸入支路數(shù)：nl=);isb=1;%input(請輸入平衡母線節(jié)

17、點(diǎn)號：isb=);pr=0.00001;%input(請輸入誤差精度：pr=);B1=1 2 0.08+0.40i 0.33333i 1 0 0.02941i; 2 3 0.10+0.40i 0.02941i 1 0 0.03448i; 1 3 0.12+0.50i 0.33333i 1 0 0.03448i; 3 4 0.30i 0 0.9 1 0B2=0 0 1.05 0 0 1; 0 0.55+0.13i 1 0 0 2; 0 0.3+0.18i 1 0 0 2; 0.5 0 1.1 1.1 0 3 ;%input(請輸入各節(jié)點(diǎn)參數(shù)形成的矩陣： B2=);Y=zeros(n);e=zer

18、os(1,n);f=zeros(1,n);V=zeros(1,n);sida=zeros(1,n);S1=zeros(nl);% % %-for i=1:nl%支路數(shù) if B1(i,6)=0%左節(jié)點(diǎn)處于1側(cè) p=B1(i,1);q=B1(i,2); else %左節(jié)點(diǎn)處于K側(cè) p=B1(i,2);q=B1(i,1); end Y(p,q)=Y(p,q)-1./(B1(i,3)*B1(i,5);%非對角元 Y(q,p)=Y(p,q); %非對角元 Y(q,q)=Y(q,q)+1./(B1(i,3)*B1(i,5)2)+B1(i,7);%對角元K側(cè) Y(p,p)=Y(p,p)+1./B1(i,3

19、)+B1(i,4);%對角元1側(cè) end%求導(dǎo)納矩陣disp(導(dǎo)納矩陣 Y=);disp(Y)%-G=real(Y);B=imag(Y);%分解出導(dǎo)納陣的實(shí)部和虛部 for i=1:n%給定各節(jié)點(diǎn)初始電壓的實(shí)部和虛部 e(i)=real(B2(i,3); f(i)=imag(B2(i,3); V(i)=B2(i,4);%PV節(jié)點(diǎn)電壓給定模值 endfor i=1:n%給定各節(jié)點(diǎn)注入功率 S(i)=B2(i,1)-B2(i,2); %i節(jié)點(diǎn)注入功率SG-SL B(i,i)=B(i,i)+B2(i,5);%i節(jié)點(diǎn)無功補(bǔ)償量 end%=P=real(S);Q=imag(S); %分解出各節(jié)點(diǎn)注入的有

20、功和無功功率ICT1=0;IT2=1;N0=2*n;N=N0+1;a=0; %迭代次數(shù)ICT1、a；不滿足收斂要求的節(jié)點(diǎn)數(shù)IT2while IT2=0 % N0=2*n 雅可比矩陣的階數(shù)；N=N0+1擴(kuò)展列 IT2=0;a=a+1; for i=1:n if i=isb%非平衡節(jié)點(diǎn) C(i)=0;D(i)=0; for j1=1:n C(i)=C(i)+G(i,j1)*e(j1)-B(i,j1)*f(j1); %(Gij*ej-Bij*fj) D(i)=D(i)+G(i,j1)*f(j1)+B(i,j1)*e(j1); %(Gij*fj+Bij*ej) end P1=C(i)*e(i)+f(i

21、)*D(i); %節(jié)點(diǎn)功率P計(jì)算ei(Gij*ej-Bij*fj)+fi(Gij*fj+Bij*ej) Q1=C(i)*f(i)-e(i)*D(i); %節(jié)點(diǎn)功率Q計(jì)算fi(Gij*ej-Bij*fj)-ei(Gij*fj+Bij*ej)%求i節(jié)點(diǎn)有功和無功功率P,Q的計(jì)算值 V2=e(i)2+f(i)2;%電壓模平方 %= 以下針對非PV節(jié)點(diǎn)來求取功率差及Jacobi矩陣元素 = if B2(i,6)=3%非PV節(jié)點(diǎn) DP=P(i)-P1;%節(jié)點(diǎn)有功功率差 DQ=Q(i)-Q1; %節(jié)點(diǎn)無功功率差 %= 以上為除平衡節(jié)點(diǎn)外其它節(jié)點(diǎn)的功率計(jì)算 =%= 求取Jacobi矩陣 = for j1=1

22、:n if j1=isb&j1=i%非平衡節(jié)點(diǎn)&非對角元 X1=-G(i,j1)*e(i)-B(i,j1)*f(i);% dP/de=-dQ/df X2=B(i,j1)*e(i)-G(i,j1)*f(i);% dP/df=dQ/de X3=X2; % X2=dp/df X3=dQ/de X4=-X1; % X1=dP/de X4=dQ/df p=2*i-1;q=2*j1-1; J(p,q)=X3;J(p,N)=DQ;m=p+1; % X3=dQ/de J(p,N)=DQ節(jié)點(diǎn)無功功率差 J(m,q)=X1;J(m,N)=DP;q=q+1; % X1=dP/de J(m,N)=DP節(jié)點(diǎn)有功功率差

23、J(p,q)=X4;J(m,q)=X2; % X4=dQ/df X2=dp/df elseif j1=i&j1=isb%非平衡節(jié)點(diǎn)&對角元 X1=-C(i)-G(i,i)*e(i)-B(i,i)*f(i); % dP/de X2=-D(i)+B(i,i)*e(i)-G(i,i)*f(i); % dP/df X3=D(i)+B(i,i)*e(i)-G(i,i)*f(i); % dQ/de X4=-C(i)+G(i,i)*e(i)+B(i,i)*f(i); % dQ/df p=2*i-1;q=2*j1-1;J(p,q)=X3;J(p,N)=DQ; %擴(kuò)展列Q m=p+1; J(m,q)=X1;q=

24、q+1;J(p,q)=X4;J(m,N)=DP; %擴(kuò)展列P J(m,q)=X2; end end else %= 下面是針對PV節(jié)點(diǎn)來求取Jacobi矩陣的元素 = DP=P(i)-P1;% PV節(jié)點(diǎn)有功誤差 DV=V(i)2-V2;% PV節(jié)點(diǎn)電壓誤差 for j1=1:n if j1=isb&j1=i%非平衡節(jié)點(diǎn)&非對角元 X1=-G(i,j1)*e(i)-B(i,j1)*f(i); % dP/de X2=B(i,j1)*e(i)-G(i,j1)*f(i); % dP/df X5=0;X6=0; p=2*i-1;q=2*j1-1;J(p,q)=X5;J(p,N)=DV; % PV節(jié)點(diǎn)電壓

25、誤差 m=p+1; J(m,q)=X1;J(m,N)=DP;q=q+1;J(p,q)=X6; % PV節(jié)點(diǎn)有功誤差 J(m,q)=X2; elseif j1=i&j1=isb %非平衡節(jié)點(diǎn)&對角元 X1=-C(i)-G(i,i)*e(i)-B(i,i)*f(i); % dP/de X2=-D(i)+B(i,i)*e(i)-G(i,i)*f(i); % dP/df X5=-2*e(i); X6=-2*f(i); p=2*i-1;q=2*j1-1;J(p,q)=X5;J(p,N)=DV; % PV節(jié)點(diǎn)電壓誤差 m=p+1; J(m,q)=X1;J(m,N)=DP;q=q+1;J(p,q)=X6;

26、% PV節(jié)點(diǎn)有功誤差 J(m,q)=X2; end end end end end%= 以上為求雅可比矩陣的各個元素及擴(kuò)展列的功率差或電壓差 = for k=3:N0 % N0=2*n （從第三行開始，第一、二行是平衡節(jié)點(diǎn)） k1=k+1;N1=N; % N=N0+1 即 N=2*n+1擴(kuò)展列P、Q 或 U for k2=k1:N1% 從k+1列的Jacobi元素到擴(kuò)展列的P、Q 或 U J(k,k2)=J(k,k2)./J(k,k); % 用K行K列對角元素去除K行K列后的非對角元素進(jìn)行規(guī)格化 end J(k,k)=1; % 對角元規(guī)格化K行K列對角元素賦1 %= 回代運(yùn)算 = if k=3

27、 % 不是第三行 k 3 k4=k-1; for k3=3:k4% 用k3行從第三行開始到當(dāng)前行的前一行k4行消去 for k2=k1:N1 % k3行后各行上三角元素 J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去運(yùn)算（當(dāng)前行k列元素消為0） end %用當(dāng)前行K2列元素減去當(dāng)前行k列元素乘以第k行K2列元素 J(k3,k)=0; %當(dāng)前行第k列元素已消為0 end if k=N0 %若已到最后一行 break; end %= 前代運(yùn)算 = for k3=k1:N0 % 從k+1行到2*n最后一行 for k2=k1:N1 % 從k+1列到擴(kuò)展列消去k+1行后各行

28、下三角元素 J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去運(yùn)算 end %用當(dāng)前行K2列元素減去當(dāng)前行k列元素乘以第k行K2列元素 J(k3,k)=0; %當(dāng)前行第k列元素已消為0 end else %是第三行k=3 %= 第三行k=3的前代運(yùn)算 = for k3=k1:N0 %從第四行到2n行（最后一行） for k2=k1:N1 %從第四列到2n+1列（即擴(kuò)展列） J(k3,k2)=J(k3,k2)-J(k3,k)*J(k,k2);%消去運(yùn)算（當(dāng)前行3列元素消為0） end %用當(dāng)前行K2列元素減去當(dāng)前行3列元素乘以第三行K2列元素 J(k3,k)=0; %當(dāng)前

29、行第3列元素已消為0 end end end%=上面是用線性變換方式高斯消去法將Jacobi矩陣化成單位矩陣= for k=3:2:N0-1 L=(k+1)./2; e(L)=e(L)-J(k,N); %修改節(jié)點(diǎn)電壓實(shí)部 k1=k+1; f(L)=f(L)-J(k1,N); %修改節(jié)點(diǎn)電壓虛部 end %-修改節(jié)點(diǎn)電壓- for k=3:N0 DET=abs(J(k,N); if DET=pr %電壓偏差量是否滿足要求 IT2=IT2+1; %不滿足要求的節(jié)點(diǎn)數(shù)加1 end end ICT2(a)=IT2; %不滿足要求的節(jié)點(diǎn)數(shù) ICT1=ICT1+1; %迭代次數(shù)end%用高斯消去法解w=-

30、J*Vdisp(迭代次數(shù)：);disp(ICT1);disp(沒有達(dá)到精度要求的個數(shù)：);disp(ICT2);for k=1:n V(k)=sqrt(e(k)2+f(k)2); %計(jì)算各節(jié)點(diǎn)電壓的模值 sida(k)=atan(f(k)./e(k)*180./pi; %計(jì)算各節(jié)點(diǎn)電壓的角度 E(k)=e(k)+f(k)*j; %將各節(jié)點(diǎn)電壓用復(fù)數(shù)表示end%= 計(jì)算各輸出量 =disp(各節(jié)點(diǎn)的實(shí)際電壓標(biāo)幺值E為(節(jié)點(diǎn)號從小到大排列)：);disp(E); %顯示各節(jié)點(diǎn)的實(shí)際電壓標(biāo)幺值E用復(fù)數(shù)表示disp(-);disp(各節(jié)點(diǎn)的電壓大小V為(節(jié)點(diǎn)號從小到大排列)：);disp(V); %顯

31、示各節(jié)點(diǎn)的電壓大小V的模值disp(-);disp(各節(jié)點(diǎn)的電壓相角sida為(節(jié)點(diǎn)號從小到大排列)：);disp(sida); %顯示各節(jié)點(diǎn)的電壓相角for p=1:n C(p)=0; for q=1:n C(p)=C(p)+conj(Y(p,q)*conj(E(q); %計(jì)算各節(jié)點(diǎn)的注入電流的共軛值 end S(p)=E(p)*C(p); %計(jì)算各節(jié)點(diǎn)的功率 S = 電壓 X 注入電流的共軛值enddisp(各節(jié)點(diǎn)的功率S為(節(jié)點(diǎn)號從小到大排列)：);disp(S); %顯示各節(jié)點(diǎn)的注入功率disp(-);disp(各條支路的首端功率Si為(順序同您輸入B1時一致)：);for i=1:n

32、l p=B1(i,1);q=B1(i,2); if B1(i,6)=0 Si(p,q)=E(p)*(conj(E(p)*conj(B1(i,4)+(conj(E(p)*B1(i,5). -conj(E(q)*conj(1./(B1(i,3)*B1(i,5); Siz(i)=Si(p,q); else Si(p,q)=E(p)*(conj(E(p)*conj(B1(i,4)+(conj(E(p)./B1(i,5). -conj(E(q)*conj(1./(B1(i,3)*B1(i,5); Siz(i)=Si(p,q); end disp(Si(p,q); SSi(p,q)=Si(p,q); ZF

33、=S(,num2str(p),num2str(q),)=,num2str(SSi(p,q); disp(ZF); disp(-);enddisp(各條支路的末端功率Sj為(順序同您輸入B1時一致)：);for i=1:nl p=B1(i,1);q=B1(i,2); if B1(i,6)=0 Sj(q,p)=E(q)*(conj(E(q)*conj(B1(i,7)+(conj(E(q)./B1(i,5). -conj(E(p)*conj(1./(B1(i,3)*B1(i,5); Sjy(i)=Sj(q,p); else Sj(q,p)=E(q)*(conj(E(q)*conj(B1(i,7)+(

34、conj(E(q)*B1(i,5). -conj(E(p)*conj(1./(B1(i,3)*B1(i,5); Sjy(i)=Sj(q,p); end disp(Sj(q,p); SSj(q,p)=Sj(q,p); ZF=S(,num2str(q),num2str(p),)=,num2str(SSj(q,p); disp(ZF); disp(-);enddisp(各條支路的功率損耗DS為(順序同您輸入B1時一致)：);for i=1:nl p=B1(i,1);q=B1(i,2); DS(i)=Si(p,q)+Sj(q,p); disp(DS(i); DDS(i)=DS(i); ZF=DS(,n

35、um2str(p),num2str(q),)=,num2str(DDS(i); disp(ZF); disp(-);endfigure(1);subplot(1,2,1);plot(V);xlabel(節(jié)點(diǎn)號);ylabel(電壓標(biāo)幺值);grid on;subplot(1,2,2);plot(sida);xlabel(節(jié)點(diǎn)號);ylabel(電壓角度);grid on;figure(2);subplot(2,2,1);P=real(S);Q=imag(S);bar(P);xlabel(節(jié)點(diǎn)號);ylabel(節(jié)點(diǎn)注入有功);grid on;subplot(2,2,2);bar(Q);xlab

36、el(節(jié)點(diǎn)號);ylabel(節(jié)點(diǎn)注入無功);grid on;subplot(2,2,3);P1=real(Siz);Q1=imag(Siz);bar(P1);xlabel(支路號);ylabel(支路首端注入有功);grid on;subplot(2,2,4);bar(Q1);xlabel(支路號);ylabel(支路首端注入無功);grid on;運(yùn)算結(jié)果和數(shù)據(jù)圖片：B1 = 1.0000 2.0000 0.0800 + 0.4000i 0 + 0.3333i 1.0000 0 0 + 0.0294i 2.0000 3.0000 0.1000 + 0.4000i 0 + 0.0294i 1.0000 0 0 + 0.0345i 1.0000 3.0000 0.1200 + 0.5000i 0 + 0.3333i 1.0000 0 0 + 0.0345i 3.0000 4.0000 0 + 0.3000i 0 0.9000 1.0000 0 B2 = 0 0 1