基于多項(xiàng)式正態(tài)變換和最大熵估計(jì)的含風(fēng)電系統(tǒng)電壓穩(wěn)定概率分析

謝應(yīng)昭，盧繼平

（重慶大學(xué) 輸配電裝備及系統(tǒng)安全與新技術(shù)國(guó)家重點(diǎn)實(shí)驗(yàn)室，重慶 400044）

0 引言

隨著風(fēng)力發(fā)電逐步向大規(guī)模、高集中開(kāi)發(fā)的方向發(fā)展，電力系統(tǒng)必然會(huì)面臨大量可再生能源的接入問(wèn)題。風(fēng)力發(fā)電的規(guī)模性接入，對(duì)電力系統(tǒng)運(yùn)行的各方面都會(huì)帶來(lái)明顯影響，電壓穩(wěn)定研究作為系統(tǒng)安全穩(wěn)定運(yùn)行的重要部分，將面臨更多的不確定因素。傳統(tǒng)電壓穩(wěn)定分析都是基于確定性模型，忽略了系統(tǒng)中的隨機(jī)注入功率的不確定性，以系統(tǒng)最危險(xiǎn)工作模式為研究對(duì)象進(jìn)行電壓穩(wěn)定性評(píng)估，評(píng)估結(jié)果過(guò)于保守，靈活性不強(qiáng)。為了彌補(bǔ)確定性方法的不足，在電力系統(tǒng)不確定問(wèn)題研究中引入概率分析方法，這些方法可分為以下3類(lèi)：蒙特卡洛MC（Monte Carlo）法、解析法以及點(diǎn)估計(jì)法。其中，蒙特卡洛法可以得到高精度的隨機(jī)變量統(tǒng)計(jì)信息，但由于計(jì)算效率低下，一般僅用作其他概率分析方法的驗(yàn)證標(biāo)準(zhǔn)；解析法通過(guò)卷積技術(shù)或輸出變量的累積量求解概率密度函數(shù)，這種方法假定輸入和輸出隨機(jī)變量成線性關(guān)系，不能反映系統(tǒng)運(yùn)行的實(shí)際情況；點(diǎn)估計(jì)法根據(jù)已知變量的概率分布求解未知變量的各階矩信息，計(jì)算量小、精度高，是一種理想的概率分析方法。上述3類(lèi)方法中，只有蒙特卡洛法可直接得出隨機(jī)變量概率分布，后面2類(lèi)方法還需要結(jié)合其他方法得出隨機(jī)變量的概率分布。

現(xiàn)有文獻(xiàn)主要集中于概率潮流計(jì)算問(wèn)題的研究［1-5］，關(guān)于電壓穩(wěn)定問(wèn)題概率分析的文獻(xiàn)較少［6-9］：文獻(xiàn)［6］將點(diǎn)估計(jì)法應(yīng)用于電壓穩(wěn)定概率分析的計(jì)算，但并未計(jì)及風(fēng)電功率的影響；文獻(xiàn)［7］研究了考慮分布式電源的靜態(tài)電壓穩(wěn)定概率問(wèn)題，研究過(guò)程忽略了輸入變量的相關(guān)性，用級(jí)數(shù)展開(kāi)求取概率分布時(shí)，存在一定誤差；文獻(xiàn)［8］采用隨機(jī)響應(yīng)面法對(duì)電壓穩(wěn)定概率問(wèn)題進(jìn)行求解，取得了較為理想的分析結(jié)果；文獻(xiàn)［9］對(duì)含風(fēng)電系統(tǒng)進(jìn)行了電壓穩(wěn)定概率分析，文中采用Nataf變換處理相關(guān)輸入變量，蒙特卡洛方法求解概率問(wèn)題，求解過(guò)程較為復(fù)雜。

本文提出了一種基于多項(xiàng)式正態(tài)變換PNT（Polynomial Normal Transformation）和最大熵估計(jì)的含風(fēng)電系統(tǒng)電壓穩(wěn)定概率分析方法。利用多項(xiàng)式正態(tài)變換方法處理輸入變量之間的相關(guān)性，結(jié)合點(diǎn)估計(jì)法將電壓穩(wěn)定概率分析問(wèn)題轉(zhuǎn)化為確定性負(fù)荷裕度非線性?xún)?yōu)化模型的求解，并計(jì)算負(fù)荷裕度的統(tǒng)計(jì)特征，采用最大熵估計(jì)方法估計(jì)輸出變量概率分布函數(shù)。文中還分析了風(fēng)速參數(shù)變化、相關(guān)系數(shù)矩陣變化、風(fēng)電接入容量變化以及風(fēng)電場(chǎng)功率因數(shù)變化對(duì)負(fù)荷裕度統(tǒng)計(jì)特征的影響。本文方法解決了點(diǎn)估計(jì)法無(wú)法處理輸入變量相關(guān)性和得出輸出變量概率分布的不足，保持了點(diǎn)估計(jì)法計(jì)算簡(jiǎn)便、精度高的優(yōu)點(diǎn)。算例分析表明，本文所提方法具有令人滿意的計(jì)算精度，且計(jì)算效率高，變換方法簡(jiǎn)潔靈活，能夠得到全面反映系統(tǒng)電壓穩(wěn)定性的概率信息。

1 多項(xiàng)式正態(tài)變換

概率問(wèn)題中包含不同概率分布的隨機(jī)變量，有些變量之間還具有相關(guān)性，分布形式的不同和變量的相關(guān)性都會(huì)對(duì)分析結(jié)果產(chǎn)生影響。因此有必要對(duì)輸入變量進(jìn)行預(yù)處理，以減小這些因素帶來(lái)的影響。多項(xiàng)式正態(tài)變換可以同時(shí)考慮多個(gè)非正態(tài)相關(guān)的隨機(jī)變量，利用統(tǒng)一的變換步驟將其變換到獨(dú)立正態(tài)空間，變換方法簡(jiǎn)潔有效，具有良好的通用性。

1.1 單變量正態(tài)變換原理

多項(xiàng)式正態(tài)變換是利用標(biāo)準(zhǔn)正態(tài)分布隨機(jī)變量z的多項(xiàng)式對(duì)非正態(tài)分布的隨機(jī)變量x進(jìn)行求解。隨機(jī)變量x的r階多項(xiàng)式正態(tài)變換表達(dá)式為：

其中，a0、a1、a2、…、ar為多項(xiàng)式變換系數(shù)。

多項(xiàng)式階數(shù)r可以取不同值，得到不同的多項(xiàng)式變換公式，其中最常用的有二階、三階和五階多項(xiàng)式。盡管高階多項(xiàng)式可利用高階矩信息提高多項(xiàng)式正態(tài)變換的近似度，但高階多項(xiàng)式對(duì)近似度的改善程度還不明確，而計(jì)算過(guò)程卻更加復(fù)雜［10］，因此階數(shù)的確定需要結(jié)合具體計(jì)算要求進(jìn)行選擇。除了階數(shù)的確定，正態(tài)變換的另一關(guān)鍵在于變換系數(shù)ai的求取。ai的求解算法主要包括積矩法、L-矩法、最小二乘法以及Fisher-Cornish級(jí)數(shù)展開(kāi)法，文獻(xiàn)［11］對(duì)這些算法做了詳細(xì)的數(shù)值分析和對(duì)比，可作為求解算法的選擇依據(jù)。結(jié)合已有文獻(xiàn)結(jié)論［10-11］，本文選擇三階多項(xiàng)式作為正態(tài)變換多項(xiàng)式，并采用L-矩法對(duì)變換系數(shù)進(jìn)行求解。

隨機(jī)變量 x 的 q 階 L-矩 λq（q=1，2，3，4）與多項(xiàng)式變換系數(shù) ai（i=0，1，2，3）滿足如式（2）所示函數(shù)關(guān)系［11］。

利用前4階L-矩得出變換系數(shù)，即可將隨機(jī)變量x轉(zhuǎn)換為標(biāo)準(zhǔn)正態(tài)變量z進(jìn)行處理。

本文采用Weibull分布對(duì)風(fēng)速進(jìn)行刻畫(huà)，可直接利用Weibull分布參數(shù)對(duì)前4階L-矩進(jìn)行計(jì)算［12］，如式（3）所示。

其中，τ3=λ3/λ2和τ4=λ4/λ2分別為 L-偏度和 L-峰度；kw和c分別為Weibull分布的形狀參數(shù)和尺度參數(shù)。通過(guò)式（3）即可完成服從Weibull分布的隨機(jī)變量多項(xiàng)式正態(tài)變換。

1.2 多變量多項(xiàng)式正態(tài)變換

多隨機(jī)變量概率分布模型的建立需同時(shí)考慮各變量的邊際分布和變量間的相關(guān)性因素。由1.1節(jié)可知，各邊際分布隨機(jī)變量均可采用單變量正態(tài)變換方法進(jìn)行變換處理，本節(jié)解決的是把原始隨機(jī)變量的相關(guān)系數(shù)矩陣變換為對(duì)應(yīng)的標(biāo)準(zhǔn)正態(tài)變量的相關(guān)系數(shù)矩陣，實(shí)現(xiàn)多變量的多項(xiàng)式正態(tài)變換。

假設(shè) X=［x1，x2，…，xn］T表示由 n 個(gè)隨機(jī)變量 xi組成的隨機(jī)向量，對(duì)應(yīng)的相關(guān)系數(shù)矩陣為Z=［z1，z2，…，zn］T表示由 n 個(gè)標(biāo)準(zhǔn)正態(tài)變量 zi組成的隨機(jī)向量，對(duì)應(yīng)的相關(guān)系數(shù)矩陣為則兩相關(guān)系數(shù)矩陣對(duì)應(yīng)元素滿足如下關(guān)系［11］：

其中，ami和 amj（m=0，1，2，3）分別為 xi和 xj的多項(xiàng)式變換系數(shù)分別為xi和xj的標(biāo)準(zhǔn)差和期望。依次求解RX中各元素對(duì)應(yīng)的等式方程，并選擇滿足不等式條件的解作為ρzizj的值，即可得到標(biāo)準(zhǔn)正態(tài)變量相關(guān)系數(shù)矩陣RZ。

結(jié)合1.1節(jié)和1.2節(jié)，可以得到由獨(dú)立標(biāo)準(zhǔn)正態(tài)隨機(jī)變量 S=［s1，s2，…，sn］T求取原始隨機(jī)變量的完整公式：

其中，Am=［am1，am2，…，amn］T（m=0，1，2，3）為由 1.1 節(jié)得出的n個(gè)隨機(jī)變量xi對(duì)應(yīng)的多項(xiàng)式變換系數(shù)；L為標(biāo)準(zhǔn)正態(tài)變量相關(guān)系數(shù)矩陣RZ經(jīng)過(guò)Cholesky分解后所得的下三角矩陣。

通過(guò)多項(xiàng)式正態(tài)變換，為點(diǎn)估計(jì)法處理含相關(guān)性的隨機(jī)變量提供數(shù)據(jù)基礎(chǔ)。

2 點(diǎn)估計(jì)法

點(diǎn)估計(jì)法是求解非線性函數(shù)Y=h（X）在隨機(jī)變量 X=［x1，x2，…，xn］T作用下的統(tǒng)計(jì)分布的一種有效方法，具有良好的計(jì)算精度和較小的計(jì)算復(fù)雜度。其原理是首先求解一個(gè)r×n的采樣點(diǎn)集F1（X）和每個(gè)采樣點(diǎn)對(duì)應(yīng)的權(quán)重系數(shù)，然后利用式（6）計(jì)算函數(shù)變量Y的前2r-1矩信息：

其中，（ μ1，μ2，…，xi，k，…，μn）和 ωi，k分別為隨機(jī)變量采樣點(diǎn)和對(duì)應(yīng)的權(quán)重系數(shù)，xi，k為X的某一分量估計(jì)值， μl（l≠i）為其他分量的自身期望；r和 n分別為隨機(jī)變量xi的估計(jì)點(diǎn)數(shù)和變量個(gè)數(shù)。研究表明當(dāng)r＞3時(shí)，采樣點(diǎn)集F1（X）和權(quán)重系數(shù)集ω可能包含復(fù)數(shù)解，不符合實(shí)際情況，因此r的數(shù)值通常取為3。

將隨機(jī)變量xi的一個(gè)估計(jì)點(diǎn)選為自身期望μi的三點(diǎn)估計(jì)法稱(chēng)為2n+1估計(jì)方案。該方案的采樣點(diǎn)集 F1（X）中有 n 個(gè)相同采樣點(diǎn)（ μ1，μ2，…，μi，…，μn），只需進(jìn)行一次計(jì)算即可完成對(duì)這n個(gè)采樣點(diǎn)的處理，具有很高的計(jì)算效率。2n+1估計(jì)方案下的位置系數(shù) ξi，k和對(duì)應(yīng)權(quán)重系數(shù) ωi，k計(jì)算公式如下［13］：

其中，i=1，2，…，n；λi，3和 λi，4分別為隨機(jī)變量 xi的偏度系數(shù)和峰度系數(shù)；σi和μi分別為xi的標(biāo)準(zhǔn)差和期望。

點(diǎn)估計(jì)法適用于輸入隨機(jī)變量相互獨(dú)立的情況，若輸入變量之間存在相關(guān)性，可利用第1節(jié)的方法進(jìn)行處理。

3 基于最大熵的概率分布估計(jì)

任意隨機(jī)變量x的熵定義為：

其中，H（x）和 f（x）分別為隨機(jī)變量 x的熵和概率密度函數(shù)。

利用最大熵原理POME（Principle Of Maximum Entropy）估計(jì)隨機(jī)變量的概率密度函數(shù)的基本思路為：給定隨機(jī)變量x的相關(guān)統(tǒng)計(jì)信息，建立約束條件下的最大熵優(yōu)化模型，在候選概率分布解集Φ［f（x）］中，選擇使H（x）最大的 f（x）作為隨機(jī)變量 x的概率密度函數(shù)最優(yōu)估計(jì)。POME數(shù)學(xué)模型如下［14-15］：

其中，φ0（x）=1、φn（x）（n=1，2，…，N）為 N+1 個(gè)已知函數(shù)，稱(chēng)為基函數(shù)集為給定的N+1個(gè)隨機(jī)變量x的相關(guān)統(tǒng)計(jì)信息。POME模型的解可用φn（x）的函數(shù)進(jìn)行表達(dá)：

其中，λn（n=0，1，…，N）為拉格朗日參數(shù)。求解模型式（8）—（10）得出拉格朗日參數(shù)，即可獲得概率密度函數(shù) f（x）的最優(yōu)估計(jì)。

本文選擇冪函數(shù)形式 φn（x）=xn（n=0，1，…，N）作為POME模型的基函數(shù)集，則POME模型變?yōu)椋?/p>

4 含風(fēng)電系統(tǒng)的電壓穩(wěn)定概率分析方法

4.1 風(fēng)電功率模型

風(fēng)速的概率分布可用雙參數(shù)Weibull分布進(jìn)行描述。風(fēng)電場(chǎng)風(fēng)速vi的概率密度函數(shù)如式（13）所示：

其中，kw和c分別為Weibull分布的形狀參數(shù)和尺度參數(shù)。

風(fēng)電機(jī)組輸出功率與風(fēng)速的函數(shù)關(guān)系為：

其中，vci、vco和vR分別為切入風(fēng)速、切出風(fēng)速和額定風(fēng)速；PR為風(fēng)電機(jī)組額定功率。

4.2 負(fù)荷方向隨機(jī)模型

電壓穩(wěn)定分析中，通常選擇各節(jié)點(diǎn)負(fù)荷作為負(fù)荷增長(zhǎng)方向，以保證各節(jié)點(diǎn)負(fù)荷功率因數(shù)不變且按同比例增加。由于系統(tǒng)測(cè)量、估計(jì)等方面的誤差，實(shí)際的節(jié)點(diǎn)負(fù)荷并不是常數(shù)，存在一定的隨機(jī)性，可選擇正態(tài)分布反映其隨機(jī)特點(diǎn)。

假定節(jié)點(diǎn)i的有功負(fù)荷增長(zhǎng)方向IPLi滿足以基態(tài)有功PLi為均值、以σLi為標(biāo)準(zhǔn)差的正態(tài)分布，則其概率密度函數(shù)如下［8］：

對(duì)應(yīng)的無(wú)功負(fù)荷增長(zhǎng)方向IQLi由式（16）表示為：

其中，QLi為節(jié)點(diǎn)i的基態(tài)無(wú)功功率。

4.3 負(fù)荷裕度非線性?xún)?yōu)化模型

負(fù)荷裕度為運(yùn)行人員提供了系統(tǒng)從當(dāng)前運(yùn)行點(diǎn)到電壓崩潰點(diǎn)距離的直觀量度，是評(píng)估系統(tǒng)電壓穩(wěn)定性的最有效指標(biāo)。選擇負(fù)荷裕度為目標(biāo)函數(shù)，同時(shí)計(jì)及系統(tǒng)等式和不等式約束，構(gòu)建計(jì)算最大負(fù)荷裕度的非線性?xún)?yōu)化模型如下［16-19］：

其中，j?i表示節(jié)點(diǎn)j與節(jié)點(diǎn)i直接相連，包括j=i的情況；SB為系統(tǒng)節(jié)點(diǎn)集合；SG為發(fā)電機(jī)節(jié)點(diǎn)集合；SR為無(wú)功源集合；SL為支路集合；PGi和QRi分別為節(jié)點(diǎn)i的電源發(fā)出的有功和無(wú)功功率；PWi和QWi分別為節(jié)點(diǎn)i的風(fēng)電輸出有功、無(wú)功功率；PLi和QLi分別為節(jié)點(diǎn)i的有功負(fù)荷和無(wú)功負(fù)荷；ρ為負(fù)荷增長(zhǎng)系數(shù)，即負(fù)荷裕度；αij=δi- δj-φij；Ui和 δi分別為節(jié)點(diǎn) i的電壓幅值和相角；Yij和φij分別為節(jié)點(diǎn)導(dǎo)納矩陣對(duì)應(yīng)位置的導(dǎo)納元素幅值和相角；Sij為節(jié)點(diǎn)i與節(jié)點(diǎn)j之間的支路潮流大小；下劃線變量和上劃線變量分別對(duì)應(yīng)各自變量的下限值和上限值。式中等式約束由含負(fù)荷裕度參數(shù)的擴(kuò)展潮流方程組成，不等式約束由系統(tǒng)靜態(tài)安全約束組成。

4.4 含風(fēng)電系統(tǒng)電壓穩(wěn)定概率分析步驟

在負(fù)荷裕度優(yōu)化模型中加入風(fēng)電輸出功率和節(jié)點(diǎn)負(fù)荷的隨機(jī)特性，并采用概率方法進(jìn)行求解，即可實(shí)現(xiàn)含風(fēng)電系統(tǒng)電壓穩(wěn)定的概率分析評(píng)估。

本文假定風(fēng)電場(chǎng)風(fēng)速和節(jié)點(diǎn)負(fù)荷的概率分布參數(shù)已知，并給出相應(yīng)的風(fēng)速相關(guān)系數(shù)矩陣，則含風(fēng)電系統(tǒng)的電壓穩(wěn)定概率分析具體步驟如下。

a.確定輸入隨機(jī)變量個(gè)數(shù)n，在獨(dú)立標(biāo)準(zhǔn)正態(tài)空間S中對(duì)n個(gè)標(biāo)準(zhǔn)正態(tài)變量進(jìn)行采樣，并求解方程（7），形成基本采樣點(diǎn)集F1（S）和對(duì)應(yīng)權(quán)重系數(shù)集ω。

b.將Weibull分布參數(shù)和風(fēng)速相關(guān)系數(shù)矩陣代入式（2）—（4），求得相應(yīng)的變換矩陣L和多項(xiàng)式變換系數(shù) Am=［am1，am2，…，amn］T（m=0，1，2，3）。

c.將基本采樣點(diǎn)集 F1（S）代入式（5），得到輸入隨機(jī)變量的采樣點(diǎn)集F1（X）。

d.將F1（X）代入負(fù)荷裕度非線性?xún)?yōu)化模型進(jìn)行求解，得出對(duì)應(yīng)的負(fù)荷裕度離散點(diǎn)集F2（ρ）。非線性?xún)?yōu)化模型采用預(yù)測(cè)-校正原對(duì)偶內(nèi)點(diǎn)法進(jìn)行求解。

e.將 ω 和 F2（ρ）代入式（6），計(jì)算負(fù)荷裕度的各階矩信息。

f.將負(fù)荷裕度各階矩信息代入模型式（12）進(jìn)行求解，即可得到負(fù)荷裕度概率分布的最優(yōu)估計(jì)。

5 算例分析

5.1 系統(tǒng)概況

本文以IEEE 30和IEEE 118節(jié)點(diǎn)系統(tǒng)為基礎(chǔ)，對(duì)含風(fēng)電場(chǎng)的電力系統(tǒng)電壓穩(wěn)定進(jìn)行概率分析，程序?qū)崿F(xiàn)平臺(tái)為MATLAB R2009a，所用計(jì)算機(jī)的CPU主頻為2.1 GHz，內(nèi)存為2 GB。風(fēng)電場(chǎng)按恒功率因數(shù)方式運(yùn)行，各風(fēng)電場(chǎng)原始相關(guān)參數(shù)如表1所示，表中接入臺(tái)數(shù)列“/”前后數(shù)字分別表示118節(jié)點(diǎn)和30節(jié)點(diǎn)的標(biāo)準(zhǔn)接入臺(tái)數(shù)，風(fēng)電場(chǎng)原始相關(guān)系數(shù)矩陣RW為：

表1 風(fēng)電場(chǎng)相關(guān)參數(shù)Table1 Parameters of wind farms

IEEE測(cè)試系統(tǒng)數(shù)據(jù)參見(jiàn)文獻(xiàn)［20］，IEEE 30節(jié)點(diǎn)系統(tǒng)中，風(fēng)電場(chǎng)分別接入節(jié)點(diǎn)9、25、28，系統(tǒng)總負(fù)荷為189.2 MW；IEEE 118節(jié)點(diǎn)系統(tǒng)中，風(fēng)電場(chǎng)分別接入節(jié)點(diǎn) 23、39、114、117，系統(tǒng)總負(fù)荷為 3668 MW。

負(fù)荷增長(zhǎng)方式為全網(wǎng)負(fù)荷同時(shí)增加，各節(jié)點(diǎn)負(fù)荷按基態(tài)功率因數(shù)等比例增長(zhǎng)。各負(fù)荷分量服從以基態(tài)負(fù)荷為均值、標(biāo)準(zhǔn)差為5%的正態(tài)分布。

5.2 基于最大熵的概率分布估計(jì)結(jié)果

應(yīng)用基于最大熵的概率分布估計(jì)方法對(duì)算例系統(tǒng)進(jìn)行含風(fēng)電系統(tǒng)的電壓穩(wěn)定裕度概率計(jì)算，得出相應(yīng)的分布曲線，并同時(shí)與40000次蒙特卡洛模擬結(jié)果進(jìn)行對(duì)比分析。

圖1、2和圖3、4分別給出了IEEE 30和 118節(jié)點(diǎn)系統(tǒng)在最大熵估計(jì)和蒙特卡洛模擬2種方法下對(duì)應(yīng)的負(fù)荷裕度概率密度分布曲線和累計(jì)概率分布曲線，圖中負(fù)荷裕度為標(biāo)幺值，后同。從圖中可以看出，由最大熵估計(jì)方法得出的分布曲線較為合理地?cái)M合了蒙特卡洛仿真結(jié)果。對(duì)比圖1和圖3中由最大熵估計(jì)方法得出的概率密度曲線可知，圖1中的密度曲線發(fā)生了畸變，存在“翹尾現(xiàn)象”，而圖3中的密度曲線變形較小，僅在其尾部有一定程度的收縮。結(jié)合5.1節(jié)的算例數(shù)據(jù)可知，IEEE 30節(jié)點(diǎn)系統(tǒng)的風(fēng)電接入比例大于IEEE 118節(jié)點(diǎn)系統(tǒng)，這表明非正態(tài)變量比例的增加會(huì)對(duì)最大熵概率分布估計(jì)結(jié)果產(chǎn)生影響。觀察圖2和圖4可知，2種方法所得累計(jì)概率分布曲線較為接近，較好地反映了負(fù)荷裕度的累積分布特性。

表2給出了服從Weibull分布的風(fēng)電場(chǎng)風(fēng)速隨機(jī)變量多項(xiàng)式正態(tài)變換系數(shù)，相應(yīng)的變換矩陣L為：

表3給出了基于最大熵的概率分布估計(jì)參數(shù)，可用對(duì)應(yīng)的估計(jì)函數(shù)描述負(fù)荷裕度的概率密度分布。表4給出了2種方法下不同算例系統(tǒng)的負(fù)荷裕度均值（標(biāo)幺值，后同）、標(biāo)準(zhǔn)差（標(biāo)幺值，后同）以及2種方法計(jì)算結(jié)果的相對(duì)誤差。由表4數(shù)據(jù)可得，均值的相對(duì)誤差都小于1%，標(biāo)準(zhǔn)差的相對(duì)誤差都小于4%，計(jì)算精度令人滿意，在工程應(yīng)用的誤差要求范圍之內(nèi)。表5給出了2種方法下的計(jì)算仿真規(guī)模和計(jì)算時(shí)間。由表5可得，最大熵估計(jì)方法的計(jì)算時(shí)間僅為蒙特卡洛方法的0.14%和0.6%，大幅節(jié)省了計(jì)算時(shí)間且能獲得較為滿意的計(jì)算結(jié)果，為電壓穩(wěn)定概率分析方法的工程實(shí)際應(yīng)用實(shí)現(xiàn)提供了可能性。

圖1 IEEE 30節(jié)點(diǎn)系統(tǒng)負(fù)荷裕度概率密度分布Fig.1 Probabilistic density distribution of load margin of IEEE 30-bus system

圖2 IEEE 30節(jié)點(diǎn)系統(tǒng)負(fù)荷裕度累計(jì)概率分布Fig.2 Cumulative probabilistic density distribution of load margin of IEEE 30-bus system

圖3 IEEE 118節(jié)點(diǎn)系統(tǒng)負(fù)荷裕度概率密度分布Fig.3 Probabilistic density distribution of load margin of IEEE 118-bus system

圖4 IEEE 118節(jié)點(diǎn)系統(tǒng)負(fù)荷裕度累計(jì)概率分布Fig.4 Cumulative probabilistic density distribution of load margin of IEEE 118-bus system

5.3 風(fēng)電參數(shù)影響分析

為了研究風(fēng)速概率密度函數(shù)變化和風(fēng)速相關(guān)系數(shù)矩陣變化對(duì)電壓穩(wěn)定裕度的影響，本節(jié)假設(shè)在不同風(fēng)速分布參數(shù)和相關(guān)系數(shù)矩陣組合下對(duì)電壓穩(wěn)定問(wèn)題進(jìn)行概率求解和結(jié)果分析，得出風(fēng)電參數(shù)變化對(duì)電壓穩(wěn)定裕度概率分布所帶來(lái)的影響。

表2 風(fēng)速隨機(jī)變量多項(xiàng)式正態(tài)變換系數(shù)Table 2 Polynomial normal transformation coefficients of wind speed random variables

表3 基于最大熵原理的概率分布估計(jì)參數(shù)Table 3 Parameters of probability distribution estimation based on maximum entropy

表4 兩算例系統(tǒng)計(jì)算結(jié)果及誤差分析Table 4 Calculative results and error analysis for two case systems

表5 兩算例系統(tǒng)仿真規(guī)模和計(jì)算時(shí)間Table 5 Simulation scale and calculation time of two case systems

本節(jié)以IEEE 118節(jié)點(diǎn)系統(tǒng)為研究對(duì)象，分別以表1的風(fēng)電場(chǎng)原始風(fēng)速參數(shù)和原始相關(guān)系數(shù)矩陣RW為標(biāo)準(zhǔn)，同方向增加或減小風(fēng)速分布參數(shù)數(shù)值和相關(guān)系數(shù)矩陣中的元素?cái)?shù)值，分別形成以下3個(gè)新風(fēng)速參數(shù)場(chǎng)景和3個(gè)新相關(guān)系數(shù)矩陣場(chǎng)景：風(fēng)速參數(shù)場(chǎng)景取名為“較小風(fēng)速場(chǎng)景”、“較大風(fēng)速場(chǎng)景1”和“較大風(fēng)速場(chǎng)景2”，各場(chǎng)景的具體參數(shù)如表6所示；相關(guān)系數(shù)矩陣場(chǎng)景取名為“零相關(guān)系數(shù)矩陣”、“較小相關(guān)系數(shù)矩陣”和“較大相關(guān)系數(shù)矩陣”，其中“零相關(guān)系數(shù)矩陣”用4階單位矩陣I4×4描述，“較小相關(guān)系數(shù)矩陣”和“較大相關(guān)系數(shù)矩陣”分別用RWS和RWL描述，各系數(shù)矩陣參數(shù)如下：

表6 不同場(chǎng)景下的風(fēng)速分布參數(shù)Table 6 Parameters of wind speed distribution for different cases

圖5為不同風(fēng)速參數(shù)場(chǎng)景下計(jì)算得出的負(fù)荷裕度概率密度曲線。觀察圖中曲線，可以看出風(fēng)速參數(shù)對(duì)概率密度分布的影響：隨著風(fēng)速參數(shù)的逐漸增大，概率密度曲線朝著扁平化的趨勢(shì)發(fā)展，曲線頂點(diǎn)對(duì)應(yīng)的概率密度數(shù)值逐漸減小，負(fù)荷裕度數(shù)值逐漸增大。這說(shuō)明風(fēng)速參數(shù)的增加有助于系統(tǒng)負(fù)荷裕度的提高，但隨著風(fēng)速參數(shù)的增加，負(fù)荷裕度隨機(jī)變量的變化也逐漸變大，風(fēng)速變化所帶來(lái)的波動(dòng)性更加明顯。

圖5 不同風(fēng)速下的負(fù)荷裕度概率密度曲線Fig.5 Probabilistic density curve of load margin for different wind speeds

表7所示為不同風(fēng)速參數(shù)和風(fēng)速相關(guān)系數(shù)矩陣組合下的負(fù)荷裕度計(jì)算結(jié)果，通過(guò)對(duì)比表中數(shù)據(jù)，驗(yàn)證了圖5所得相關(guān)結(jié)論的正確性。由上述分析可知，在選擇風(fēng)能資源時(shí)，并不能只考慮風(fēng)資源本身的變化因素，還要考慮所接入系統(tǒng)的電壓穩(wěn)定承受能力等影響因素，這樣才能將風(fēng)能資源和所接系統(tǒng)進(jìn)行合理的匹配。

表7中還給出風(fēng)速相關(guān)系數(shù)矩陣變化時(shí)的負(fù)荷裕度模型計(jì)算結(jié)果。觀察表中數(shù)據(jù)可知，相關(guān)系數(shù)變化會(huì)對(duì)負(fù)荷裕度計(jì)算結(jié)果產(chǎn)生影響：隨著風(fēng)速相關(guān)性增加，有些負(fù)荷裕度均值逐漸減小，有些負(fù)荷裕度均值先減后增，相關(guān)系數(shù)的變化引起負(fù)荷裕度均值的變化；隨著相關(guān)系數(shù)增大，負(fù)荷裕度標(biāo)準(zhǔn)差逐漸增大，表明其波動(dòng)更加明顯。因此在負(fù)荷裕度概率分析時(shí)，不能忽略風(fēng)速相關(guān)性，必須加以考慮，否則分析結(jié)果與實(shí)際情況偏差較大，不能正確反映負(fù)荷裕度變化的概率特性。本文方法所得結(jié)果和蒙特卡洛模擬所得結(jié)果的相對(duì)誤差如表7所示。從表中可以看出，采用本文方法所得的負(fù)荷裕度結(jié)果與蒙特卡洛模擬所得結(jié)果較為接近，結(jié)果均值誤差較小，標(biāo)準(zhǔn)差誤差偏大，但均在可接受范圍之內(nèi)，可認(rèn)為本文方法所得結(jié)果正確有效。綜上可知，采用基于多項(xiàng)式正態(tài)變換的點(diǎn)估計(jì)法能較好地處理風(fēng)速相關(guān)性問(wèn)題，計(jì)算結(jié)果反映了負(fù)荷裕度變化趨勢(shì)，相對(duì)蒙特卡洛模擬方法，計(jì)算時(shí)間縮短為原來(lái)的1%以下，適合于對(duì)電力系統(tǒng)概率問(wèn)題的分析求解。

5.4 風(fēng)電場(chǎng)接入容量影響分析

隨著風(fēng)電場(chǎng)的發(fā)展，風(fēng)電容量的接入比例會(huì)越來(lái)越大，本節(jié)就不同的風(fēng)電接入比例對(duì)負(fù)荷裕度概率模型進(jìn)行求解，觀察容量變化帶來(lái)的影響。標(biāo)準(zhǔn)風(fēng)電接入容量為230MW，考慮不同的標(biāo)準(zhǔn)接入容量倍數(shù)，得到不同風(fēng)電接入比例下的負(fù)荷裕度均值和標(biāo)準(zhǔn)差，具體數(shù)據(jù)如表8所示。可以看出，風(fēng)電比例越高，負(fù)荷裕度的均值和標(biāo)準(zhǔn)差均有所增加，這是因?yàn)橐环矫骘L(fēng)電的接入提供了電源支持，使得負(fù)荷裕度水平有所提高，但是由于風(fēng)電固有的隨機(jī)特性，容量的增加也讓波動(dòng)更加明顯，標(biāo)準(zhǔn)差隨之增大。因此，如何有效地減小風(fēng)電波動(dòng)性帶來(lái)的影響是大規(guī)模利用風(fēng)電必須考慮的問(wèn)題，否則大量的風(fēng)電接入會(huì)給系統(tǒng)運(yùn)行的穩(wěn)定性造成不利的影響。

表7 不同風(fēng)電參數(shù)組合下的負(fù)荷裕度模型計(jì)算結(jié)果Table 7 Calculative results of load margin for different wind power parameter combinations

表8 不同風(fēng)電接入容量下的負(fù)荷裕度模型計(jì)算結(jié)果Table 8 Calculative results of load margin for different wind power penetration capacities

5.5 風(fēng)電場(chǎng)功率因數(shù)影響分析

風(fēng)電場(chǎng)發(fā)電控制通常采用恒功率策略，不同的功率因數(shù)設(shè)定值對(duì)負(fù)荷裕度有所影響。本節(jié)采用最大熵估計(jì)法快速求解不同功率因數(shù)下的負(fù)荷裕度模型，計(jì)算結(jié)果如表9所示。表中負(fù)功率因數(shù)表示風(fēng)電場(chǎng)需要消耗無(wú)功。

表9 不同功率因數(shù)下的負(fù)荷裕度計(jì)算結(jié)果Table 9 Calculative results of load margin for different power factors

由表9可以看出，隨著功率因數(shù)的增加，負(fù)荷裕度均值逐漸增加，標(biāo)準(zhǔn)差逐漸減小。這是因?yàn)殡S著功率因數(shù)的增加，整個(gè)風(fēng)電場(chǎng)的無(wú)功消耗逐漸降低，無(wú)功輸出逐漸增加。由于風(fēng)電場(chǎng)能夠向系統(tǒng)提供無(wú)功，因此系統(tǒng)的電壓穩(wěn)定情況得以改善，電壓穩(wěn)定裕度有所增加，能夠承受更大的系統(tǒng)變化，在同樣的有功波動(dòng)情況下，能使電壓穩(wěn)定裕度的均值更大，而標(biāo)準(zhǔn)差更小，電壓穩(wěn)定特性更加優(yōu)良。表中還給出了不同相關(guān)系數(shù)矩陣條件下，功率因數(shù)變化時(shí)的負(fù)荷裕度計(jì)算結(jié)果，不同條件下的計(jì)算結(jié)果變化趨勢(shì)一致，驗(yàn)證了前述的原因分析。

6 結(jié)論

本文提出了一種基于多項(xiàng)式正態(tài)變換和最大熵估計(jì)的含風(fēng)電系統(tǒng)電壓穩(wěn)定概率分析方法。該方法結(jié)合多項(xiàng)式正態(tài)變換、點(diǎn)估計(jì)法以及最大熵估計(jì)，能夠?qū)紤]相關(guān)性的含風(fēng)電系統(tǒng)電壓穩(wěn)定問(wèn)題進(jìn)行完整的概率分析，得出全面的概率信息。本文將風(fēng)電輸入功率和節(jié)點(diǎn)負(fù)荷功率視為隨機(jī)變量，以蒙特卡洛模擬作為評(píng)價(jià)標(biāo)準(zhǔn)，驗(yàn)證了本文算法的計(jì)算精度，討論了不同風(fēng)速參數(shù)、相關(guān)系數(shù)矩陣變化、不同風(fēng)電接入容量以及風(fēng)電場(chǎng)功率因數(shù)變化對(duì)電壓穩(wěn)定裕度的影響。算例分析結(jié)果表明，基于多項(xiàng)式正態(tài)變換和最大熵估計(jì)的電壓穩(wěn)定概率分析方法能夠獲得正確的負(fù)荷裕度概率分布，且與蒙特卡洛模擬相比，具有更高的計(jì)算效率，能夠滿足負(fù)荷裕度概率信息及時(shí)獲取的要求，為運(yùn)行人員控制決策提供了重要參考。