計及隨機變量相關(guān)性的多點線性化概率潮流計算

葉 希，王彥灃，黃 楊，廖幫昆，歐陽雪彤，文云峰

（1.國網(wǎng)四川省電力公司，成都 610041；2.輸變電新技術(shù)教育部工程研究中心（湖南大學(xué)），長沙 410082）

在“碳達峰/碳中和”發(fā)展目標和新型電力系統(tǒng)構(gòu)建的背景下，風(fēng)電、光伏等新能源發(fā)電裝機容量快速增長，源?網(wǎng)?荷各環(huán)節(jié)電力電子化程度顯著加深。以同步電源為主的電力系統(tǒng)正逐漸演進為“高比例新能源”并網(wǎng)和“高比例電力電子裝備”饋入的“雙高”電力系統(tǒng)[1]。在“雙高”電力系統(tǒng)中，源荷雙側(cè)的功率不確定性顯著增強，電力電量平衡以及容量充裕度的概念與方法將由目前確定性的思路向概率性思路轉(zhuǎn)變[2]。

概率潮流計算能夠計及不確定性因素并分析其對系統(tǒng)運行特性的影響，可為“雙高”電力系統(tǒng)安全穩(wěn)定分析提供更為可靠的理論支撐[3?4]。自1974年Borkowska[5]首次提出概率潮流概念后，經(jīng)過國內(nèi)外學(xué)者近50年的大量研究，形成了多種計算方法，大致可以歸納為3 類：模擬法[6?7]、解析法[8?10]和近似法[11?13]。

模擬法以蒙特卡洛模擬法MCSM（Monte Carlo simulation method）為代表，通過簡單隨機抽樣SRS（simple random sampling）法對輸入隨機變量進行抽樣得到樣本，并對樣本中每組數(shù)據(jù)進行獨立的確定性潮流計算，最后統(tǒng)計輸出結(jié)果得到各節(jié)點狀態(tài)與支路潮流的概率分布。該方法的優(yōu)點是簡單易行，適用于任意節(jié)點規(guī)模的系統(tǒng)，且在樣本規(guī)模足夠大的情況下具有較高的精度；缺點是當(dāng)樣本規(guī)模以及節(jié)點規(guī)模較大時，大量的潮流計算相當(dāng)耗時，效率低下，難以應(yīng)用于快速分析或?qū)崟r監(jiān)測評估領(lǐng)域，常用于評判其他計算方法的優(yōu)劣。

為了提升概率潮流的計算效率，國內(nèi)外學(xué)者進行了大量研究。例如：解析法中的半不變量法[10]，即結(jié)合了數(shù)值特征分析手段的概率潮流分析方法；文獻[14]采用線性化的潮流方程并與半不變量法相結(jié)合，降低了計算時間；文獻[15]提出了一種利用最大熵原理改進的C 型Gram?Charlier 即CGC（C?type Gram?Charlier）級數(shù)展開法，提升了傳統(tǒng)CGC 方法的計算效率，但該方法中尚未考慮輸入隨機變量之間的相關(guān)性，且采用Weibull 分布描述風(fēng)電場出力概率模型，對于分布特征較為復(fù)雜的概率模型可能存在擬合效果較差的情況。點估計法是近似法中的代表，其基本原理與MCSM 相同，但樣本規(guī)模相對于MCSM 更少，結(jié)合數(shù)值統(tǒng)計方法可使效率更高。文獻[16]提出了一種可用于計算含服從非正態(tài)分布負荷的概率潮流問題的改進型兩點估計法，其精度高于傳統(tǒng)的兩點估計法，但點估計法的缺點是當(dāng)系統(tǒng)節(jié)點規(guī)模十分龐大時，潮流計算次數(shù)也隨之增加，可能出現(xiàn)計算效率較低的情況。

由于MCSM的計算精度較高且原理較為簡單，近年來國內(nèi)外學(xué)者針對其抽樣技術(shù)以及潮流求解方法進行優(yōu)化改進，以解決計算耗時長的問題。文獻[6]利用最大似然估計建立混合高斯模型，在新能源機組出力擾動及負荷需求靈活多變的情況下對其概率分布進行精確描述，建立源荷雙側(cè)統(tǒng)一概率模型；文獻[17]提出了擴展拉丁超立方采樣技術(shù)，改進了MCSM的抽樣手段，在減少抽樣規(guī)模的同時保證了采樣覆蓋效果；文獻[18]基于Copula理論構(gòu)建多維隨機變量聯(lián)合概率分布以獲得具有相關(guān)性的數(shù)據(jù)樣本，但其仍基于傳統(tǒng)蒙特卡羅模擬法進行概率潮流計算，未能進一步有效地提高潮流計算效率。

針對蒙特卡洛模擬法計算概率潮流效率低下以及輸入隨機變量概率分布與相關(guān)性建模不足的問題，本文提出一種計及隨機變量相關(guān)性的多點線性化概率潮流計算方法。對于有實際監(jiān)測數(shù)據(jù)的新能源場站或負荷，采用核密度估計法建立其概率分布模型，而無監(jiān)測數(shù)據(jù)時則采用參數(shù)估計法獲得其概率分布。通過構(gòu)建Copula 函數(shù)描述輸入隨機變量間的相關(guān)性，獲取具有相關(guān)性的隨機變量樣本。在MCSM 基礎(chǔ)上將非線性潮流方程轉(zhuǎn)換為線性化方程，并采用多點線性化技術(shù)減少線性化帶來的截斷誤差，以適應(yīng)“雙高”電力系統(tǒng)中考慮強不確定性因素的計算需求。在IEEE?30 節(jié)點系統(tǒng)上進行算例分析，驗證所提方法的準確性和有效性。

1 源荷雙側(cè)概率分布模型的建立

由于新能源出力不確定性以及負荷隨機擾動性，“雙高”電力系統(tǒng)中節(jié)點注入功率具有較強的時變特征，難以應(yīng)用單一的標準分布模型準確擬合其概率密度函數(shù)PDF（probability density function）以及累積分布函數(shù)CDF（cumulative distribution func?tion）[10]，而源荷雙側(cè)概率分布模型的擬合精確程度直接影響到概率潮流計算結(jié)果的精確性，因此需要運用更加準確的擬合方法對源荷概率模型進行刻畫。核密度估計KDE（kernel density estimation）[19]是一種非參數(shù)估計方法，完全根據(jù)數(shù)據(jù)本身的特點和性質(zhì)來擬合分布，能夠比參數(shù)估計法得到更好的模型。本文將它應(yīng)用于源荷雙側(cè)有功概率分布模型的建立。

設(shè)P1，P2，…，Pn是取自某一新能源機組輸出功率P的n個樣本，則該機組所在節(jié)點的輸出功率的概率密度函數(shù)f（P）的核密度估計定義為

核函數(shù)通常選用以0 為中心的單峰對稱概率密度函數(shù)，結(jié)合大部分常規(guī)風(fēng)電、光伏機組輸出功率分布特征，可選擇高斯核函數(shù)，即

基于高斯核函數(shù)的核密度估計函數(shù)為

當(dāng)核函數(shù)K定義在[?1，1]上對稱且滿足核函數(shù)約束條件[19]，則當(dāng)h→0、nh→+∞時，有

特別地，當(dāng)使用高斯核函數(shù)進行核密度估計時，則最優(yōu)窗寬為

圖1 風(fēng)電場實際輸出功率的核密度估計Fig.1 Kernel density estimation of real power output from wind farm

2 輸入隨機變量間的相關(guān)性處理

大部分概率潮流計算方法一般將源荷雙側(cè)隨機變量視為相互獨立[15]。然而實際系統(tǒng)運行中，由于氣候、地理等因素影響，新能源機組之間或負荷之間存在著一定相關(guān)性。若忽略這些輸入隨機變量之間的相關(guān)性，將會給概率潮流計算結(jié)果帶來誤差。Copula理論[18]能夠把隨機變量之間的相關(guān)關(guān)系與變量的邊際分布分開進行研究，不僅能完整描述輸入隨機變量間的相互關(guān)系，而且不受輸入隨機變量的邊緣分布類型及是否線性相關(guān)的限制，導(dǎo)出的相關(guān)性度量在線性變換和嚴格單增變換下都不會改變。為此，本文采用Copula理論描述源荷雙側(cè)輸入隨機變量之間的相關(guān)性。

本文采用橢球Copula 函數(shù)中的正態(tài)Copula 函數(shù)描述多維隨機變量間的相關(guān)性[18]。以構(gòu)建二維高斯Copula函數(shù)為例，相關(guān)系數(shù)為ρvi,vj的隨機變量vi與vj的Copula函數(shù)C的表達式為

式中：Φ?1和Φρvi,vj分別為標準正態(tài)分布的累積分布逆函數(shù)以及二維正態(tài)聯(lián)合累積分布函數(shù)；ui與uj為隨機變量vi與vj的累積量。多維隨機變量間的Copula 函數(shù)可通過參數(shù)辨識求解或結(jié)合各隨機變量邊緣累積分布函數(shù)進行直接求解，詳細求解過程可參考文獻[18]。

步驟3根據(jù)等概率轉(zhuǎn)換原則，生成滿足給定相關(guān)性條件的m個具有相關(guān)性的隨機變量數(shù)據(jù)樣本Vm×N，其等概率轉(zhuǎn)換公式為

3 計及隨機變量相關(guān)性的多點線性化概率潮流計算

基于蒙特卡洛模擬的傳統(tǒng)概率潮流計算方法，將規(guī)模龐大的數(shù)據(jù)樣本分別進行多次獨立潮流計算，存在計算耗時、效率低下的難題，無法有效應(yīng)用于實際大規(guī)模電力系統(tǒng)[7]。為此，本文提出了一種計及隨機變量相關(guān)性的多點線性化概率潮流計算方法，通過引入潮流方程線性化處理技術(shù)，顯著降低了蒙特卡洛模擬過程中的非線性潮流計算量，避免了單點線性化潮流計算無法適應(yīng)節(jié)點注入功率過于偏離基準運行點的情況。對于“雙高”電力系統(tǒng)新能源出力高度不確定性和負荷側(cè)功率日內(nèi)大幅度變化的情況，該方法能夠在兼顧計算精度的同時大幅提升計算效率。

3.1 潮流方程的線性化

牛頓?拉夫遜法潮流方程的節(jié)點功率方程以及支路潮流方程簡化為

式中：W為由各節(jié)點注入功率構(gòu)成的矩陣；X為節(jié)點狀態(tài)矩陣；Z為支路潮流功率矩陣；f(X)和g(X)分別為節(jié)點電壓方程與支路潮流方程。

將節(jié)點注入功率W視為系統(tǒng)的輸入隨機變量，則潮流計算結(jié)果X與Z可視為輸出隨機變量，三者均可寫為基準值+隨機擾動值的形式，即

式中：W0、X0和Z0分別為節(jié)點注入功率、節(jié)點狀態(tài)和支路潮流的基準值；ΔW、ΔX和ΔZ為節(jié)點注入功率、節(jié)點狀態(tài)和支路潮流的隨機擾動值。

將潮流方程在系統(tǒng)運行的基準點進行泰勒展開，并忽略高次項，則可得潮流方程的線性化表達式[14]為

3.2 潮流方程的多點線性化

上述潮流方程的線性化處理屬于單點線性化，在系統(tǒng)基準運行點進行泰勒展開，在各節(jié)點的注入功率擾動范圍不大時，具有較好的精度。但高比例新能源并網(wǎng)后源側(cè)出力不確定性增強，使得并網(wǎng)節(jié)點的注入功率處于大幅度變化狀態(tài)，新能源輸出功率偏離給定期望值較大時，會導(dǎo)致單點泰勒展開的潮流方程與實際方程產(chǎn)生較大偏離[14]。此時單點線性化的潮流方程已無法滿足計算精度的需求，需要對其進行優(yōu)化改進。

此時，若節(jié)點i風(fēng)電場輸出功率的某一樣本數(shù)據(jù)的擾動處于第k1個區(qū)域（k1=1，2，…，n）之內(nèi)，且其輸出功率全體樣本數(shù)據(jù)的均值為Pimean，處于第k2個區(qū)域（k2=1，2，…，n）之內(nèi)，則可得到由該節(jié)點注入功率擾動引起的節(jié)點狀態(tài)與支路潮流擾動量的表達式為

應(yīng)指出，若系統(tǒng)中同時存在多個需要進行多點線性化處理的節(jié)點，在對其中單個節(jié)點進行線性化處理時，其余節(jié)點注入功率按各自的基準值處理。設(shè)系統(tǒng)中節(jié)點i與j需要進行多點線性化處理，則先對節(jié)點i注入功率進行分段劃分，結(jié)合式（13）得到節(jié)點i引起的ΔXi與ΔZi，此時節(jié)點j處注入功率取其基準值；而后對節(jié)點j進行處理時，節(jié)點i處注入功率取其基準值，得到節(jié)點i引起的ΔXj與ΔZj。將上述節(jié)點產(chǎn)生的擾動量疊加即可得到多點線性化處理下節(jié)點i、j注入功率對系統(tǒng)各節(jié)點狀態(tài)和各支路潮流的擾動量。綜上所述，多個節(jié)點的多點線性化處理原則為：各節(jié)點的線性化處理過程相互獨立且依次進行，在針對單個節(jié)點的處理過程中不同時考慮其他節(jié)點的擾動性。

當(dāng)系統(tǒng)中存在較多的功率處于大范圍擾動的節(jié)點，則在多點線性化處理過程中需要進行多次潮流計算，將會影響計算效率。針對該情況，可適當(dāng)擴大各節(jié)點的注入功率劃分的區(qū)間寬度，通過減少分段數(shù)量來減少潮流計算次數(shù)，且對于系統(tǒng)中功率擾動范圍相對較小的節(jié)點，可進一步減少其節(jié)點注入功率的分段數(shù)量，或?qū)ζ洳捎脝吸c線性化處理。顯而易見，上述方法將提升計算效率，但計算誤差略有增大。

3.3 算法流程

本文所提出的概率潮流計算方法，主旨是對蒙特卡洛模擬法中的大量非線性潮流計算進行改進并轉(zhuǎn)化為線性化潮流計算，進一步針對單點線性化計算所存在的誤差問題進行優(yōu)化，利用多點線性化處理降低線性化誤差。該算法主要流程如圖2所示。

步驟1源荷雙側(cè)概率分布模型的建立。對于有監(jiān)測數(shù)據(jù)的新能源場站或負荷，通過核密度估計得到各節(jié)點注入功率的累積分布函數(shù)F及逆函數(shù)F?1；對于無監(jiān)測數(shù)據(jù)的新能源場站或負荷，則采用Weibull分布或正態(tài)分布等概率分布模型。

步驟3統(tǒng)計各節(jié)點注入功率數(shù)據(jù)，確定系統(tǒng)節(jié)點注入功率基準值W0，基于基準值進行潮流計算，得到節(jié)點狀態(tài)基準X0、支路潮流基準Z0以及雅可比矩陣J0與靈敏度矩陣S0。對功率擾動范圍大的節(jié)點進行區(qū)域劃分，求解各區(qū)域?qū)?yīng)的靈敏度矩陣Sk。式（11）中ΔW由上一步獲得的樣本Vn×N減去W0得到。

步驟4對于功率擾動范圍大的節(jié)點，先確定其節(jié)點注入功率數(shù)據(jù)所在區(qū)域，結(jié)合式（13）進行多點線性化計算；其余節(jié)點則結(jié)合式（12）進行單點線性化計算。兩者疊加得到輸出隨機變量擾動值ΔXn×N與ΔZn×N。

步驟5將X0、Z0與ΔXn×N、ΔZn×N按式（11）進行疊加，得到輸出隨機變量數(shù)據(jù)Xn×N與Zn×N，最終獲得節(jié)點電壓以及支路潮流的概率分布。

4 算例分析

4.1 測試系統(tǒng)

基于Matlab 平臺搭建IEEE?30 節(jié)點系統(tǒng)，對所提出的概率潮流計算方法可行性與有效性進行驗證。測試系統(tǒng)的網(wǎng)絡(luò)結(jié)構(gòu)以及線路參數(shù)等，詳見文獻[20]。在節(jié)點6、22、27、28 處各接入30 臺額定出力為2.5 MW 的風(fēng)電機組，監(jiān)測數(shù)據(jù)來自中國東南沿海某4 個風(fēng)電場2020 年1 月1 日—12 月31 日的8 784組有功出力數(shù)據(jù)。風(fēng)電機組控制手段為恒功率因數(shù)控制。假設(shè)測試系統(tǒng)中各負荷服從正態(tài)分布，其期望值為原始系統(tǒng)數(shù)據(jù)的負荷基準值，標準差為期望值的30%。對測試系統(tǒng)接入的風(fēng)電場輸出功率概率模型進行核密度估計擬合。以節(jié)點6與節(jié)點22 風(fēng)電場輸出功率PW1與PW2概率分布為例，其PDF擬合結(jié)果如圖3所示。

圖3 核密度估計與混合高斯模型PDF 曲線比較Fig.3 Comparison of PDF curve between kernel density estimation and Gaussian mixture model

圖3 中，將監(jiān)測數(shù)據(jù)實際概率密度直方圖與核密度估計所得的概率密度函數(shù)以及由最大似然估計法確定的混合高斯模型三者進行對比，可以看出，核密度估計法更好地刻畫了數(shù)據(jù)真實的分布情況，相比于混合高斯模型更為精確；且通過改變窗寬可任意調(diào)節(jié)擬合曲線的平滑度，相對于混合高斯模型有更好的泛用性。

經(jīng)計算得到PW1與PW2的監(jiān)測數(shù)據(jù)相關(guān)系數(shù)為0.817 4。對2 個風(fēng)電場的監(jiān)測數(shù)據(jù)進行等距抽樣獲得800 組離散數(shù)據(jù)，繪制出數(shù)據(jù)直方圖與散點圖，如圖4所示。圖中柱狀圖表示風(fēng)電場實際輸出功率的概率密度直方圖。分析圖4 中散點分布情況可看出，2個風(fēng)電場監(jiān)測數(shù)據(jù)之間具有相關(guān)性，且屬于高度相關(guān)。

圖4 風(fēng)電場實際輸出功率散點圖與直方圖Fig.4 Scatter plot and histogram of real power output from wind farm

結(jié)合抽樣技術(shù)獲得節(jié)點6與節(jié)點22的1 000組滿足相關(guān)性條件的風(fēng)電場輸出功率樣本數(shù)據(jù)PWS1與PWS2，具體數(shù)據(jù)分布情況如圖5所示。

圖5 風(fēng)電場模擬輸出功率散點圖與直方圖Fig.5 Scatter plot and histogram of simulated power output from wind farm

圖5 中柱狀圖表示風(fēng)電場實際輸出功率的概率密度直方圖。分析圖5中數(shù)據(jù)分布并與圖4進行對比可看出，基于Copula 理論所獲得的節(jié)點6與節(jié)點22 風(fēng)電場樣本數(shù)據(jù)間也具有與原始數(shù)據(jù)基本一致的相關(guān)性。進一步對樣本數(shù)據(jù)進行統(tǒng)計，得到其相關(guān)系數(shù)為0.838 8，與實際監(jiān)測數(shù)據(jù)的相關(guān)系數(shù)接近，驗證了基于Copula理論生成具有相關(guān)性隨機變量樣本數(shù)據(jù)方法的有效性。

4.2 仿真計算與誤差分析

采用本文所提出的計算方法對測試系統(tǒng)進行仿真計算。首先對各節(jié)點風(fēng)電場功率進行分段劃分。設(shè)定每段功率區(qū)間寬度為5 MW，共劃分15個區(qū)域，基準容量取100 MV·A，功率區(qū)間上、下限取標幺值，分別為[0，0.05），[0.05，0.10），…，[0.70，0.75]。以60 000次MCSM計算結(jié)果作為參考值，進行結(jié)果比較。

以節(jié)點15 電壓幅值、支路24?25 有功功率、支路6?28無功功率累積分布函數(shù)為例，2種方法所得曲線對比結(jié)果如圖6所示，圖中縱坐標RCDF表示概率的累積量。分析圖6中曲線可知，本文所提出的方法與MCSM具有十分接近的計算精度，兩者所得到的曲線結(jié)果十分接近。

圖6 節(jié)點電壓幅值或支路潮流CDF 曲線比較Fig.6 Comparison of CDF curve of node voltage magnitude or branch flow power

圖7所示為將MCSM、單點線性化（方法1）和多點線性化（方法2）3 種方法所求得的節(jié)點3 電壓PDF 與CDF 曲線同時進行對比的曲線。由圖7 可知，單點線性化方法由于僅在系統(tǒng)運行基準點進行泰勒展開線性化處理，當(dāng)輸入隨機變量相對于基準點偏離程度增大時，所求得曲線偏離程度也隨之增大。而得益于對非線性潮流方程的多點泰勒展開，本文所提方法求得的曲線與MCSM結(jié)果基本吻合，有效地降低了由于新能源機組出力擾動范圍過大而引發(fā)的截斷誤差，提升了計算精度。

圖7 不同方法下節(jié)點電壓幅值CDF 與PDF 曲線比較Fig.7 Comparison of CDF and PDF curves of node voltage magnitude for different methods

進一步，采用方差和的根均值A(chǔ)RMS（average root mean square）[21]來衡量不同方法下概率潮流計算的誤差。ARMS定義為

式中：CMCi與Ci分別為MCSM與待分析計算方法所求累積分布函數(shù)上第i個數(shù)據(jù)點的值；N為數(shù)據(jù)規(guī)模。ARMS的值越小說明計算結(jié)果誤差越小。表1中列舉了部分節(jié)點電壓以及支路潮流的ARMS，并將方法1 與方法2 的精確度進行了比較。分析表1數(shù)據(jù)可見，所提出的多點線性化計算方法在高比例新能源并網(wǎng)節(jié)點功率大幅度變化的場景下，仍能保持較高的準確度，如：使用方法2得到的節(jié)點3電壓幅值U3的ARMS 為0.16，相較方法1 的1.82 提升了91%；節(jié)點14、19 電壓幅值與支路12?25、6?28 有功功率等計算精度提升幅度均在80%以上；多個節(jié)點電壓以及支路功率的ARMS 均能控制在0.5%以內(nèi)?？梢姡疚乃岱椒ㄏ啾扔趩吸c線性化方法，精度有大幅度提升。

表1 部分節(jié)點電壓幅值或支路功率的ARMSTab.1 ARMS of part node voltage magnitude or branch power

將60 000 次MCSM 所需時間與相同抽樣規(guī)模的方法1 和方法2 的進行對比，結(jié)果如表2 所示。分析表2 中數(shù)據(jù)可知，相對于MCSM 的計算時間684.83 s，本文所提方法的計算所需時間縮短至3.20 s，僅為MCSM 的0.4%，有效地提升了計算效率。相對于單點線性化的計算方法，本文所提方法在高比例新能源出力強不確定性場景下，能更有效地降低線性化所帶來的誤差，而代價則是在多點線性化處理過程中需要額外進行少量非線性潮流計算，以確定各區(qū)域所對應(yīng)的靈敏度矩陣，計算時間略有增加，由方法1的2.63 s增加至3.20 s。

表2 不同方法的計算時間對比Tab.2 Comparison of computational time for different methods

5 結(jié) 論

針對蒙特卡洛模擬法概率潮流計算效率低下以及輸入隨機變量相關(guān)性考慮不足的問題，本文提出了一種考慮隨機變量相關(guān)性的多點線性化概率潮流計算方法。該方法主要有以下特點：

（1）基于核密度估計建立輸入隨機變量的概率分布模型，對其分布特征刻畫更為準確；

（2）構(gòu)建Copula 函數(shù)描述多維隨機變量間的相關(guān)性，為潮流計算生成更加符合系統(tǒng)運行實際情況的樣本數(shù)據(jù)；

（3）對傳統(tǒng)蒙特卡洛模擬法中的非線性潮流計算方法進行改進，結(jié)合多點線性化技術(shù)對潮流方程進行處理，大幅度提高計算效率的同時有效減小了線性化誤差。

IEEE?30節(jié)點系統(tǒng)算例測試驗證了該方法相對于MCSM的計算效率有大幅提高，同時與單點線性化方法相比，計算精度也有較大提升，可適用于“雙高”電力系統(tǒng)考慮源荷雙側(cè)強不確定性因素下的安全運行與規(guī)劃調(diào)度等分析計算工作，為規(guī)劃調(diào)度人員提供了更全面和更準確的輔助決策信息，具有一定的工程實際應(yīng)用價值。