第一篇：基于MATLAB的電力系統(tǒng)潮流計(jì)算設(shè)計(jì)

關(guān)鍵詞：電力系統(tǒng)分析；潮流計(jì)算；matlab仿真

中圖分類(lèi)號(hào)：tm744 文獻(xiàn)標(biāo)識(shí)碼：a 文章編號(hào)：1006-4311（2016）21-0185-03

0 引言

潮流計(jì)算是電力系統(tǒng)穩(wěn)態(tài)運(yùn)行中的基本計(jì)算方法中的一種計(jì)算方法，也是電力系統(tǒng)穩(wěn)態(tài)運(yùn)行中最重要的運(yùn)算。潮流計(jì)算是保證電力系統(tǒng)安全、經(jīng)濟(jì)運(yùn)行的根本。在新電網(wǎng)建設(shè)的初期規(guī)劃中，有了潮流計(jì)算，可規(guī)劃出電源的容量及其接入點(diǎn)，可計(jì)算出無(wú)功補(bǔ)償?shù)娜萘?，選擇合適的補(bǔ)償方式，以滿(mǎn)足在電網(wǎng)潮流的控制、調(diào)壓、調(diào)相、調(diào)峰的交換要求。潮流計(jì)算可以選擇電力系統(tǒng)的運(yùn)行方式，便于定期對(duì)電力系統(tǒng)中的元件進(jìn)行檢修。

潮流計(jì)算的過(guò)程

1.1 原始資料

①系統(tǒng)圖：兩個(gè)發(fā)電廠(chǎng)分別通過(guò)變壓器和輸電線(xiàn)路與四個(gè)變電所相連。（圖1）

②發(fā)電廠(chǎng)資料：

③變電所資料：

1）變電所1、2、3、4低壓母線(xiàn)的電壓等級(jí)分別為：10kv，35kv，10kv，35kv。

3）每個(gè)變電所的功率因數(shù)均為cosφ=0.9。

④輸電線(xiàn)路資料：

發(fā)電廠(chǎng)和變電所之間的輸電線(xiàn)路的電壓等級(jí)及長(zhǎng)度標(biāo)于圖中，單位長(zhǎng)度的電阻為0.17ω，單位長(zhǎng)度的電抗為0.402ω，單位長(zhǎng)度的電納為2.78*10-6s。

1.2 基本要求

①對(duì)給定的網(wǎng)絡(luò)查找潮流計(jì)算所需的各元件等值參數(shù)，畫(huà)出等值電路圖。

②輸入各支路數(shù)據(jù)，各節(jié)點(diǎn)數(shù)據(jù)，利用simulink搭建仿真模型等方法，進(jìn)行在變電所的某一負(fù)荷情況下的潮流計(jì)算及仿真，并對(duì)計(jì)算結(jié)果進(jìn)行分析。

③如果各母線(xiàn)電壓不滿(mǎn)足要求，進(jìn)行電壓的調(diào)整。（變電所低壓母線(xiàn)電壓10kv要求調(diào)整范圍在9.5-10.5之間；電壓35kv要求調(diào)整范圍在35-36之間）。

④利用matlab軟件，進(jìn)行上述各種情況潮流的計(jì)算及仿真。

1.3 節(jié)點(diǎn)設(shè)置及分析

由上述系統(tǒng)圖可知，該系統(tǒng)圖為雙端供電網(wǎng)絡(luò)。將母線(xiàn)1，2設(shè)為節(jié)點(diǎn)1，10，將變電所1、2、3、4的高低壓側(cè)分別設(shè)為節(jié)點(diǎn)2、3、4、5、6、7、8、9。并且，將節(jié)點(diǎn)1設(shè)為平衡節(jié)點(diǎn)，將節(jié)點(diǎn)10設(shè)為pv節(jié)點(diǎn)，其余節(jié)點(diǎn)設(shè)為pq節(jié)點(diǎn)。

1.4 參數(shù)求取

將參數(shù)整理如表

1、表2所示。

1.5 進(jìn)行潮流計(jì)算

圖2為仿真模型圖。

從潮流計(jì)算的結(jié)果可得到，系統(tǒng)的各個(gè)節(jié)點(diǎn)電壓的標(biāo)幺值可歸納為表3。

由matlab編程調(diào)節(jié)后，可得到表4的發(fā)電廠(chǎng)電壓和變壓器分接頭電壓得標(biāo)幺值。

在得到了上述調(diào)節(jié)后的電壓標(biāo)幺值，對(duì)電機(jī)模型和變壓器模型進(jìn)行更改。表5為調(diào)節(jié)前后各節(jié)點(diǎn)的電壓標(biāo)幺值。

由題意可知，變電所低壓母線(xiàn)電壓10kv要求調(diào)整范圍在9.5-10.5之間；電壓35kv要求調(diào)整范圍在35-36之間。因此我們可以看出，經(jīng)過(guò)調(diào)節(jié)后，節(jié)點(diǎn)3、5、7、9點(diǎn)電壓已經(jīng)滿(mǎn)足了系統(tǒng)的要求。表6是電壓調(diào)節(jié)前后對(duì)線(xiàn)路損耗進(jìn)行分析的記錄。

由表6的電壓調(diào)節(jié)前后功率損耗對(duì)比，可以看出有功功率隨著變壓器分接頭變比的增大而逐漸增大，使得變壓器的低壓側(cè)的電壓處于允許范圍內(nèi)，符合其要求。

表7為調(diào)節(jié)后的各支路電壓首末端的功率整理表

表8為各節(jié)點(diǎn)功率s的標(biāo)幺值。

1.6 對(duì)比

由上面的三種方法簡(jiǎn)單地比較，我們可以看出，在同一個(gè)電力系統(tǒng)中，用不同的方法進(jìn)行潮流計(jì)算，所得到的結(jié)果是大致相同的。

結(jié)束語(yǔ)

第二篇：用matlab電力系統(tǒng)潮流計(jì)算

題目：潮流計(jì)算與matlab

教學(xué)單位 電氣信息學(xué)院

姓 名 學(xué) 號(hào)

年 級(jí)

專(zhuān) 業(yè) 電氣工程及其自動(dòng)化

指導(dǎo)教師

職 稱(chēng) 副教授

摘 要

電力系統(tǒng)穩(wěn)態(tài)分析包括潮流計(jì)算和靜態(tài)安全分析。本文主要運(yùn)用的事潮流計(jì)算，潮流計(jì)算是電力網(wǎng)絡(luò)設(shè)計(jì)與運(yùn)行中最基本的運(yùn)算，對(duì)電力網(wǎng)絡(luò)的各種設(shè)計(jì)方案及各種運(yùn)行方式進(jìn)行潮流計(jì)算，可以得到各種電網(wǎng)各節(jié)點(diǎn)的電壓，并求得網(wǎng)絡(luò)的潮流及網(wǎng)絡(luò)中的各元件的電力損耗，進(jìn)而求得電能損耗。本位就是運(yùn)用潮流計(jì)算具體分析，并有MATLAB仿真。

關(guān)鍵詞： 電力系統(tǒng) 潮流計(jì)算 MATLAB

Abstract Electric power system steady flow calculation and analysis of the static safety analysis.This paper, by means of the calculation, flow calculation is the trend of the power network design and operation of the most basic operations of electric power network, various design scheme and the operation ways to tide computation, can get all kinds of each node of the power grid voltage and seek the trend of the network and the network of the components of the power loss, and getting electric power.The standard is to use the power flow calculation and analysis, the specific have MATLAB simulation.Key words: Power system;Flow calculation;MATLAB simulation

目 錄 任務(wù)提出與方案論證....................................................................................................................................2 2 總體設(shè)計(jì)........................................................................................................................................................3 2.1潮流計(jì)算等值電路.............................................................................................................................3 2.2建立電力系統(tǒng)模型.............................................................................................................................3 2.3模型的調(diào)試與運(yùn)行.............................................................................................................................3 3 詳細(xì)設(shè)計(jì)........................................................................................................................................................4 3.1 計(jì)算前提............................................................................................................................................4 3.2手工計(jì)算.............................................................................................................................................7 4設(shè)計(jì)圖及源程序...........................................................................................................................................11 4.1MATLAB仿真.......................................................................................................................................11 4.2潮流計(jì)算源程序...............................................................................................................................11 5 總結(jié).............................................................................................................................................................31 參考文獻(xiàn)..........................................................................................................................................................32 任務(wù)提出與方案論證

潮流計(jì)算是在給定電力系統(tǒng)網(wǎng)絡(luò)結(jié)構(gòu)、參數(shù)和決定系統(tǒng)運(yùn)行狀態(tài)的邊界條件的情況下確定系統(tǒng)穩(wěn)態(tài)運(yùn)行狀態(tài)的一種基本方法,是電力系統(tǒng)規(guī)劃和運(yùn)營(yíng)中不可缺少的一個(gè)重要組成部分?？梢哉f(shuō),它是電力系統(tǒng)分析中最基本、最重要的計(jì)算,是系統(tǒng)安全、經(jīng)濟(jì)分析和實(shí)時(shí)控制與調(diào)度的基礎(chǔ)。常規(guī)潮流計(jì)算的任務(wù)是根據(jù)給定的運(yùn)行條件和網(wǎng)路結(jié)構(gòu)確定整個(gè)系統(tǒng)的運(yùn)行狀態(tài)，如各母線(xiàn)上的電壓（幅值及相角）、網(wǎng)絡(luò)中的功率分布以及功率損耗等。潮流計(jì)算的結(jié)果是電力系統(tǒng)穩(wěn)定計(jì)算和故障分析的基礎(chǔ)。在電力系統(tǒng)運(yùn)行方式和規(guī)劃方案的研究中，都需要進(jìn)行潮流計(jì)算以比較運(yùn)行方式或規(guī)劃供電方案的可行性、可靠性和經(jīng)濟(jì)性。同時(shí)，為了實(shí)時(shí)監(jiān)控電力系統(tǒng)的運(yùn)行狀態(tài)，也需要進(jìn)行大量而快速的潮流計(jì)算。因此，潮流計(jì)算是電力系統(tǒng)中應(yīng)用最廣泛、最基本和最重要的一種電氣運(yùn)算。在系統(tǒng)規(guī)劃設(shè)計(jì)和安排系統(tǒng)的運(yùn)行方式時(shí)，采用離線(xiàn)潮流計(jì)算；在電力系統(tǒng)運(yùn)行狀態(tài)的實(shí)時(shí)監(jiān)控中，則采用在線(xiàn)潮流計(jì)算。是電力系統(tǒng)研究人員長(zhǎng)期研究的一個(gè)課題。它既是對(duì)電力系統(tǒng)規(guī)劃設(shè)計(jì)和運(yùn)行方式的合理性、可靠性及經(jīng)濟(jì)性進(jìn)行定量分析的依據(jù)，又是電力系統(tǒng)靜態(tài)和暫態(tài)穩(wěn)定計(jì)算的基礎(chǔ)。

潮流計(jì)算經(jīng)歷了一個(gè)由手工到應(yīng)用數(shù)字電子計(jì)算機(jī)的發(fā)展過(guò)程，現(xiàn)在的潮流算法都以計(jì)算機(jī)的應(yīng)用為前提用計(jì)算機(jī)進(jìn)行潮流計(jì)算主要步驟在于編制計(jì)算機(jī)程序，這是一項(xiàng)非常復(fù)雜的工作。對(duì)系統(tǒng)進(jìn)行潮流分析，本文利用 MATLAB中的SimpowerSystems工具箱設(shè)計(jì)電力系統(tǒng)，在simulink 環(huán)境下，不僅可以仿真系統(tǒng)的動(dòng)態(tài)過(guò)程，還可以對(duì)系統(tǒng)進(jìn)行穩(wěn)態(tài)潮流分析。

總體設(shè)計(jì)

SimpowerSystems使用Simulink環(huán)境，可以將該系統(tǒng)中的發(fā)電機(jī)、變壓器，線(xiàn)路等模型聯(lián)結(jié)起來(lái)，形成電力系統(tǒng)仿真模擬圖。在加人測(cè)量模塊，并對(duì)各元件的參數(shù)進(jìn)行設(shè)置后，用measurement和sink中的儀器可以觀(guān)察各元件的電壓、電流、功率的大小。

2.1潮流計(jì)算等值電路

10MWYN,d114?63MWVA15MWGGGG120MW10kV??p0?15.7kW??ps?73kW?I%?0.5?0?Vs%?10.5YN,d1116MWVA4?63MW“?xd?0.134?x2?0.161?x?0.06?0?cos?N?0.8510kV??p0?11kW??ps?50kW?I%?0.55?0?Vs%?10.5YN,d112?10MWVA35kV32km25MW110kV80km25MW110kV70km110kVYN,d112?20MWVA20MWGGG4?15MW”?xd?0.136?x2?0.16?x?0.073?0?cos?N?0.8??p0?18.6kW??ps?89kW?I%?0.530MW?0?Vs%?10.510kV??p0?15.7kW??ps?73kW?I%?0.5GG?0V%?10.5s?YN,d112?16MWVA63MWVA??p0?44kW??ps?121kW10kV?I%?0.35?0?Vs%?10.5GG35MWYN,Y,d112?10MVA10kVGGG3?12MW1?50MW“?xd?0.128?x2?0.154?x?0.054?0?cos?N?0.852?25MW”?xd?0.128?x2?0.157?x?0.0591?0?cos?N?0.8?x?0.136?x2?0.161?x?0.075?0?cos?N?0.8“d

2.2建立電力系統(tǒng)模型

在Simulink中按照電力系統(tǒng)原型選擇元件進(jìn)行建模。所建立的模型和建立的方法在詳細(xì)設(shè)計(jì)中詳述。

在電力系統(tǒng)模型的建立工程中主要涉及到的是：元器件的選擇及其參數(shù)的設(shè)置；發(fā)電機(jī)選型；變壓器選擇；線(xiàn)路的選擇；負(fù)荷模型的選擇；母線(xiàn)選擇。

2.3模型的調(diào)試與運(yùn)行

建立系統(tǒng)模型，并設(shè)置好參數(shù)以后，就可以在Simulink環(huán)境下進(jìn)行仿真運(yùn)行。運(yùn)行的具體結(jié)果和分析也在詳細(xì)設(shè)計(jì)中詳述。

30km35kV0km31??p0?44kW??ps?121kW?I%?0.35?0?Vs%?10.5??p0?13.2kW??ps?63kW?I0%?0.55?V%?10.5s(1?2)80MW?Vs(2?3)%?6.55MW??Vs(1?3)%?17.53 詳細(xì)設(shè)計(jì)

3.1 計(jì)算前提

首先是發(fā)電機(jī)的參數(shù)計(jì)算，先對(duì)5個(gè)發(fā)電廠(chǎng)簡(jiǎn)化為5臺(tái)發(fā)電機(jī)來(lái)計(jì)算。發(fā)電機(jī)G1：

P1?4?15?60MWQ1?60?tan(arccos0.8)?45MVar發(fā)電機(jī)G2：

P2?4?63?252MWQ2?252?tan(arccos0.85)?156MVar發(fā)電機(jī)G3：

P3?3?12?36MWQ3?36?tan(arccos0.8)?27MVar發(fā)電機(jī)G4：

P4?1?50?50MWQ4?50?tan(arccos0.85)?31MVar發(fā)電機(jī)G5：

P5?2?25?50MWQ5?50?tan(arccos0.8)?37.5MVar

其次是變電站的參數(shù)計(jì)算，我們還是對(duì)7個(gè)變電站簡(jiǎn)化為7臺(tái)變壓器來(lái)計(jì)算。變壓器T1：

2?ps?VN73?11023RT1??10??103?3.450?232SN(16?10)2Vs%?VN10.5?1102XT1??10??10?79.406?SN16?103?S01??p0?j變壓器T2：（雙并聯(lián)）

I0%?SN?(0.0157?j0.0800)MVA 100RT2XT221?ps?VN189?11023???10???103?1.346?2322SN2(20?10)21Vs%?VN110.5?1102???10???10?31.7625? 32SN220?10?S02?2?(?p0?j變壓器T3：（四并聯(lián)）

I0%?SN)?(0.0372?j0.2000)MVA100 1?ps?VN21121?11023RT3???10???103?0.092?2324SN4(63?10)XT31Vs%?VN2110.5?1102???10???10?5.042? 4SN463?103I0%?SN)?(0.1760?j0.8820)MVA100?S03?4?(?p0?j變壓器T4：（雙并聯(lián)）

1RT1?1.7250?21 XT4?XT1?39.7030?2?S04?2?S01?(0.0314?j0.1600)MVART4?變壓器T5：

RT5?4RT3?0.3680?XT5?4XT3?20.168??S05?1?S03?(0.0440?j0.2205)MVA4163?3523???10?0.386? 322(10?10)

變壓器T6：（兩個(gè)三繞組變壓器并聯(lián)）

RT6?1?RT6?2?RT6?31?[Vs(1?2)%?Vs(1?3)%?Vs(2?3)%]?10.7521Vs2%??[Vs(1?2)%?Vs(2?3)%?Vs(1?3)%]??0.25

21Vs3%??[Vs(1?3)%?Vs(2?3)%?Vs(1?2)%]?6.752Vs1%?21Vs1%?VNXT6?1???10?6.584?2SNXT6?2XT6?321Vs2%?VN???10??0.153?2SN21Vs3%?VN???10?4.134? 2SN?S06?2?(?P06?j變壓器T7：（雙并聯(lián)）

I0%?10)?(0.0264?j0.1100)MVA 100 RT7XT721?ps?VN150?3523???10???103?0.306?2322SN2(10?10)21Vs%?VN110.5?352???10???10?6.431?2SN210?103

?S07?2?(?p0?jI0%?SN)?(0.0220?j0.1100)MVA100再次是傳輸線(xiàn)參數(shù)計(jì)算，5條傳輸線(xiàn)的具體計(jì)算如下。

根據(jù)教材查得r0?0.21?/km x0?0.4?/km b0?2.8?10S/km ?6線(xiàn)路L1：

線(xiàn)路L2：

線(xiàn)路L3：（雙回路）

線(xiàn)路L4：

線(xiàn)路L5：（雙回路）RL1?r0?l1?0.21?40?8.4?XL1?x0?l1?0.4?40?16?B?6L1?b0?l1?2.8?10?40?1.12?10?4S ?Q1L1??2BL1V2N??0.6776MVarRL2?r0?l2?0.21?130?27.3?XL2?x0?l2?0.4?130?52?B?6L2?b0?l2?2.8?10?130?3.64?10?4S ?Q1L2??2BL2V2N??2.2022MVarR?12?rl1L30?3?2?0.21?70?7.35?X11L3?2?x0?l3?2?0.4?70?14? BL3?2?b?40?l3?2?2.8?10?6?70?3.92?10S?Q1L3??2B2L3VN??2.3716MVarRL4?r0?l4?0.21?60?12.6?XL4?x0?l4?0.4?60?24?BL4?b0?l4?2.8?10?6?60?1.68?10?4S ?Q12L4??2BL4VN??1.0164MVar

11RL5??r0?l5??0.21?20?2.1?2211XL5??x0?l5??0.4?20?4? 22BL5?2?b0?l5?2?2.8?10?6?20?1.12?10?4S1?QL5??BL3VN2??0.0686MVar23.2手工計(jì)算

FLR1：

P2102?ST1?2(RT1?jXT1)?(3.450?j74.406)?(0.0285?j0.6562)MVA2VN110Sa?10MW??ST1??S01?j?QL1?(10.0442?j0.1142)MVAP2?Q210.04422?0.11422?SL1?(RL1?jXL1)?(8.4?j16)?(0.070?j0.1334)MVAVN21102?ST2P2?Q2402?452?(RT2?jXT2)?(1.346?j31.7625)?(0.4032?j9.5156)MVAVN21102FLR2Sb?SG1?20??ST2?60?j45?20?0.4032?j9.5156?(39.5968?j35.4844)MVASc?Sb?Sa?25?jQL1??SL1?(4.4826?j35.9144)MVA：

?ST3P2?Q22522?1562?(RT3?jXT3)?(0.092?j5.042)?(0.6679?j36.6024)MVA22VN110P?Q4.4931?34.1048(R?jX)?(27.3?j52)?(2.67?j5.0854)MVAL2L222VN1102222Sc'?(4.4931?j34.1048)MVA?SL2?FLRSd?SG2?Sc'?120??ST3??S03?jQL2??SL2?(132.9792?j149.229)MVA3：

?ST4P2?Q262?272?(RT4?jXT4)?(1.725?j39.703)?(0.1091?j2.5101)MVAVN21102P2?Q2133.59552?149.99562?(RL3?jXL3)?(7.35?j14)?(24.51?j46.682)MVA22VN110'Sd?(133.5955?j149.9956)MVA?SL3'Se?SG3?Sd?30?25??ST4??S04?jQL3??SL3?(89.945?j130.0151)MVAFLR4： ?ST5P2?Q2502?312?(RT5?jXT5)?(0.368?j20.168)?(0.1052?j5.7687)MVA22VN110P2?Q292.74872?133.99372?(RL4?jXL4)?(12.6?j24)?(27.654?j52.674)MVA22VN110Se'?(92.7481?j133.9937)MVA?SL4Sf?SG4?Se'?80??ST5??S05?jQL4??SL4?(34.9449?j107.3469)MVAFLR5: 152?ST7?2?(0.306?j6.431)?(0.0562?j1.1812)MVA35Sh?15??ST7??S07?j?QL5?(15.0782?j0.3422)MVA15.07822?0.34222?SL5??(2.1?j4)?(0.3899?j0.743)MVA352Si?Sh??SL5??S06?j?QL5?5?(20.4945?j1.1266)MVA 152?37.52?ST6?3??(0.386?j4.34)?(0.514?j5.7793)MVA35220.65052?0.54512?ST6?2??(0.386?j0.153)?(0.1345?j0.0533)MVA23526.3362?98.73692?ST6?1??(0.386?j6.584)?(3.2905?j56.1256)MVA352Sg?Sf??ST6?1?SG5??ST6?2??ST6?3?Si?35?(25.5114?j194.12)MVA計(jì)算每一個(gè)FLR的功率分布和電壓分布計(jì)算如下： FLR1：

?VT2?PR?QX40?1.346?45?31.7625??12.8970kVVN115 Vb?115??VT2?102.1030kVPR?QX10.0442?8.4?0.1442?16?VL1???0.8489kVVb102.1030Va?Vb??VL1?101.2541kVFLR2：

功率分布：

SL2?ZZ?Z*T3*L2*Sd?T3(Vb?VN)Z?ZL2****?VN?(0.092?j5.042)?(132.9792?j149.229)?1418.6727.392?j57.042T3?(4.8812?j13.8097)MVAST3?ZZ?Z*L2*L2*Sd?T3(Vb?VN)Z?ZL2?VN?(27.3?j52)?(132.9792?j149.229)?1418.6727.392?j57.042T3?(108.687?j122.62)MVA 電壓分布：

Sc1?SL2??SL2?(4.8812?j13.8097)?(2.67?j5.0854)?(7.5512?j8.7243)MVA7.5512?27.3?8.7243?52??2.424kV102.1030Vd?Vb??VL2?102.103?(?2.424)?104.527kV?VL2?功率分布：

FLR3：

SL3?ZZ?Z*T4*L3*Se?T4(VG3?Vd)Z?ZL3****?VN?(1.725?j39.703)?(89.945?j130.0151)?1037.9279.075?j53.73T4?(59.444?j16.846)MVAST4?ZZ?Z*L3*L3*Se?T4(Vb?VN)Z?ZL3?VN?(7.35?j14)?(89.945?j130.0151)?1037.9279.075?j53.73T4?(31.811?j60.1256)MVA 電壓分布：

Se1?SL3??SL3?(59.444?j19.846)?(24.51?j46.682)?(83.954?j26.836)MVA83.954?7.35?26.836?14?9.404kV105.5643Ve?Vd??VL3?96.16kV?VL3?功率分布：

FLR4：

SL4?ZZ?Z*T5*L4*Sf?T5(VG3?Vd)Z?ZL4****?VN=(0.368?j20.168)?(34.9449?j107.3469)?1037.92712.968?j44.168T5?(20.843?j19.689)MVAST4?ZZ?Z*L4*L4*Sf?T5(VG3?Vd)Z?ZL4?VN=(12.6?j24)?(34.9449?j107.3469)?1037.92712.968?j44.168T5?(1.398?j44.389)MVA 電壓分布：

Se1?SL3??SL3?(59.444?j16.846)?(24.51?j46.682)?(83.954?j63.528)MVA83.954?12.6?63.528?24?24.464kV105.5643Ve?Vd??VL3?81.10kV?VL4?FLR5: 這里我們先將f點(diǎn)和發(fā)電機(jī)G5當(dāng)做電源，經(jīng)過(guò)ZT61和ZT63構(gòu)成兩端供電網(wǎng)絡(luò)以g點(diǎn)作為運(yùn)算負(fù)荷進(jìn)行計(jì)算。ST6?ST4(0.386?j4.134)?(20.2656?j70.9293)?(22.0938?37)?35?(3.900?j25.1175)MVA0.772?j10.718(0.386?j6.584)?(20.2656?j70.9293)?(22.0938?37)?35??(16.5061?j91.7905)MVA0.772?j10.718電壓分布：

ST631?ST63??ST63?(16.6421?j97.5698)MVA16.6421?0.386?97.5698?4.134?10.9186kV37Vg?37??VT63?26.0814V?VT63?20.2656?0.386?70.9293?(?0.153)??0.1162kV

26.0814Vi?Vg??VT62?26.1976?VT62?20.4945?2.1?1.1266?4?1.815kV26.1976Vh?Vi??VL5?24.3826?VL5?

4設(shè)計(jì)圖及源程序

4.1MATLAB仿真

相關(guān)的原始數(shù)據(jù)輸入格式如下：

1、B1是支路參數(shù)矩陣，第一列和第二列是節(jié)點(diǎn)編號(hào)。節(jié)點(diǎn)編號(hào)由小到大編寫(xiě)。

2、對(duì)于含有變壓器的支路，第一列為低壓側(cè)節(jié)點(diǎn)編號(hào)，第二列為高壓側(cè)節(jié)點(diǎn)編號(hào)，將變壓器的串聯(lián)阻抗置于低壓側(cè)處理，第三列為支路的串列阻抗參數(shù)，第四列為支路的對(duì)地導(dǎo)納參數(shù)，第五烈為含變壓器支路的變壓器的變比，第六列為變壓器是否是否含有變壓器的參數(shù)，其中“1”為含有變壓器，“0”為不含有變壓器。

3、B2為節(jié)點(diǎn)參數(shù)矩陣，其中第一列為節(jié)點(diǎn)注入發(fā)電功率參數(shù)；第二列為節(jié)點(diǎn)負(fù)荷功率參數(shù)；第三列為節(jié)點(diǎn)電壓參數(shù)；第六列為節(jié)點(diǎn)類(lèi)型參數(shù)，其中“1”為平衡節(jié)點(diǎn)，“2”為PQ節(jié)點(diǎn)，“3”為PV節(jié)點(diǎn)參數(shù)。

4、X為節(jié)點(diǎn)號(hào)和對(duì)地參數(shù)矩陣。其中第一列為節(jié)點(diǎn)編號(hào)，第二列為節(jié)點(diǎn)對(duì)地參數(shù)。

4.2潮流計(jì)算源程序

%本程序的功能是用牛頓——拉夫遜法進(jìn)行11節(jié)點(diǎn)潮流計(jì)算 clear;n=11;%input('請(qǐng)輸入節(jié)點(diǎn)數(shù)：n=');nl=11;%input('請(qǐng)輸入支路數(shù)：nl=');isb=1;%input('請(qǐng)輸入平衡母線(xiàn)節(jié)點(diǎn)號(hào)：isb=');pr=0.00001;%input('請(qǐng)輸入誤差精度：pr=');B1=[1

0.03512+0.08306i

0.13455i

0;

0.0068+0.18375i

0

1.02381

1;

0.05620+0.13289i

0.05382i

0;

0.00811+0.24549i

0

1.02381

1;

0.05620+0.13289i

0.05382i

0;

0.04215+0.09967i

0.04037i

0;

0.0068+0.18375i

0

1.02381

1;

0.02810+0.06645i

0.10764i

0;

0.05620+0.13289i

0.05382i

0;0.00811+0.24549i

0

1;

0.03512+0.08306i

0.13455i

0] B2=[0

0

1.1

1.1

0

1;

0

0

0

0

2;

0

0.343+0.21256i

0

0

2;

0

0

0

0

2;

0

0.204+0.12638i

0

0

2;

0

0

0

0

2;

0

0.306+0.18962i

0

0

2;

0

0

0

0

2;

0.5

0

1.1

1.1

0

3;

0

0.343+0.21256i

0

0

2;

0

0

0

0

2];% B1矩陣：

1、支路首端號(hào)；

2、末端號(hào)；

3、支路阻抗；

4、支路對(duì)地電納 %

5、支路的變比；

6、支路首端處于K側(cè)為1，1側(cè)為0 % 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)所接的無(wú)功補(bǔ)償設(shè)備的容量 %

6、節(jié)點(diǎn)分類(lèi)標(biāo)號(hào)：1為平衡節(jié)點(diǎn)（應(yīng)為1號(hào)節(jié)點(diǎn)）；2為PQ節(jié)點(diǎn)； %

3為PV節(jié)點(diǎn)；

%input('請(qǐng)輸入各節(jié)點(diǎn)參數(shù)形成的矩陣： B2=');Y=zeros(n);e=zeros(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));%非對(duì)角元

Y(q,p)=Y(p,q);

%非對(duì)角元

Y(q,q)=Y(q,q)+1./(B1(i,3)*B1(i,5)^2)+B1(i,4)./2;%對(duì)角元K側(cè)

Y(p,p)=Y(p,p)+1./B1(i,3)+B1(i,4)./2;

%對(duì)角元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)電壓給定模值

end for 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)無(wú)功補(bǔ)償量

end %===================== P=real(S);Q=imag(S);

%分解出各節(jié)點(diǎn)注入的有功和無(wú)功功率 ICT1=0;IT2=1;N0=2*n;N=N0+1;a=0;%迭代次數(shù)ICT1、a；不滿(mǎn)足收斂要求的節(jié)點(diǎn)數(shù)IT2 while 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)*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)有功和無(wú)功功率P',Q'的計(jì)算值

V2=e(i)^2+f(i)^2;

%電壓模平方

%========= 以下針對(duì)非PV節(jié)點(diǎn)來(lái)求取功率差及Jacobi矩陣元素 =========

if B2(i,6)~=3

%非PV節(jié)點(diǎn)

DP=P(i)-P1;

%節(jié)點(diǎn)有功功率差

DQ=Q(i)-Q1;

%節(jié)點(diǎn)無(wú)功功率差

%=============== 以上為除平衡節(jié)點(diǎn)外其它節(jié)點(diǎn)的功率計(jì)算 ================= %================= 求取Jacobi矩陣 ===================

for j1=1:n

if j1~=isb&j1~=i

%非平衡節(jié)點(diǎn)&非對(duì)角元

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)無(wú)功功率差

J(m,q)=X1;J(m,N)=DP;q=q+1;

% X1=dP/de J(m,N)=DP節(jié)點(diǎn)有功功率差

J(p,q)=X4;J(m,q)=X2;

% X4=dQ/df X2=dp/df

elseif j1==i&j1~=isb %非平衡節(jié)點(diǎn)&對(duì)角元

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=q+1;J(p,q)=X4;J(m,N)=DP;%擴(kuò)展列△P

J(m,q)=X2;

end

end

else

%=============== 下面是針對(duì)PV節(jié)點(diǎn)來(lái)求取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)&非對(duì)角元

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)電壓誤差

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)&對(duì)角元

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;

% PV節(jié)點(diǎn)有功誤差

J(m,q)=X2;

end

end

end

end

end %========= 以上為求雅可比矩陣的各個(gè)元素及擴(kuò)展列的功率差或電壓差 =====================

for k=3:N0

% N0=2*n（從第三行開(kāi)始，第一、二行是平衡節(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)+eps);% 用K行K列對(duì)角元素去除K行K列后的非對(duì)角元素進(jìn)行規(guī)格化

end

J(k,k)=1;

% 對(duì)角元規(guī)格化K行K列對(duì)角元素賦1

%==================== 回代運(yùn)算

=======================================

if k~=3

% 不是第三行

k > 3

k4=k-1;

for k3=3:k4

% 用k3行從第三行開(kāi)始到當(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行K14 列元素

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行后各行下三角元素

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)前行第3列元素已消為0

end

end

end %====上面是用線(xiàn)性變換方式高斯消去法將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

%電壓偏差量是否滿(mǎn)足要求

IT2=IT2+1;%不滿(mǎn)足要求的節(jié)點(diǎn)數(shù)加1

end

end

ICT2(a)=IT2;

%不滿(mǎn)足要求的節(jié)點(diǎn)數(shù)

ICT1=ICT1+1;

%迭代次數(shù) end %用高斯消去法解”w=-J*V“

disp('迭代次數(shù)：');disp(ICT1);disp('沒(méi)有達(dá)到精度要求的個(gè)數(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)號(hào)從小到大排列)：');disp(E);

%顯示各節(jié)點(diǎn)的實(shí)際電壓標(biāo)幺值E用復(fù)數(shù)表示 disp('----------------------');disp('各節(jié)點(diǎn)的電壓大小V為(節(jié)點(diǎn)號(hào)從小到大排列)：');disp(V);

%顯示各節(jié)點(diǎn)的電壓大小V的模值 disp('----------------------');disp('各節(jié)點(diǎn)的電壓相角sida為(節(jié)點(diǎn)號(hào)從小到大排列)：');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 注入電流的共軛值 end disp('各節(jié)點(diǎn)的功率S為(節(jié)點(diǎn)號(hào)從小到大排列)：');disp(S);

%顯示各節(jié)點(diǎn)的注入功率

disp('----------------------');disp('各條支路的首端功率Si為(順序同您輸入B1時(shí)一致)：');for i=1:nl

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)./2)+(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)./2)+(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=['S(',num2str(p),',',num2str(q),')=',num2str(SSi(p,q))];

disp(ZF);

disp('----------------------');end disp('各條支路的末端功率Sj為(順序同您輸入B1時(shí)一致)：');

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,4)./2)+(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,4)./2)+(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('----------------------');end disp('各條支路的功率損耗DS為(順序同您輸入B1時(shí)一致)：');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(',num2str(p),',',num2str(q),')=',num2str(DDS(i))];

disp(ZF);

disp('----------------------');end

%本程序的功能是用牛頓——拉夫遜法進(jìn)行10節(jié)點(diǎn)潮流計(jì)算 %本程序的功能是用牛頓——拉夫遜法進(jìn)行潮流計(jì)算 clear;n=10;%input('請(qǐng)輸入節(jié)點(diǎn)數(shù)：n=');nl=10;%input('請(qǐng)輸入支路數(shù)：nl=');isb=1;%input('請(qǐng)輸入平衡母線(xiàn)節(jié)點(diǎn)號(hào)：isb=');pr=0.00001;%input('請(qǐng)輸入誤差精度：pr=');B1=[1

0.03512+0.08306i

0.13455i

0;

0.0068+0.18375i

0

1.02381

1;

0.05620+0.13289i

0.05382i

0;

0.00811+0.24549i

0

1.02381

1;

0.05620+0.13289i

0.05382i

0;

0.04215+0.09967i

0.04037i

0;

0.0068+0.18375i

0

1.02381

1;

0.02810+0.06645i

0.10764i

0;0.00811+0.24549i

0

1;

0.03512+0.08306i

0.13455i

0] B2=[0

0

1.1

1.1

0

1;

0

0

0

0

2;

0

0.343+0.21256i

0

0

2;

0

0

0

0

2;

0

0.204+0.12638i

0

0

2;

0

0

0

0

2;

0

0.306+0.18962i

0

0

2;

0

0

0

0

2;

0.5

0

1.1

1.1

0

3;

0

0.343+0.21256i

0

0

2];% B1矩陣：

1、支路首端號(hào)；

2、末端號(hào)；

3、支路阻抗；

4、支路對(duì)地電納 %

5、支路的變比；

6、支路首端處于K側(cè)為1，1側(cè)為0 % 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)所接的無(wú)功補(bǔ)償設(shè)備的容量 %

6、節(jié)點(diǎn)分類(lèi)標(biāo)號(hào)：1為平衡節(jié)點(diǎn)（應(yīng)為1號(hào)節(jié)點(diǎn)）；2為PQ節(jié)點(diǎn)； %

3為PV節(jié)點(diǎn)；

%input('請(qǐng)輸入各節(jié)點(diǎn)參數(shù)形成的矩陣： B2=');Y=zeros(n);e=zeros(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));%非對(duì)角元

Y(q,p)=Y(p,q);

%非對(duì)角元

Y(q,q)=Y(q,q)+1./(B1(i,3)*B1(i,5)^2)+B1(i,4)./2;%對(duì)角元K側(cè)

Y(p,p)=Y(p,p)+1./B1(i,3)+B1(i,4)./2;

%對(duì)角元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)電壓給定模值

end for 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)無(wú)功補(bǔ)償量

end %===================== P=real(S);Q=imag(S);

%分解出各節(jié)點(diǎn)注入的有功和無(wú)功功率

ICT1=0;IT2=1;N0=2*n;N=N0+1;a=0;%迭代次數(shù)ICT1、a；不滿(mǎn)足收斂要求的節(jié)點(diǎn)數(shù)IT2 while 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)*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)有功和無(wú)功功率P',Q'的計(jì)算值

V2=e(i)^2+f(i)^2;

%電壓模平方

%========= 以下針對(duì)非PV節(jié)點(diǎn)來(lái)求取功率差及Jacobi矩陣元素 =========

if B2(i,6)~=3

%非PV節(jié)點(diǎn)

DP=P(i)-P1;

%節(jié)點(diǎn)有功功率差

DQ=Q(i)-Q1;

%節(jié)點(diǎn)無(wú)功功率差

%=============== 以上為除平衡節(jié)點(diǎn)外其它節(jié)點(diǎn)的功率計(jì)算 ================= %================= 求取Jacobi矩陣 ===================

for j1=1:n

if j1~=isb&j1~=i

%非平衡節(jié)點(diǎn)&非對(duì)角元

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)無(wú)功功率差

J(m,q)=X1;J(m,N)=DP;q=q+1;

% X1=dP/de J(m,N)=DP節(jié)點(diǎn)有功功率差

J(p,q)=X4;J(m,q)=X2;

% X4=dQ/df X2=dp/df

elseif j1==i&j1~=isb %非平衡節(jié)點(diǎn)&對(duì)角元

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=q+1;J(p,q)=X4;J(m,N)=DP;%擴(kuò)展列△P

J(m,q)=X2;

end

end

else

%=============== 下面是針對(duì)PV節(jié)點(diǎn)來(lái)求取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)&非對(duì)角元

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)電壓誤差

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)&對(duì)角元

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;

% PV節(jié)點(diǎn)有功誤差

J(m,q)=X2;

end

end

end

end

end %========= 以上為求雅可比矩陣的各個(gè)元素及擴(kuò)展列的功率差或電壓差 =====================

for k=3:N0

% N0=2*n（從第三行開(kāi)始，第一、二行是平衡節(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列對(duì)角元素去除K行K列后的非對(duì)角元素進(jìn)行規(guī)格化

end

J(k,k)=1;

% 對(duì)角元規(guī)格化K行K列對(duì)角元素賦1

%==================== 回代運(yùn)算

=======================================

if k~=3

% 不是第三行

k > 3

k4=k-1;

for k3=3:k4

% 用k3行從第三行開(kāi)始到當(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行后各行下三角元素

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)前行第3列元素已消為0

end

end

end %====上面是用線(xiàn)性變換方式高斯消去法將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

%電壓偏差量是否滿(mǎn)足要求

IT2=IT2+1;%不滿(mǎn)足要求的節(jié)點(diǎn)數(shù)加1

end

end

ICT2(a)=IT2;

%不滿(mǎn)足要求的節(jié)點(diǎn)數(shù)

ICT1=ICT1+1;

%迭代次數(shù) end %用高斯消去法解”w=-J*V“ disp('迭代次數(shù)：');disp(ICT1);disp('沒(méi)有達(dá)到精度要求的個(gè)數(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)號(hào)從小到大排列)：');disp(E);

%顯示各節(jié)點(diǎn)的實(shí)際電壓標(biāo)幺值E用復(fù)數(shù)表示 disp('----------------------');disp('各節(jié)點(diǎn)的電壓大小V為(節(jié)點(diǎn)號(hào)從小到大排列)：');disp(V);

%顯示各節(jié)點(diǎn)的電壓大小V的模值 disp('----------------------');disp('各節(jié)點(diǎn)的電壓相角sida為(節(jié)點(diǎn)號(hào)從小到大排列)：');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 注入電流的共軛值 end disp('各節(jié)點(diǎn)的功率S為(節(jié)點(diǎn)號(hào)從小到大排列)：');disp(S);

%顯示各節(jié)點(diǎn)的注入功率

disp('----------------------');disp('各條支路的首端功率Si為(順序同您輸入B1時(shí)一致)：');for i=1:nl

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)./2)+(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)./2)+(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=['S(',num2str(p),',',num2str(q),')=',num2str(SSi(p,q))];

disp(ZF);

disp('----------------------');end disp('各條支路的末端功率Sj為(順序同您輸入B1時(shí)一致)：');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,4)./2)+(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,4)./2)+(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('----------------------');end disp('各條支路的功率損耗DS為(順序同您輸入B1時(shí)一致)：');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(',num2str(p),',',num2str(q),')=',num2str(DDS(i))];

disp(ZF);

disp('----------------------');end

%本程序的功能是用牛頓——拉夫遜法進(jìn)行12節(jié)點(diǎn)潮流計(jì)算 %本程序的功能是用牛頓——拉夫遜法進(jìn)行潮流計(jì)算 clear;n=12;%input('請(qǐng)輸入節(jié)點(diǎn)數(shù)：n=');nl=12;%input('請(qǐng)輸入支路數(shù)：nl=');isb=1;%input('請(qǐng)輸入平衡母線(xiàn)節(jié)點(diǎn)號(hào)：isb=');pr=0.00001;%input('請(qǐng)輸入誤差精度：pr=');B1=[1

0.03512+0.08306i

0.13455i

0;

0.0068+0.18375i

0

1.02381

1;

0.05620+0.13289i

0.05382i

0;

0.00811+0.24549i

0

1.02381

1;

0.05620+0.13289i

0.05382i

0;

0.04215+0.09967i

0.04037i

0;

0.0068+0.18375i

0

1.02381

1;

0.02810+0.06645i

0.10764i

0;

0.05620+0.13289i

0.05382i

0;0.00811+0.24549i

0

1;

0.03512+0.08306i

0.13455i

0;

0.03512+0.08306i

0.13455i

0] B2=[0

0

1.1

1.1

0

1;

0

0

0

0

2;

0

0.343+0.21256i

0

0

2;

0

0

0

0

2;

0

0.204+0.12638i

0

0

2;

0

0

0

0

2;

0

0.306+0.18962i

0

0

2;

0

0

0

0

2;

0.5

0

1.1

1.1

0

3;

0

0.343+0.21256i

0

0

2;

0

0

0

0

2;

0

0

0

0

2];% B1矩陣：

1、支路首端號(hào)；

2、末端號(hào)；

3、支路阻抗；

4、支路對(duì)地電納 %

5、支路的變比；

6、支路首端處于K側(cè)為1，1側(cè)為0 % 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)所接的無(wú)功補(bǔ)償設(shè)備的容量 %

6、節(jié)點(diǎn)分類(lèi)標(biāo)號(hào)：1為平衡節(jié)點(diǎn)（應(yīng)為1號(hào)節(jié)點(diǎn)）；2為PQ節(jié)點(diǎn)； %

3為PV節(jié)點(diǎn)；

%input('請(qǐng)輸入各節(jié)點(diǎn)參數(shù)形成的矩陣： B2=');Y=zeros(n);e=zeros(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));%非對(duì)角元

Y(q,p)=Y(p,q);

%非對(duì)角元

Y(q,q)=Y(q,q)+1./(B1(i,3)*B1(i,5)^2)+B1(i,4)./2;%對(duì)角元K側(cè)

Y(p,p)=Y(p,p)+1./B1(i,3)+B1(i,4)./2;

%對(duì)角元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)電壓給定模值

end for 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)無(wú)功補(bǔ)償量

end %===================== P=real(S);Q=imag(S);

%分解出各節(jié)點(diǎn)注入的有功和無(wú)功功率 ICT1=0;IT2=1;N0=2*n;N=N0+1;a=0;%迭代次數(shù)ICT1、a；不滿(mǎn)足收斂要求的節(jié)點(diǎn)數(shù)IT2 while 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)*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)有功和無(wú)功功率P',Q'的計(jì)算值

V2=e(i)^2+f(i)^2;

%電壓模平方

%========= 以下針對(duì)非PV節(jié)點(diǎn)來(lái)求取功率差及Jacobi矩陣元素 =========

if B2(i,6)~=3

%非PV節(jié)點(diǎn)

DP=P(i)-P1;

%節(jié)點(diǎn)有功功率差

DQ=Q(i)-Q1;

%節(jié)點(diǎn)無(wú)功功率差

%=============== 以上為除平衡節(jié)點(diǎn)外其它節(jié)點(diǎn)的功率計(jì)算 ================= %================= 求取Jacobi矩陣 ===================

for j1=1:n

if j1~=isb&j1~=i

%非平衡節(jié)點(diǎn)&非對(duì)角元

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)無(wú)功功率差

J(m,q)=X1;J(m,N)=DP;q=q+1;

% X1=dP/de J(m,N)=DP節(jié)點(diǎn)有功功率差

J(p,q)=X4;J(m,q)=X2;

% X4=dQ/df X2=dp/df

elseif j1==i&j1~=isb %非平衡節(jié)點(diǎn)&對(duì)角元

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=q+1;J(p,q)=X4;J(m,N)=DP;%擴(kuò)展列△P

J(m,q)=X2;

end

end

else

%=============== 下面是針對(duì)PV節(jié)點(diǎn)來(lái)求取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)&非對(duì)角元

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)電壓誤差

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)&對(duì)角元

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;

% PV節(jié)點(diǎn)有功誤差

J(m,q)=X2;

end

end

end

end

end %========= 以上為求雅可比矩陣的各個(gè)元素及擴(kuò)展列的功率差或電壓差 =====================

for k=3:N0

% N0=2*n（從第三行開(kāi)始，第一、二行是平衡節(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)+eps);% 用K行K列對(duì)角元素去除K行K列后的非對(duì)角元素進(jìn)行規(guī)格化

end

J(k,k)=1;

% 對(duì)角元規(guī)格化K行K列對(duì)角元素賦1

%==================== 回代運(yùn)算

=======================================

if k~=3

% 不是第三行

k > 3

k4=k-1;

for k3=3:k4

% 用k3行從第三行開(kāi)始到當(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行后各行下三角元素

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)前行第3列元素已消為0

end

end

end %====上面是用線(xiàn)性變換方式高斯消去法將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

%電壓偏差量是否滿(mǎn)足要求

IT2=IT2+1;%不滿(mǎn)足要求的節(jié)點(diǎn)數(shù)加1

end

end

ICT2(a)=IT2;

%不滿(mǎn)足要求的節(jié)點(diǎn)數(shù)

ICT1=ICT1+1;

%迭代次數(shù) end %用高斯消去法解”w=-J*V" disp('迭代次數(shù)：');disp(ICT1);disp('沒(méi)有達(dá)到精度要求的個(gè)數(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)號(hào)從小到大排列)：');disp(E);

%顯示各節(jié)點(diǎn)的實(shí)際電壓標(biāo)幺值E用復(fù)數(shù)表示 disp('----------------------');disp('各節(jié)點(diǎn)的電壓大小V為(節(jié)點(diǎn)號(hào)從小到大排列)：');disp(V);

%顯示各節(jié)點(diǎn)的電壓大小V的模值 disp('----------------------');disp('各節(jié)點(diǎn)的電壓相角sida為(節(jié)點(diǎn)號(hào)從小到大排列)：');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 注入電流的共軛值 end disp('各節(jié)點(diǎn)的功率S為(節(jié)點(diǎn)號(hào)從小到大排列)：');disp(S);

%顯示各節(jié)點(diǎn)的注入功率

disp('----------------------');disp('各條支路的首端功率Si為(順序同您輸入B1時(shí)一致)：');for i=1:nl

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)./2)+(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)./2)+(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=['S(',num2str(p),',',num2str(q),')=',num2str(SSi(p,q))];

disp(ZF);

disp('----------------------');end disp('各條支路的末端功率Sj為(順序同您輸入B1時(shí)一致)：');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,4)./2)+(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,4)./2)+(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('----------------------');end disp('各條支路的功率損耗DS為(順序同您輸入B1時(shí)一致)：');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(',num2str(p),',',num2str(q),')=',num2str(DDS(i))];

disp(ZF);

disp('----------------------');end

如果源程序的運(yùn)行結(jié)果需要作圖可用下面的程序

figure(1);subplot(1,2,1);plot(V);xlabel('節(jié)點(diǎn)號(hào)');ylabel('電壓標(biāo)幺值');grid on;subplot(1,2,2);plot(sida);xlabel('節(jié)點(diǎn)號(hào)');ylabel('電壓角度');grid on;figure(2);subplot(2,2,1);P=real(S);Q=imag(S);bar(P);xlabel('節(jié)點(diǎn)號(hào)');ylabel('節(jié)點(diǎn)注入有功');grid on;subplot(2,2,2);bar(Q);xlabel('節(jié)點(diǎn)號(hào)');ylabel('節(jié)點(diǎn)注入無(wú)功');grid on;subplot(2,2,3);P1=real(Siz);Q1=imag(Siz);bar(P1);xlabel('支路號(hào)');ylabel('支路首端注入有功');grid on;subplot(2,2,4);bar(Q1);xlabel('支路號(hào)');ylabel('支路首端注入無(wú)功');grid on;

總結(jié)

通過(guò)本次課程設(shè)計(jì)讓我有復(fù)習(xí)了一次潮流計(jì)算的相關(guān)知識(shí)，跟家清晰了什么事潮流計(jì)算以及潮流計(jì)算的在電力系統(tǒng)的重要性。電力系統(tǒng)的穩(wěn)定運(yùn)行狀況即是正常運(yùn)行狀況，是指電力系統(tǒng)在穩(wěn)定運(yùn)行條件下電壓、功率的分布，也稱(chēng)為潮流分布。電力系統(tǒng)分析的潮流計(jì)算是電力系統(tǒng)分析的一個(gè)重要的部分。通過(guò)對(duì)電力系統(tǒng)潮流分布的分析和計(jì)算，可進(jìn)一步對(duì)系統(tǒng)運(yùn)行的安全性，經(jīng)濟(jì)性進(jìn)行分析、評(píng)估，提出改進(jìn)措施。同時(shí)潮流分布也是電力系統(tǒng)規(guī)劃設(shè)計(jì)的一項(xiàng)基礎(chǔ)工作。

整個(gè)計(jì)算過(guò)程的模型建立并不是十分復(fù)雜，但計(jì)算過(guò)程十分繁瑣、計(jì)算量相當(dāng)?shù)拇?，而且由于枝?jié)太多很容易算錯(cuò)。不過(guò)在計(jì)算潮流計(jì)算的過(guò)程中卻對(duì)以往學(xué)過(guò)的電力系統(tǒng)分析的相關(guān)知識(shí)進(jìn)行了一次較為深入的復(fù)習(xí)。而且整個(gè)計(jì)算對(duì)計(jì)算量的要求很大，鍛煉了我們的計(jì)算能力。而且對(duì)細(xì)節(jié)的把握也得到了鍛煉，做題的精細(xì)程度得到了提高。

參考文獻(xiàn)

[1]何仰贊, 溫增銀《電力系統(tǒng)分析》（第三版）［M].華中科技大學(xué)，2002 [2]http://baike.baidu.com/view/627420.[3]王守相，劉玉田 電力系統(tǒng)潮流計(jì)算研究現(xiàn)狀--《山東電力技術(shù)》1996年05期

[4]何仰贊，溫增銀.電力系統(tǒng)分析（上冊(cè)）第三版［M].湖北：華中科技大學(xué)出版社，2002 [5] 劉同娟.MATLAB在電路分析中的應(yīng)用.電氣電子教學(xué)學(xué)報(bào).2002 [6] 西安交通大學(xué)等.電力系統(tǒng)計(jì)算[M].北京：水利電力出版社，1993.12 [7] 李光琦.電力系統(tǒng)暫態(tài)分析[M].北京： 水利電力出版社，2002.5 [8]何仰贊，溫增銀.電力系統(tǒng)分析（下冊(cè)）第三版［M].湖北：華中科技大學(xué)出版社，2002 [9]韋化，李濱，杭乃善，等．大規(guī)模水一火電力系統(tǒng)最優(yōu)潮流的現(xiàn)代內(nèi)點(diǎn)算法實(shí)現(xiàn)[J]．中國(guó)電機(jī)工程學(xué)報(bào)，2003．23(6)：13一l8．

[10]Chen Luo—nan，Suzuki Hideki，Katou Kazuo．Mean fieldtheory for optimal power flow[J]．IEEE Transactions OilPower Systems，1997，12(4)：1481·1486

第三篇：基于MATLAB的電力系統(tǒng)潮流計(jì)算

基于MATLAB的電力系統(tǒng)潮流計(jì)算

%簡(jiǎn)單潮流計(jì)算的小程序，相關(guān)的原始數(shù)據(jù)數(shù)據(jù)數(shù)據(jù)輸入格式如下：

%B1是支路參數(shù)矩陣，第一列和第二列是節(jié)點(diǎn)編號(hào)。節(jié)點(diǎn)編號(hào)由小到大編寫(xiě) %對(duì)于含有變壓器的支路，第一列為低壓側(cè)節(jié)點(diǎn)編號(hào)，第二列為高壓側(cè)節(jié)點(diǎn) %編號(hào)，將變壓器的串聯(lián)阻抗置于低壓側(cè)處理。%第三列為支路的串列阻抗參數(shù)。%第四列為支路的對(duì)地導(dǎo)納參數(shù)。

%第五烈為含變壓器支路的變壓器的變比

%第六列為變壓器是否是否含有變壓器的參數(shù)，其中“1”為含有變壓器，%“0”為不含有變壓器。

%B2為節(jié)點(diǎn)參數(shù)矩陣，其中第一列為節(jié)點(diǎn)注入發(fā)電功率參數(shù)；第二列為節(jié)點(diǎn) %負(fù)荷功率參數(shù)；第三列為節(jié)點(diǎn)電壓參數(shù)；第六列為節(jié)點(diǎn)類(lèi)型參數(shù)，其中 %“1”為平衡節(jié)點(diǎn)，“2”為PQ節(jié)點(diǎn)，“3”為PV節(jié)點(diǎn)參數(shù)。

%X為節(jié)點(diǎn)號(hào)和對(duì)地參數(shù)矩陣。其中第一列為節(jié)點(diǎn)編號(hào)，第二列為節(jié)點(diǎn)對(duì)地 %參數(shù)。

n=input('請(qǐng)輸入節(jié)點(diǎn)數(shù):n=');n1=input('請(qǐng)輸入支路數(shù):n1=');isb=input('請(qǐng)輸入平衡節(jié)點(diǎn)號(hào):isb=');pr=input('請(qǐng)輸入誤差精度:pr=');B1=input('請(qǐng)輸入支路參數(shù):B1=');B2=input('請(qǐng)輸入節(jié)點(diǎn)參數(shù):B2=');X=input('節(jié)點(diǎn)號(hào)和對(duì)地參數(shù):X=');Y=zeros(n);Times=1;%置迭代次數(shù)為初始值 %創(chuàng)建節(jié)點(diǎn)導(dǎo)納矩陣 for i=1:n1 if B1(i,6)==0 %不含變壓器的支路 p=B1(i,1);q=B1(i,2);Y(p,q)=Y(p,q)-1/B1(i,3);Y(q,p)=Y(p,q);Y(p,p)=Y(p,p)+1/B1(i,3)+0.5*B1(i,4);Y(q,q)=Y(q,q)+1/B1(i,3)+0.5*B1(i,4);else %含有變壓器的支路 p=B1(i,1);q=B1(i,2);Y(p,q)=Y(p,q)-1/(B1(i,3)*B1(i,5));Y(q,p)=Y(p,q);Y(p,p)=Y(p,p)+1/B1(i,3);Y(q,q)=Y(q,q)+1/(B1(i,5)^2*B1(i,3));end end Y OrgS=zeros(2*n-2,1);DetaS=zeros(2*n-2,1);%將OrgS、DetaS初始化

%創(chuàng)建OrgS，用于存儲(chǔ)初始功率參數(shù) h=0;j=0;for i=1:n %對(duì)PQ節(jié)點(diǎn)的處理 if i~=isb&B2(i,6)==2 h=h+1;for j=1:n OrgS(2*h-1,1)=OrgS(2*h-1,1)+real(B2(i,3))*(real(Y(i,j))*real(B2(j,3))-imag(Y(i,j))*imag(B2(j,3)))+imag(B2(i,3))*(real(Y(i,j))*imag(B2(j,3))+imag(Y(i,j))*real(B2(j,3)));OrgS(2*h,1)=OrgS(2*h,1)+imag(B2(i,3))*(real(Y(i,j))*real(B2(j,3))-imag(Y(i,j))*imag(B2(j,3)))-real(B2(i,3))*(real(Y(i,j))*imag(B2(j,3))+imag(Y(i,j))*real(B2(j,3)));end end end for i=1:n %對(duì)PV節(jié)點(diǎn)的處理，注意這時(shí)不可再將h初始化為0 if i~=isb&B2(i,6)==3 h=h+1;for j=1:n OrgS(2*h-1,1)=OrgS(2*h-1,1)+real(B2(i,3))*(real(Y(i,j))*real(B2(j,3))-imag(Y(i,j))*imag(B2(j,3)))+imag(B2(i,3))*(real(Y(i,j))*imag(B2(j,3))+imag(Y(i,j))*real(B2(j,3)));OrgS(2*h,1)=OrgS(2*h,1)+imag(B2(i,3))*(real(Y(i,j))*real(B2(j,3))-imag(Y(i,j))*imag(B2(j,3)))-real(B2(i,3))*(real(Y(i,j))*imag(B2(j,3))+imag(Y(i,j))*real(B2(j,3)));end end end OrgS %創(chuàng)建PVU 用于存儲(chǔ)PV節(jié)點(diǎn)的初始電壓 PVU=zeros(n-h-1,1);t=0;for i=1:n if B2(i,6)==3 t=t+1;PVU(t,1)=B2(i,3);end end PVU %創(chuàng)建DetaS，用于存儲(chǔ)有功功率、無(wú)功功率和電壓幅值的不平衡量 h=0;for i=1:n %對(duì)PQ節(jié)點(diǎn)的處理 if i~=isb&B2(i,6)==2 h=h+1;DetaS(2*h-1,1)=real(B2(i,2))-OrgS(2*h-1,1);DetaS(2*h,1)=imag(B2(i,2))-OrgS(2*h,1);end end t=0;for i=1:n %對(duì)PV節(jié)點(diǎn)的處理，注意這時(shí)不可再將h初始化為0 if i~=isb&B2(i,6)==3 h=h+1;t=t+1;DetaS(2*h-1,1)=real(B2(i,2))-OrgS(2*h-1,1);DetaS(2*h,1)=real(PVU(t,1))^2+imag(PVU(t,1))^2-real(B2(i,3))^2-imag(B2(i,3))^2;end end DetaS %創(chuàng)建I，用于存儲(chǔ)節(jié)點(diǎn)電流參數(shù) i=zeros(n-1,1);h=0;for i=1:n if i~=isb h=h+1;I(h,1)=(OrgS(2*h-1,1)-OrgS(2*h,1)*sqrt(-1))/conj(B2(i,3));end end I %創(chuàng)建Jacbi(雅可比矩陣)Jacbi=zeros(2*n-2);h=0;k=0;for i=1:n %對(duì)PQ節(jié)點(diǎn)的處理 if B2(i,6)==2 h=h+1;for j=1:n if j~=isb k=k+1;if i==j %對(duì)角元素的處理

Jacbi(2*h-1,2*k-1)=-imag(Y(i,j))*real(B2(i,3))+real(Y(i,j))*imag(B2(i,3))+imag(I(h,1));Jacbi(2*h-1,2*k)=real(Y(i,j))*real(B2(i,3))+imag(Y(i,j))*imag(B2(i,3))+real(I(h,1));Jacbi(2*h,2*k-1)=-Jacbi(2*h-1,2*k)+2*real(I(h,1));Jacbi(2*h,2*k)=Jacbi(2*h-1,2*k-1)-2*imag(I(h,1));else %非對(duì)角元素的處理

Jacbi(2*h-1,2*k-1)=-imag(Y(i,j))*real(B2(i,3))+real(Y(i,j))*imag(B2(i,3));Jacbi(2*h-1,2*k)=real(Y(i,j))*real(B2(i,3))+imag(Y(i,j))*imag(B2(i,3));Jacbi(2*h,2*k-1)=-Jacbi(2*h-1,2*k);Jacbi(2*h,2*k)=Jacbi(2*h-1,2*k-1);end if k==(n-1)%將用于內(nèi)循環(huán)的指針置于初始值，以確保雅可比矩陣換行

k=0;end end end end end k=0;for i=1:n %對(duì)PV節(jié)點(diǎn)的處理 if B2(i,6)==3 h=h+1;for j=1:n if j~=isb k=k+1;if i==j %對(duì)角元素的處理

Jacbi(2*h-1,2*k-1)=-imag(Y(i,j))*real(B2(i,3))+real(Y(i,j))*imag(B2(i,3))+imag(I(h,1));Jacbi(2*h-1,2*k)=real(Y(i,j))*real(B2(i,3))+imag(Y(i,j))*imag(B2(i,3))+real(I(h,1));Jacbi(2*h,2*k-1)=2*imag(B2(i,3));Jacbi(2*h,2*k)=2*real(B2(i,3));else %非對(duì)角元素的處理

Jacbi(2*h-1,2*k-1)=-imag(Y(i,j))*real(B2(i,3))+real(Y(i,j))*imag(B2(i,3));Jacbi(2*h-1,2*k)=real(Y(i,j))*real(B2(i,3))+imag(Y(i,j))*imag(B2(i,3));Jacbi(2*h,2*k-1)=0;Jacbi(2*h,2*k)=0;end if k==(n-1)%將用于內(nèi)循環(huán)的指針置于初始值，以確保雅可比矩陣換行

k=0;end end end end end Jacbi %求解修正方程，獲取節(jié)點(diǎn)電壓的不平衡量 DetaU=zeros(2*n-2,1);DetaU=inv(Jacbi)*DetaS;DetaU %修正節(jié)點(diǎn)電壓 j=0;for i=1:n %對(duì)PQ節(jié)點(diǎn)處理 if B2(i,6)==2 j=j+1;B2(i,3)=B2(i,3)+DetaU(2*j,1)+DetaU(2*j-1,1)*sqrt(-1);end end for i=1:n %對(duì)PV節(jié)點(diǎn)的處理 if B2(i,6)==3 j=j+1;B2(i,3)=B2(i,3)+DetaU(2*j,1)+DetaU(2*j-1,1)*sqrt(-1);end end B2 %開(kāi)始循環(huán)********************************************************************** while abs(max(DetaU))>pr OrgS=zeros(2*n-2,1);%！!初始功率參數(shù)在迭代過(guò)程中是不累加的，所以在這里必須將其初始化為零矩陣 h=0;j=0;for i=1:n if i~=isb&B2(i,6)==2 h=h+1;for j=1:n OrgS(2*h-1,1)=OrgS(2*h-1,1)+real(B2(i,3))*(real(Y(i,j))*real(B2(j,3))-imag(Y(i,j))*imag(B2(j,3)))+imag(B2(i,3))*(real(Y(i,j))*imag(B2(j,3))+imag(Y(i,j))*real(B2(j,3)));OrgS(2*h,1)=OrgS(2*h,1)+imag(B2(i,3))*(real(Y(i,j))*real(B2(j,3))-imag(Y(i,j))*imag(B2(j,3)))-real(B2(i,3))*(real(Y(i,j))*imag(B2(j,3))+imag(Y(i,j))*real(B2(j,3)));end end end for i=1:n if i~=isb&B2(i,6)==3 h=h+1;for j=1:n OrgS(2*h-1,1)=OrgS(2*h-1,1)+real(B2(i,3))*(real(Y(i,j))*real(B2(j,3))-imag(Y(i,j))*imag(B2(j,3)))+imag(B2(i,3))*(real(Y(i,j))*imag(B2(j,3))+imag(Y(i,j))*real(B2(j,3)));OrgS(2*h,1)=OrgS(2*h,1)+imag(B2(i,3))*(real(Y(i,j))*real(B2(j,3))-imag(Y(i,j))*imag(B2(j,3)))-real(B2(i,3))*(real(Y(i,j))*imag(B2(j,3))+imag(Y(i,j))*real(B2(j,3)));end end end OrgS %創(chuàng)建DetaS h=0;for i=1:n if i~=isb&B2(i,6)==2 h=h+1;DetaS(2*h-1,1)=real(B2(i,2))-OrgS(2*h-1,1);DetaS(2*h,1)=imag(B2(i,2))-OrgS(2*h,1);end end t=0;for i=1:n if i~=isb&B2(i,6)==3 h=h+1;t=t+1;DetaS(2*h-1,1)=real(B2(i,2))-OrgS(2*h-1,1);DetaS(2*h,1)=real(PVU(t,1))^2+imag(PVU(t,1))^2-real(B2(i,3))^2-imag(B2(i,3))^2;end end DetaS %創(chuàng)建I i=zeros(n-1,1);h=0;for i=1:n if i~=isb h=h+1;I(h,1)=(OrgS(2*h-1,1)-OrgS(2*h,1)*sqrt(-1))/conj(B2(i,3));end end I %創(chuàng)建Jacbi Jacbi=zeros(2*n-2);h=0;k=0;for i=1:n if B2(i,6)==2 h=h+1;for j=1:n if j~=isb k=k+1;if i==j Jacbi(2*h-1,2*k-1)=-imag(Y(i,j))*real(B2(i,3))+real(Y(i,j))*imag(B2(i,3))+imag(I(h,1));Jacbi(2*h-1,2*k)=real(Y(i,j))*real(B2(i,3))+imag(Y(i,j))*imag(B2(i,3))+real(I(h,1));Jacbi(2*h,2*k-1)=-Jacbi(2*h-1,2*k)+2*real(I(h,1));Jacbi(2*h,2*k)=Jacbi(2*h-1,2*k-1)-2*imag(I(h,1));else Jacbi(2*h-1,2*k-1)=-imag(Y(i,j))*real(B2(i,3))+real(Y(i,j))*imag(B2(i,3));Jacbi(2*h-1,2*k)=real(Y(i,j))*real(B2(i,3))+imag(Y(i,j))*imag(B2(i,3));Jacbi(2*h,2*k-1)=-Jacbi(2*h-1,2*k);Jacbi(2*h,2*k)=Jacbi(2*h-1,2*k-1);end if k==(n-1)k=0;end end end end end k=0;for i=1:n if B2(i,6)==3 h=h+1;for j=1:n if j~=isb k=k+1;if i==j Jacbi(2*h-1,2*k-1)=-imag(Y(i,j))*real(B2(i,3))+real(Y(i,j))*imag(B2(i,3))+imag(I(h,1));Jacbi(2*h-1,2*k)=real(Y(i,j))*real(B2(i,3))+imag(Y(i,j))*imag(B2(i,3))+real(I(h,1));Jacbi(2*h,2*k-1)=2*imag(B2(i,3));Jacbi(2*h,2*k)=2*real(B2(i,3));else Jacbi(2*h-1,2*k-1)=-imag(Y(i,j))*real(B2(i,3))+real(Y(i,j))*imag(B2(i,3));Jacbi(2*h-1,2*k)=real(Y(i,j))*real(B2(i,3))+imag(Y(i,j))*imag(B2(i,3));Jacbi(2*h,2*k-1)=0;Jacbi(2*h,2*k)=0;end if k==(n-1)k=0;end end end end end Jacbi DetaU=zeros(2*n-2,1);DetaU=inv(Jacbi)*DetaS;DetaU %修正節(jié)點(diǎn)電壓 j=0;for i=1:n if B2(i,6)==2 j=j+1;B2(i,3)=B2(i,3)+DetaU(2*j,1)+DetaU(2*j-1,1)*sqrt(-1);end end for i=1:n if B2(i,6)==3 j=j+1;B2(i,3)=B2(i,3)+DetaU(2*j,1)+DetaU(2*j-1,1)*sqrt(-1);end end B2 Times=Times+1;%迭代次數(shù)加1 end Times 一個(gè)原始數(shù)據(jù)的例子 節(jié)點(diǎn)數(shù) 5 支路數(shù) 5平衡節(jié)點(diǎn)編號(hào) 5 精度pr 0.000001 B1(支路參數(shù)矩陣)[1 2 0.04+0.25i 0.5i 1 0;1 3 0.1+0.35i 0 1 0;2 3 0.08+0.30i 0.5i 1 0;4 2 0.015i 0 1.05 1;5 3 0.03i 0 1.05 1] B2(節(jié)點(diǎn)參數(shù)矩陣)[0-1.6-0.8i 1 0 0 2;0-2-1i 1 0 0 2;0-3.7-1.3i 1 0 0 2;0 5+0i 1.05 1.05 0 3;0 0 1.05 1.05 0 1] X(節(jié)點(diǎn)號(hào)和對(duì)地參數(shù))[1 0;2 0;3 0;4 0;5 0]

電力系統(tǒng)潮流計(jì)算

——9結(jié)點(diǎn)算例-PQ法

原始數(shù)據(jù)錄入data.txt文檔:

標(biāo)號(hào),起始結(jié)點(diǎn),終止結(jié)點(diǎn),支路電阻參數(shù),支路電抗參數(shù),支路對(duì)地導(dǎo)納參數(shù) 1,2,5,0.0,0.063,0.0, 2,5,9,0.019,0.072,0.075, 3,6,9,0.012,0.101,0.105, 4,3,6,0.0,0.059,0.0, 5,6,8,0.039,0.17,0.179, 6,4,8,0.017,0.092,0.079, 7,5,7,0.032,0.161,0.153, 8,4,7,0.01,0.085,0.088, 9,1,4,0.0,0.058,0.0, 潮流程序chaoliu.txt文檔: #include #include #define N 9 /*總結(jié)點(diǎn)數(shù)*/ #define M 6 /*PQ結(jié)點(diǎn)數(shù)*/ #define K 9 /*線(xiàn)路數(shù)*/ #define eps 1e-4 void guass(int n,int m,float c[],float b[][N],float x[])/*高斯函數(shù)*/ { float a[N][N],y[N];int i,j,k;for(i=0;i=0;i--){ for(j=i+1;jHeadnode-1;j=t->Endnode-1;lr=t->R;lx=t->X;lb1=t->b;lg=lr/(lr*lr+lx*lx);lb=-lx/(lr*lr+lx*lx);g[i][i]+=lg;g[j][j]+=lg;b[i][i]+=lb+lb1;b[j][j]+=lb+lb1;h[i][j]=h[j][i]=-lb1;g[i][j]-=lg;g[j][i]-=lg;b[i][j]-=lb;b[j][i]-=lb;} getch();printf(“n=====jie dian dao na ju zhen=====n”);for(i=0;i

;

float x1[N-1],x2

;for(i=1;imax)max=fabs(dpu[i]);} if(max>=eps){ for(i=0;i

guass(1,N-1,y1,B,x1);for(i=1;imax)max=fabs(dqu[i]);} if(max>=eps){ for(i=0;i

guass(N-M,M,y2,B,x2);for(i=N-M;i

else { kp=0;if(kq==0)val(u,g,b,r,ku,kr,h);else goto top;} } } void val(float u[N],float g[N][N],float b[N][N],float r[N],int ku, int kr,float h[N][N]){ float ps=0,pv1=0,pv2=0;float qs=0,qv1=0,qv2=0;float p[N][N]={0};float q[N][N]={0};float s[N][N];float dp[N][N]={0};float dq[N][N]={0};float ds[N][N];float dSp=0,dSq=0;int i,j;FILE *fp1;printf(”n=====ping heng jie dian gong lv =====n“);getch();for(i=0;i

printf(”n=======shu ju bao cun=====n“);fp1=fopen(”jieguo.txt“,”w+“);{ fprintf(fp1,”xian lu cao liu:n“);for(i=0;i

第四篇：基于MATLAB的電力系統(tǒng)潮流計(jì)算_張寧

DOI :10.13207

/j.cnki.jnwafu.20

4.12.028

第 32卷 第 12期 402 年 12月

西北農(nóng)林科技大學(xué)學(xué)報(bào)(自然科學(xué)版)

ruoJ.of

No

rthw

est

Sci-Tech

Univ.o

f

Ag

ri.and

F

or.(Nat.Sci.Ed.)loV.32 No.12

Dec.204

基于 MATLAB的電力系統(tǒng)潮流計(jì)算

1張 寧,江紅梅,張 渭2

(1西北農(nóng)林科技大學(xué)水利與建筑工程學(xué)院,陜西楊凌 71210

;2哈密市高級(jí)中學(xué)計(jì)算機(jī)室 ,新疆哈密 83900)

[摘 要] 隨著計(jì)算機(jī)語(yǔ)言技術(shù)的不斷發(fā)展和成熟 ,基于 MAT

LA

B的潮流計(jì)算研究近年來(lái)得到了長(zhǎng)足的發(fā) 展。針對(duì)這一現(xiàn)狀,以 P-Q分解法為例,分析了 BA

SIC , FO RT RA

N和 MA

T

LAB高級(jí)語(yǔ)言潮流計(jì)算的異同 ,指出了 其優(yōu)缺點(diǎn),并針對(duì)潮流計(jì)算模型結(jié)構(gòu)的特點(diǎn) ,提出了基于 MA

T LAB的潮流算法。

[關(guān)鍵詞] 電力系統(tǒng);潮流計(jì)算;MA

T LAB [中圖分類(lèi)號(hào)] T

M715;T M744

[文獻(xiàn)標(biāo)識(shí)碼] A

[文章編號(hào)] 1671-9387(204)12-0124-03

潮流計(jì)算是電力系統(tǒng)規(guī)劃、運(yùn)行的基本研究方

法,隨著現(xiàn)代電力系統(tǒng)大系統(tǒng)、強(qiáng)非線(xiàn)性與多元件特 點(diǎn)的日益突出,其計(jì)算量與計(jì)算復(fù)雜度急劇增

[ 1]

加。在處理潮流計(jì)算時(shí),其計(jì)算機(jī)軟件的速度已

[ 2] 無(wú)法滿(mǎn)足大電網(wǎng)模擬和實(shí)時(shí)控制的仿真要求 ,而 高效的潮流問(wèn)題相關(guān)軟件的研究已成為大規(guī)模電力

3,4] 系統(tǒng)仿真計(jì)算的關(guān)鍵[。隨著計(jì)算機(jī)技術(shù)的不斷 發(fā)展和成熟,對(duì) MA

T

LA B潮流計(jì)算的研究為快速、詳細(xì)地解決大電網(wǎng)仿真技術(shù)問(wèn)題開(kāi)辟了新思 路。針對(duì)這一現(xiàn)狀,本文以某電力網(wǎng)絡(luò)為例,分 析了 BASIC、FO

RT RA

N和 M A

T

LAB 高級(jí)語(yǔ)言潮 流計(jì)算的異同,并提出了基于 MA

T LAB 的潮流算 法,以期為電力系統(tǒng)潮流計(jì)算提供參考。[ 5,6] 的有效方法。該方法把非線(xiàn)性方程線(xiàn)性化,由于線(xiàn)

性方程的系數(shù)矩陣結(jié)構(gòu)上是稀疏的非對(duì)稱(chēng)矩陣,結(jié) 合稀疏矩陣技術(shù)可使計(jì)算機(jī)內(nèi)存占用量大大減少 , 計(jì)算速度大大加快;P-Q分解法[

8]是在 N

ew

ton jQ ,(i =1 , 2 ,3 ,…, n)(1)

ji j∑ =P j=1 U＊

i nj=1

計(jì)算機(jī)潮流計(jì)算的基本要求是:計(jì)算方法具有

[ 11 , 12 ]

∑·

ji jZ I = 式中,P iQ 分別為節(jié)點(diǎn) i向網(wǎng)絡(luò)注入的有功功率，·i和無(wú)功功率;jU為節(jié)點(diǎn) j的電壓相量 ·＊ P-Q

i,(i =1 ,2 ,3 ,…, n)(2)

ji·

＊

Ui

一定的可靠性和收斂性;盡量選用占計(jì)算機(jī)內(nèi)

存較少的存儲(chǔ)方式;在可靠收斂的前提下,選用計(jì)算 速度較快的方法;人機(jī)交互環(huán)境便于數(shù)據(jù)輸入、校核 和修改,且具有一定的靈活性。

BAS

IC是一種解決數(shù)學(xué)問(wèn)題的語(yǔ)言,由于其取

;Uii

為節(jié)點(diǎn)

消了編譯工作及連接過(guò)程,語(yǔ)言的執(zhí)行時(shí)間相應(yīng)較

·

AS

IC語(yǔ)句的標(biāo)識(shí)符可為常數(shù)和變量,矩陣也 的電壓共軛相量;jj的電流相量;I為節(jié)點(diǎn) Y為節(jié) 慢。Bji點(diǎn)導(dǎo)納矩陣;Z為節(jié)點(diǎn)阻抗矩陣。式(1)和式(2)各

ji

[ 13 ]

可規(guī)定為變量進(jìn)行計(jì)算。

有 n個(gè)非線(xiàn)性復(fù)數(shù)方程,對(duì)其作不同的應(yīng)用和處 FO RT RAN 是 1957年發(fā)明的用于科學(xué)計(jì)算的

語(yǔ)言理,就形成了不同的潮流計(jì)算方法。其中, N ew

ton個(gè)工程數(shù)學(xué)運(yùn)算函數(shù),可實(shí)現(xiàn)潮流計(jì)算中的矩陣求.70141E

+38

[ 15 ]。要達(dá)到題目所要求的計(jì)算精積、求逆、稀疏矩陣形成、復(fù)數(shù)運(yùn)算以及初等數(shù)學(xué)運(yùn) 算等[ 1]。MA

T

LA B語(yǔ)言允許用戶(hù)以數(shù)學(xué)形式的語(yǔ) 度，且兼顧矩陣程序設(shè)計(jì)的難易程度 ,MA

T LA

B則成為

言編寫(xiě)程序 ,其比 BASIC 語(yǔ)言和 FO RT RA

N等更

首選潮流計(jì)算的計(jì)算機(jī)語(yǔ)言。對(duì)于大矩陣的潮流計(jì) 為接近書(shū)寫(xiě)的數(shù)學(xué)表達(dá)格式 ,且程序易調(diào)試。在計(jì)

算 ,M A

T LAB提供的 M文件可將輸入矩陣按格式

算要求相同的情況下 ,使用 MA

T LAB編程 ,工作量

先寫(xiě)入一個(gè)文本文件,在編程時(shí) ,按文本文件名在命

將會(huì)大為減少。

令窗調(diào)用 ,大矩陣就被輸入到內(nèi)存中了。MA

T LAB 3

基于 MA

T LAB的電力系統(tǒng)潮流計(jì) 也能直接創(chuàng)建復(fù)數(shù)矩陣,這兩點(diǎn)也是 BA

SIC

,(2)求系數(shù)矩陣 B

′, B ″的逆陣,計(jì)算各節(jié)點(diǎn)電壓

FO RT

RA

N 算 的相位角 語(yǔ)言不可比擬的。Δδ。相位角 iΔiδ計(jì)算公式為

:

P-Q分解法計(jì)算電力系統(tǒng)的潮流分布 ,其

試用

-(k)

i步驟為:

-(B ′)(ΔP

(k)/ U

(k))=(U

式中 , ΔP

(k), U

(k)為已知量,相位角 Δδ的計(jì)算主(1)形成系數(shù)矩陣 B ′, B ″(設(shè)該矩陣為

4×4矩

i要(k)Δδ)是矩陣的求積運(yùn)陣)。

算

。AMA

T L

B和 BA

SIC語(yǔ)言都用命令 inv

(B)或 a

1a

2…

a 14 IN

V

(B)來(lái)實(shí)現(xiàn)矩陣的求逆功能,B

′=B

″= a

a

…

MA

T LA

B和 … … … …

a 24 BAS

IC語(yǔ)言求積運(yùn)算均與數(shù)學(xué)書(shū)寫(xiě)格式一致 ,但 a

a

…

BAS

IC語(yǔ)言必須對(duì)矩陣進(jìn)行 MA

T說(shuō)明 ,而且執(zhí)行

a 44

在程序設(shè)計(jì)上 ,用 MA

T LA

B語(yǔ)言編寫(xiě)程序形

時(shí)間相應(yīng)較長(zhǎng) ,計(jì)算精度較低。FO RT

RAN 語(yǔ)言雖

成的系數(shù)矩陣 B

′和B

″,遠(yuǎn)比 BASIC

,FO

RT RA

N簡(jiǎn)

然在參考資料上提供矩陣求逆、求積程序[ 1618

] ,但

單,矩陣輸入、輸出與數(shù)學(xué)書(shū)寫(xiě)格式相似。即在命令

其按線(xiàn)性代數(shù)的矩陣求逆、求積步驟編程[ ] ,編寫(xiě)窗口輸入

: 的程序至少需用 1個(gè)三重循環(huán)語(yǔ)句 , 所以與

>>B

… a

(3)相比計(jì)算平衡節(jié)點(diǎn)功率和線(xiàn)路功率MA

T LA

B ,工作量大而效率低。

=[

a

1a

… a 2 2a

a 1平衡節(jié)點(diǎn)功率為 …

a

4]

再在命令窗口輸入 >>B

=,則窗口將顯示出 B

′, B ″矩陣。

… a

… … … a

4

2￣S

s

=

+Q

線(xiàn)路功率為

·

＊·

U s is

i =P 1 ＊

Y

U

i =1

n

∑

￣S

ij

= U

iI

=P

ji

ji

+Q

ij

平衡節(jié)點(diǎn)功率和線(xiàn)路功率的計(jì)算屬?gòu)?fù)數(shù)運(yùn)算。MA

T LA

B和 FO RT

RAN 語(yǔ)言都提供復(fù)數(shù)功能語(yǔ) 在 BASIC

, FO RT RAN 語(yǔ)言編程時(shí) , 必須進(jìn)行

句[ 15 ],MA

T LAB 以數(shù)學(xué)上的復(fù)數(shù)書(shū)寫(xiě)格式編寫(xiě)程 矩陣說(shuō)明 , BASIC 用 DIM B(4, 4)語(yǔ)句說(shuō)明矩陣，序 ,方便而且不容易出錯(cuò);而 FO RT

RA

N必須對(duì)復(fù)

F ORT

RA

N用 DOU

BLE

PRECISION

B(1 ∶4

, 1 ∶4

數(shù)變量進(jìn)行)

CMP LX語(yǔ)句的復(fù)數(shù)類(lèi)型說(shuō)明[

],顯

語(yǔ)句同時(shí)說(shuō)明矩陣和矩陣精度;在建立矩陣時(shí), 得較為繁瑣。BA

SIC不提供復(fù)數(shù)功能語(yǔ)句 ,一般將實(shí)

BASIC

,F

ORT RA

N語(yǔ)言用 REA

D語(yǔ)句和 DA

TA 語(yǔ) 數(shù)、虛數(shù)分開(kāi)計(jì)算 ,最后用輸出語(yǔ)句寫(xiě)在一起 ,所以 句輸入矩陣;FO RT

RAN 語(yǔ)言用 WRIT E語(yǔ)句輸出 一般不用 BASIC 語(yǔ)言編寫(xiě)復(fù)數(shù)運(yùn)算程序

。矩陣,而 BASIC 語(yǔ)言只能用二重循環(huán)語(yǔ)句輸出二維

另外 ,MA

T LA

B可以提供潮流計(jì)算中稀疏矩陣

· ·＊

126

西北農(nóng)林科技大學(xué)學(xué)報(bào)(自然科學(xué)版)

第 32卷 的建立命令 spconvetr,從而將外部數(shù)據(jù)轉(zhuǎn)化為稀疏 矩陣,而且 MA

T

LA B函數(shù) lu也可直接實(shí)現(xiàn) L

R三 角分解[了

方便。

一步。矩陣輸入、輸出格式簡(jiǎn)單 ,與數(shù)學(xué)書(shū)寫(xiě)格式相

似;以雙精度類(lèi)型進(jìn)行數(shù)據(jù)的存儲(chǔ)和運(yùn)算 ,數(shù)據(jù)精確

;潮流計(jì)算中復(fù)雜矩陣的輸入問(wèn)題可通過(guò)創(chuàng)建 M文件來(lái)解決;MA

T LA B稱(chēng)為矩陣實(shí)驗(yàn)室 ,其能進(jìn) 行潮流計(jì)算中的各種矩陣運(yùn)算,包括求逆、求積和矩 19 , 20

],從而為現(xiàn)代電力系統(tǒng)潮流計(jì)算提供度高 討

論

陣 L

R分解等 ,其程序的編寫(xiě)也因 M A

T LAB提供 了許多功能函數(shù)而變得簡(jiǎn)單易行。另外 , MA

T LAB

稀疏矩陣技術(shù)的引入 ,使電力系統(tǒng)潮流計(jì)算由傳統(tǒng) 方法轉(zhuǎn)變?yōu)閮?yōu)化算法成為可能。通過(guò)與 BA

SIC, FO RT RA

N語(yǔ)言的比較,基于 M A

T LAB的電力系統(tǒng)潮流計(jì)算使計(jì)算機(jī)在計(jì)算、分 析、研究復(fù)雜的電力系統(tǒng)潮流分布問(wèn)題上又前進(jìn)了

[參考文獻(xiàn)] [ 1 ]

薛

巍 ,舒繼武 ,王心豐 ,等.電力系統(tǒng)并行算法的研究進(jìn)展 [ J ]

.清華大學(xué)學(xué)報(bào)(自然科學(xué)版), 2002

, 42

(9):1192

-195

.[ 2 ]

F

ernaod

L A

.Computalnoi

complexit

y in

pow

er

systems

[ J ]

.IE

E

T

rans

on

PAS, 1976 , 95(4):1028

-1037

.[ 3 ]

周濟(jì) ,羅應(yīng)立 ,張建華 ,等.基于 MAT

LA B的次同步諧振的特征根分析 [ J ]

.現(xiàn)代電力 , 1999 ,(4):67-71.[ 4 ]

王克英 ,穆

鋼 ,韓學(xué)山 ,等.使潮流方程直接可解的 PMV配置方案研究 [ J ]

.中國(guó)電機(jī)工程學(xué)報(bào) , 1999 , 10

(2):67-69.[ 5]

盧天成 ,韓

或 ,李東風(fēng).MA

T LAB在鳴禽語(yǔ)圖分析中的應(yīng)用 [

J ].遼寧師范大學(xué)學(xué)報(bào)(自然科學(xué)版), 2001

, 24

(4):56-58.[ 6]

宋維君.基于 MA

T LAB的時(shí)域分析法實(shí)現(xiàn) [ J ]

.遼寧師專(zhuān)學(xué)報(bào)(自然科學(xué)版), 2002, 1(4):72

-74.[ 9]

邱關(guān)源.電路(下冊(cè))[ M]

.北京 :高等教育出版[ ,7] 0 周孝信等..國(guó)電力百科全書(shū) ·電力系統(tǒng)卷 [ M].第 2版.北京 :中國(guó)電力出版社 , 2001

.31

.社 199

.67-78 [ 10

]

趙晉泉 ,侯志儉 ,吳際舜.牛頓最優(yōu)潮流算法中離散控制量的新處理方法 [

J ].電力系統(tǒng)自動(dòng)化 , 1999 , 23

(3):32

[ 8]

陳

珩.電力系統(tǒng)穩(wěn)態(tài)分析 [ M]

.北京 :高等教育出版社 ,20

.139

-193

.-34

.[ 11

]

[英]勞頓 M A,塞 M G.電氣工程師技術(shù)手冊(cè) [

M]

.北京 :機(jī)械工業(yè)出版社 , 1992

.826

-827

.[ 12

]

樓順天 ,陳生潭 ,雷虎明.MA

T LAB5

.X程序設(shè)計(jì)語(yǔ)言 [ M]

.西安 :西安電子科技大學(xué)出版社 , 2001

.57-78;163

-162

.[ 14

]

潭浩強(qiáng) ,田淑清.F

ORT

RA

N語(yǔ)言 ———F ORT

RA

N 77結(jié)構(gòu)化程序設(shè)計(jì) [

M].北京 :清華大學(xué)出版[ 13

]

王沫然.MA

T LAB

6.0與科學(xué)計(jì)算[

M]

.北京 :電子工業(yè)出版社 , 2001

.1

-5;42

.社 , 1995.58;98.[ 15 ]

劉文遠(yuǎn) ,李興成 ,趙

莉.BA

SIC語(yǔ)言與農(nóng)業(yè)實(shí)用程序 [ 041

.[ 16 ]

潭浩強(qiáng) ,崔武子 ,田淑清.FO RT

RA

N程序設(shè)計(jì)(二級(jí))教程[

社 , 1998.105

-123

.[ 17 ]

鮑有文 ,周海燕 ,趙重敏 ,等.FO RT

RAN

77程序設(shè)計(jì)試題匯編 [ M]

.北京 :清華大學(xué)出版[ 20

]

鄧建中 ,葛仁杰 ,程正興.計(jì)算方法[ M]

.西安 :西安交通大學(xué)出版社 , 1998

.21

社 ,.68;107

.-31

1998;201

.[ 18 ]

劉萬(wàn)春 ,朱玉文 ,龔圓明.FO RT

RA

N程序設(shè)計(jì)[

M].北京 :國(guó)防工業(yè)出版社 , 2003

.102

-124

.[ 19 ]

同濟(jì)大學(xué)教研室.線(xiàn)性代數(shù)[

M].北京 :高等教育出版社 , 1990

.35

;70-71

.M].北京 :清華大學(xué)出版

M]

.沈陽(yáng) :遼寧科學(xué)技術(shù)出版社 , 1987.139

-Pow

er

flow

computaino

of

the

electir

pow

er

system

based

on

MA

TLA

B

1GNAHZ

Ning, JIANG

Hong

-mei

,ZHANG

Wei

(1

Colge

of

Water

Resocru

and

A

rchlarueti

E

ngiren

, Northwest

A

&F University, Yanglin, Shanxi

71210

, China

;

Computer

Rom

of

High

Schol

in

Hami, Ham i, Xinj

iang

83900, China)

Acartsb

:Wit

h

the

devolpment

of

the

computer

langseu

in

rectn

years

,t

he

reshcra

of

M A

T LAB

power wolf

algohmtirs

got

substanli

devolpment

.T

he

paper

, using

P

-Q

method

, analy

zed

the

dif

erensc

andimsilarseti

of

BA

SIC

, FO RT RAN

and

MA

T LAB

advaecn

langseu, and

pointed

out

t hei

r

advatn

ages

and M A

T LAB

pow

er

flow

algoir

thm

s.yeK

words

:electir

pow

er

sy

stem

;power

f low

computaino

;MA

T LAB w

eaknes

, based

on

the

charsietc

of

the

model

const

ructnoi

in

pow

er

flow

algohmtirs

, bring

up

t he

第五篇：電力系統(tǒng)潮流計(jì)算

南 京 理 工 大 學(xué)

《電力系統(tǒng)穩(wěn)態(tài)分析》

課程報(bào)告

姓名

XX

學(xué) 號(hào)： 5*** 自動(dòng)化學(xué)院 電氣工程

基于牛頓-拉夫遜法的潮流計(jì)算例題編程報(bào)學(xué)院(系): 專(zhuān)

業(yè): 題

目: 任課教師 碩士導(dǎo)師 告

楊偉 XX

2015年6月10號(hào)

基于牛頓-拉夫遜法的潮流計(jì)算例題編程報(bào)告

摘要：電力系統(tǒng)潮流計(jì)算的目的在于：確定電力系統(tǒng)的運(yùn)行方式、檢查系統(tǒng)中各元件是否過(guò)壓或者過(guò)載、為電力系統(tǒng)繼電保護(hù)的整定提供依據(jù)、為電力系統(tǒng)的穩(wěn)定計(jì)算提供初值、為電力系統(tǒng)規(guī)劃和經(jīng)濟(jì)運(yùn)行提供分析的基礎(chǔ)。潮流計(jì)算的計(jì)算機(jī)算法包含高斯—賽德?tīng)柕?、牛頓-拉夫遜法和P—Q分解法等，其中牛拉法計(jì)算原理較簡(jiǎn)單、計(jì)算過(guò)程也不復(fù)雜，而且由于人們引入泰勒級(jí)數(shù)和非線(xiàn)性代數(shù)方程等在算法里從而進(jìn)一步提高了算法的收斂性和計(jì)算速度。同時(shí)基于MATLAB的計(jì)算機(jī)算法以雙精度類(lèi)型進(jìn)行數(shù)據(jù)的存儲(chǔ)和運(yùn)算, 數(shù)據(jù)精確度高，能進(jìn)行潮流計(jì)算中的各種矩陣運(yùn)算，使得傳統(tǒng)潮流計(jì)算方法更加優(yōu)化。

一 研究?jī)?nèi)容

通過(guò)一道例題來(lái)認(rèn)真分析牛頓-拉夫遜法的原理和方法（采用極坐標(biāo)形式的牛拉法），同時(shí)掌握潮流計(jì)算計(jì)算機(jī)算法的相關(guān)知識(shí)，能看懂并初步使用MATLAB軟件進(jìn)行編程，培養(yǎng)自己電力系統(tǒng)潮流計(jì)算機(jī)算法編程能力。

例題如下：用牛頓-拉夫遜法計(jì)算下圖所示系統(tǒng)的潮流分布，其中系統(tǒng)中5為平衡節(jié)點(diǎn)，節(jié)點(diǎn)5電壓保持U=1.05為定值，其他四個(gè)節(jié)點(diǎn)分別為PQ節(jié)點(diǎn)，給定的注入功率如圖所示。計(jì)算精度要求各節(jié)點(diǎn)電壓修正量不大于10-6。

二 牛頓-拉夫遜法潮流計(jì)算 1 基本原理

牛頓法是取近似解x(k)之后，在這個(gè)基礎(chǔ)上，找到比x(k)更接近的方程的根，一步步地迭代，找到盡可能接近方程根的近似根。牛頓迭代法其最大優(yōu)點(diǎn)是在方程f(x)=0的單根附近時(shí)誤差將呈平方減少，而且該法還可以用來(lái)求方程的重根、復(fù)根。電力系統(tǒng)潮流計(jì)算，一般來(lái)說(shuō)，各個(gè)母線(xiàn)所供負(fù)荷的功率是已知的，各個(gè)節(jié)點(diǎn)的電壓是未知的（平衡節(jié)點(diǎn)外）可以根據(jù)網(wǎng)絡(luò)結(jié)構(gòu)形成節(jié)點(diǎn)導(dǎo)納矩陣，然后由節(jié)點(diǎn)導(dǎo)納矩陣列寫(xiě)功率方程，由于功率方程里功率是已知的，電壓的幅值和相角是未知的，這樣潮流計(jì)算的問(wèn)題就轉(zhuǎn)化為求解非線(xiàn)性方程組的問(wèn)題了。為了便于用迭代法解方程組，需要將上述功率方程改寫(xiě)成功率平衡方程，并對(duì)功率平衡方程求偏導(dǎo)，得出對(duì)應(yīng)的雅可比矩陣，給未知節(jié)點(diǎn)賦電壓初值，將初值帶入功率平衡方程，得到功率不平衡量，這樣由功率不平衡量、雅可比矩陣、節(jié)點(diǎn)電壓不平衡量（未知的）構(gòu)成了誤差方程，解誤差方程，得到節(jié)點(diǎn)電壓不平衡量，節(jié)點(diǎn)電壓加上節(jié)點(diǎn)電壓不平衡量構(gòu)成節(jié)點(diǎn)電壓新的初值，將新的初值帶入原來(lái)的功率平衡方程，并重新形成雅可比矩陣，然后計(jì)算新的電壓不平衡量，這樣不斷迭代，不斷修正，一般迭代三到五次就能收斂。2 基本步驟和設(shè)計(jì)流程圖

形成了雅克比矩陣并建立了修正方程式，運(yùn)用牛頓-拉夫遜法計(jì)算潮流的核心問(wèn)題已經(jīng)解決，已有可能列出基本計(jì)算步驟并編制流程圖。由課本總結(jié)基本步驟如下：

1)形成節(jié)點(diǎn)導(dǎo)納矩陣Y；

2)設(shè)各節(jié)點(diǎn)電壓的初值，如果是直角坐標(biāo)的話(huà)設(shè)電壓的實(shí)部e和虛部f；如果是極坐標(biāo)的話(huà)則設(shè)電壓的幅值U和相角a；

3)將各個(gè)節(jié)點(diǎn)電壓的初值代入公式求修正方程中的不平衡量以及修正方程的系數(shù)矩陣的雅克比矩陣；

4)解修正方程式，求各節(jié)點(diǎn)電壓的變化量，即修正量； 5)計(jì)算各個(gè)節(jié)點(diǎn)電壓的新值，即修正后的值；

6)利用新值從第（3）步開(kāi)始進(jìn)入下一次迭代，直至達(dá)到精度退出循環(huán)； 7)計(jì)算平衡節(jié)點(diǎn)的功率和線(xiàn)路功率，輸出最后計(jì)算結(jié)果； ① 公式推導(dǎo)

② 流程圖

三

matlab編程代碼

clear;

% 如圖所示1，2，3，4為PQ節(jié)點(diǎn),5為平衡節(jié)點(diǎn)

y=0;

% 輸入原始數(shù)據(jù)，求節(jié)點(diǎn)導(dǎo)納矩陣

y(1,2)=1/(0.07+0.21j);

y(4,5)=0;y(1,3)=1/(0.06+0.18j);

y(1,4)=1/(0.05+0.10j);

y(1,5)=1/(0.04+0.12j);

y(2,3)=1/(0.05+0.10j);

y(2,5)=1/(0.08+0.24j);

y(3,4)=1/(0.06+0.18j);

for i=1:5

for j=i:5

y(j,i)=y(i,j);

end

end

Y=0;

% 求節(jié)點(diǎn)導(dǎo)納矩陣中互導(dǎo)納

for i=1:5

for j=1:5

if i~=j

Y(i,j)=-y(i,j);

end

end

end

% 求節(jié)點(diǎn)導(dǎo)納矩陣中自導(dǎo)納

for i=1:5

Y(i,i)=sum(y(i,:));

end

Y

% Y為導(dǎo)納矩陣

G=real(Y);

B=imag(Y);% 輸入原始節(jié)點(diǎn)的給定注入功率

S(1)=0.3+0.3j;

S(2)=-0.5-0.15j;

S(3)=-0.6-0.25j;

S(4)=-0.7-0.2j;

S(5)=0;

P=real(S);

Q=imag(S);

% 賦初值，U為節(jié)點(diǎn)電壓的幅值，a為節(jié)點(diǎn)電壓的相位角

U=ones(1,5);

U(5)=1.05;

a=zeros(1,5);

x1=ones(8,1);

x2=ones(8,1);

k=0;

while max(x2)>1e-6

for i=1:4

for j=1:4

H(i,j)=0;

N(i,j)=0;

M(i,j)=0;

L(i,j)=0;

oP(i)=0;

oQ(i)=0;

end

end

% 求有功、無(wú)功功率不平衡量

for i=1:4

for j=1:5

oP(i)=oP(i)-U(i)*U(j)*(G(i,j)*cos(a(i)-a(j))+B(i,j)*sin(a(i)-a(j)));

oQ(i)=oQ(i)-U(i)*U(j)*(G(i,j)*sin(a(i)-a(j))-B(i,j)*cos(a(i)-a(j)));

end

oP(i)=oP(i)+P(i);

oQ(i)=oQ(i)+Q(i);

end

x2=[oP,oQ]';

% x2為不平衡量列向量

% 求雅克比矩陣

% 當(dāng)i~=j時(shí)，求H,N,M,L

for i=1:4

for j=1:4

if i~=j

H(i,j)=-U(i)*U(j)*(G(i,j)*sin(a(i)-a(j))-B(i,j)*cos(a(i)-a(j)));

N(i,j)=-U(i)*U(j)*(G(i,j)*cos(a(i)-a(j))+B(i,j)*sin(a(i)-a(j)));

L(i,j)=H(i,j);

M(i,j)=-N(i,j);

end

end

end

% 當(dāng)i=j時(shí)，求H,N,M,L

for i=1:4

for j=1:5

if i~=j H(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)));N(i,i)=N(i,i)-U(i)*U(j)*(G(i,j)*cos(a(i)-a(j))+B(i,j)*sin(a(i)-a(j)));

M(i,i)=M(i,i)-U(i)*U(j)*(G(i,j)*cos(a(i)-a(j))+B(i,j)*sin(a(i)-a(j)));

L(i,i)=L(i,i)-U(i)*U(j)*(G(i,j)*sin(a(i)-a(j))-B(i,j)*cos(a(i)-a(j)))

end

end

N(i,i)=N(i,i)-2*(U(i))^2*G(i,i);

L(i,i)=L(i,i)+2*(U(i))^2*B(i,i);

end

J=[H,N;M,L]

% J為雅克比矩陣

x1=-((inv(J))*x2);

% x1為所求△x的列向量

% 求節(jié)點(diǎn)電壓新值，準(zhǔn)備下一次迭代

for i=1:4

oa(i)=x1(i);

oU(i)=x1(i+4)*U(i);

end

for i=1:4

a(i)=a(i)+oa(i);

U(i)=U(i)+oU(i);

end

k=k+1;

end

k,U,a

% 求節(jié)點(diǎn)注入功率

i=5;

for j=1:5

P(i)=U(i)*U(j)*(G(i,j)*cos(a(i)-a(j))+B(i,j)*sin(a(i)-a(j)))+P(i);

Q(i)=U(i)*U(j)*(G(i,j)*sin(a(i)-a(j))-B(i,j)*cos(a(i)-a(j)))+Q(i);

end

S(5)=P(5)+Q(5)*sqrt(-1);

S

% 求節(jié)點(diǎn)注入電流

I=Y*U'

四

運(yùn)行結(jié)果

節(jié)點(diǎn)導(dǎo)納矩陣

經(jīng)過(guò)五次迭代后的雅克比矩陣

迭代次數(shù)以及節(jié)點(diǎn)電壓的幅值和相角（弧度數(shù)）

節(jié)點(diǎn)注入功率和電流

五 結(jié)果分析

在這次學(xué)習(xí)和實(shí)際操作過(guò)程里：首先，對(duì)電力系統(tǒng)分析中潮流計(jì)算的部分特別是潮流計(jì)算的計(jì)算機(jī)算法中的牛頓-拉夫遜法進(jìn)行深入的研讀，弄明白了其原理、計(jì)算過(guò)程、公式推導(dǎo)以及設(shè)計(jì)流程。牛頓-拉夫遜法是求解非線(xiàn)性方程的迭代過(guò)程，其計(jì)算公式為?F?J?X，式中J為所求函數(shù)的雅可比矩陣；?X為需要求的修正值；?F為不平衡的列向量。利用x(*)=x(k+1)+?X(k+1)進(jìn)行多次迭代，通過(guò)迭代判據(jù)得到所需要的精度值即準(zhǔn)確值x(*)。六 結(jié)論

通過(guò)這個(gè)任務(wù)，自己在matlab編程，潮流計(jì)算，word文檔的編輯功能等方面均有提高，但也暴漏出一些問(wèn)題：理論知識(shí)儲(chǔ)備不足，對(duì)matlab的性能和特點(diǎn)還不能有一個(gè)全面的把握，對(duì)word軟件也不是很熟練，相信通過(guò)以后的學(xué)習(xí)能彌補(bǔ)這些不足，達(dá)到一個(gè)新的層次。