考慮風(fēng)速時(shí)空相關(guān)性的動(dòng)態(tài)概率潮流計(jì)算

余爽，翁程琳，張程，臧海祥

(1.國網(wǎng)江蘇省電力有限公司檢修分公司，江蘇 南京 211100；2.河海大學(xué) 能源與電氣學(xué)院，江蘇 南京 211100；3.南京航空航天大學(xué)，江蘇 南京 210016)

隨著人類環(huán)保意識(shí)的增強(qiáng)，以及國家“雙碳戰(zhàn)略”的提出，化石燃料的使用將不斷減少，風(fēng)電等新能源發(fā)電形式的應(yīng)用將不斷增強(qiáng)[1-2]。然而風(fēng)電的大量接入會(huì)對(duì)電網(wǎng)產(chǎn)生諸多影響，如潮流越限、諧波等。風(fēng)力機(jī)出力與風(fēng)速等自然因素息息相關(guān)，因此存在天然的不確定性，這對(duì)電網(wǎng)的安全穩(wěn)定運(yùn)行也帶來了新的挑戰(zhàn)。

概率潮流(probabilistic power flow, PPF)[3-4]是處理風(fēng)電并網(wǎng)不確定性的有效方法之一，其求解方法主要包括模擬法[5]、解析法[6-7]和近似法[8]?；诤?jiǎn)單隨機(jī)采樣的蒙特卡羅模擬法(simple random sampling Monte Carlo simulation, SRS-MCS)具有精度高的優(yōu)點(diǎn)，但其計(jì)算效率較低，因此一般將其作為與其他方法比較的評(píng)判標(biāo)準(zhǔn)?；诶〕⒎讲蓸?Latin hypercube sampling, LHS)的MCS在保證一定精度的同時(shí)提高了PPF的計(jì)算效率，且不受輸入變量概率分布類型的約束，具有實(shí)現(xiàn)簡(jiǎn)單、精度高的優(yōu)點(diǎn)。解析法和近似法相較模擬法具有更高的計(jì)算效率，然而解析法在線性化過程中存在誤差，近似法忽略了高階矩中包含的信息，兩者的計(jì)算精度劣于模擬法。

傳統(tǒng)PPF基于單時(shí)間斷面的潮流分布特性，而實(shí)際電力系統(tǒng)是動(dòng)態(tài)變化的，因此近年來學(xué)者提出了動(dòng)態(tài)概率潮流(dynamic probability power flow, DPPF)的概念。文獻(xiàn)[8]提出一種基于奇異值分解和點(diǎn)估計(jì)法的隨機(jī)性處理方法，文中假設(shè)風(fēng)速誤差服從正態(tài)分布，然而實(shí)際風(fēng)速誤差會(huì)存在較大的峰度和偏度，使用正態(tài)分布描述預(yù)測(cè)誤差可能對(duì)DPPF的計(jì)算結(jié)果產(chǎn)生較大的誤差。文獻(xiàn)[9]提出含電動(dòng)汽車、風(fēng)力發(fā)電機(jī)和光伏組件的DPPF。文獻(xiàn)[10]考慮電動(dòng)汽車充電負(fù)荷的配電網(wǎng)DPPF，并假設(shè)預(yù)測(cè)誤差服從正態(tài)分布。文獻(xiàn)[11-12]提出含風(fēng)電、電動(dòng)汽車、負(fù)荷的DPPF，在求解過程中假設(shè)輸入變量相互獨(dú)立，然而實(shí)際電力系統(tǒng)中輸入變量在時(shí)間和空間上具有一定的相關(guān)性，因而忽略時(shí)空相關(guān)性會(huì)對(duì)DPPF的結(jié)果帶來誤差。

基于此，本文首先建立風(fēng)速的動(dòng)態(tài)概率模型，采用非參數(shù)核密度估計(jì)根據(jù)風(fēng)速誤差歷史數(shù)據(jù)擬合得到誤差分布，并考慮風(fēng)速誤差的時(shí)空相關(guān)性。接著，通過等概率變換理論和Nataf變換理論得到具有時(shí)空相關(guān)性的風(fēng)速誤差樣本。利用LHS方法進(jìn)行DPPF計(jì)算，得到節(jié)點(diǎn)電壓幅值和支路潮流的樣本，根據(jù)樣本得到輸出變量的數(shù)字特征、動(dòng)態(tài)概率分布和樣本函數(shù)曲線。最后，將本文方法所得結(jié)果與傳統(tǒng)方法進(jìn)行對(duì)比，驗(yàn)證本文所提方法的正確性和實(shí)用性。

1 風(fēng)速和風(fēng)電場(chǎng)的概率模型

1.1 風(fēng)速動(dòng)態(tài)概率模型

實(shí)際風(fēng)速與是預(yù)測(cè)風(fēng)速之間存在一定誤差，可描述為

vt=Vt+xt.

(1)

式中：vt為實(shí)際風(fēng)速，t=1,2,3,…,T，T為時(shí)段數(shù)；Vt為預(yù)測(cè)風(fēng)速；xt為預(yù)測(cè)誤差。

實(shí)際風(fēng)速在時(shí)間和空間上具有一定相關(guān)性。通過特定的風(fēng)速預(yù)測(cè)方法，可以得到具有時(shí)空相關(guān)性的風(fēng)速預(yù)測(cè)序列[13]。由于風(fēng)速預(yù)測(cè)必然存在一定誤差，因此本文在預(yù)測(cè)值的基礎(chǔ)上，考慮預(yù)測(cè)誤差的時(shí)空相關(guān)性。一般而言，誤差過程xt的自相關(guān)系數(shù)函數(shù)僅是時(shí)間差的函數(shù)，因而xt可以用自相關(guān)系數(shù)平穩(wěn)過程進(jìn)行刻畫。

自相關(guān)系數(shù)平穩(wěn)過程{Xt}具有以下特征：

a){Xt}的一、二階矩都存在；

b){Xt}的自相關(guān)函數(shù)ρ(t1,t2)僅依賴于t1與t2的時(shí)間間隔，即

(2)

式中：C(t1,t2)為協(xié)方差函數(shù)；D(t1)、D(t2)為方差函數(shù)。

由于誤差過程的分布類型和分布參數(shù)需要根據(jù)實(shí)際情況和歷史數(shù)據(jù)進(jìn)行分析，因此本文采用非參數(shù)核密度估計(jì)直接根據(jù)風(fēng)速誤差歷史數(shù)據(jù)得到誤差分布[14]。這種方法避免了假設(shè)分布與實(shí)際分布不符以及參數(shù)估計(jì)帶來的誤差。

假設(shè)Y1,Y2, …,YN是N個(gè)風(fēng)速誤差樣本數(shù)據(jù)，則風(fēng)速誤差概率密度函數(shù)的核估計(jì)

(3)

式中：h為帶寬；K(·)為核函數(shù)。

1.2 風(fēng)電場(chǎng)出力概率模型

風(fēng)電場(chǎng)由安裝在同一位置的幾十甚至上百臺(tái)風(fēng)電機(jī)組組成。風(fēng)電場(chǎng)內(nèi)風(fēng)電機(jī)組的風(fēng)力機(jī)類型和功率特性不盡相同，并且受尾流效應(yīng)和場(chǎng)內(nèi)電氣損耗的影響，風(fēng)電場(chǎng)的功率特性曲線不能簡(jiǎn)單地由單臺(tái)風(fēng)力機(jī)的功率特性曲線疊加得到。為此，文獻(xiàn)[14-16]采用最小二乘法，使用歷史風(fēng)速和風(fēng)電場(chǎng)出力作為樣本，擬合得到風(fēng)電場(chǎng)的功率特性曲線，并用多種分布擬合風(fēng)電場(chǎng)功率的誤差分布?；诖?，本文采用概率模型將風(fēng)電場(chǎng)的功率特性描述為如下分段函數(shù)形式：

(4)

式中：Pw和Pr分別為風(fēng)電場(chǎng)的輸出功率和額定功率；v為場(chǎng)內(nèi)平均風(fēng)速；a、b、c、d為功率特性曲線的系數(shù)；vci、vr和vco分別為風(fēng)電機(jī)組的切入風(fēng)速、額定風(fēng)速和切出風(fēng)速；Pε(v)為風(fēng)電場(chǎng)實(shí)際功率誤差，是與風(fēng)速有關(guān)的隨機(jī)變量[15]。

風(fēng)電場(chǎng)實(shí)際功率誤差服從零均值的正態(tài)分布，其概率密度函數(shù)

(5)

式中α為風(fēng)電場(chǎng)功率誤差的標(biāo)準(zhǔn)差與實(shí)際風(fēng)速之間的比例系數(shù)。

2 基于LHS的DPPF模型

2.1 預(yù)測(cè)誤差時(shí)空相關(guān)性的處理方法

地理位置接近的風(fēng)電場(chǎng)風(fēng)速在時(shí)間和空間上具有較強(qiáng)的相關(guān)性。通過風(fēng)速預(yù)測(cè)得到的各風(fēng)電場(chǎng)的風(fēng)速序列本身存在時(shí)空相關(guān)性，同樣的，風(fēng)速預(yù)測(cè)誤差也存在相關(guān)性。因此，本節(jié)處理風(fēng)速預(yù)測(cè)相對(duì)誤差的時(shí)空相關(guān)性。

設(shè)K座風(fēng)電場(chǎng)，T個(gè)時(shí)段風(fēng)速相對(duì)誤差的矩陣為

(6)

式中：xij為風(fēng)速相對(duì)誤差，xij的概率密度函數(shù)fij(xij)和概率分布函數(shù)Fij(xij)通過1.1節(jié)非參數(shù)核密度估計(jì)獲得。

根據(jù)等概率轉(zhuǎn)換原則，引入相關(guān)標(biāo)準(zhǔn)正態(tài)分布的K×T的隨機(jī)變量矩陣Y，其中第i行第j列的元素

yij=Φ-1(Fij(xij)),

(7)

式中Φ-1(·)為標(biāo)準(zhǔn)正態(tài)分布的逆概率分布函數(shù)。

工程中常用相關(guān)系數(shù)矩陣來描述隨機(jī)變量間的相關(guān)性，各時(shí)段的風(fēng)速相對(duì)誤差的空間相關(guān)系數(shù)矩陣

(8)

式中ρij,t為第i與第j座風(fēng)電場(chǎng)風(fēng)速在第t時(shí)段的相關(guān)系數(shù)。

各座風(fēng)電場(chǎng)的時(shí)間相關(guān)系數(shù)矩陣

(9)

式中k=1,2,…,K，ρk,ij為第k座風(fēng)電場(chǎng)風(fēng)速在第i和第j個(gè)時(shí)段的相關(guān)系數(shù)。

根據(jù)Nataf變換理論，推導(dǎo)出標(biāo)準(zhǔn)正態(tài)分布矩陣Y的時(shí)空相關(guān)系數(shù)矩陣滿足：

(10)

式中：μit和σit分別為第i座風(fēng)電場(chǎng)風(fēng)速在第t個(gè)時(shí)段的期望和標(biāo)準(zhǔn)差；ρij,t和ρk,ij分別為時(shí)空相關(guān)系數(shù)矩陣的時(shí)間和空間分量，φ(·)為標(biāo)準(zhǔn)正態(tài)分布的聯(lián)合概率密度函數(shù)。

求解式(10)即可得等效時(shí)空相關(guān)系數(shù)矩陣CK,t和Ck,T。關(guān)于式(10)的計(jì)算方法，文獻(xiàn)[17]給出了10種典型分布的經(jīng)驗(yàn)公式，文獻(xiàn)[18]采用二分法進(jìn)行求解，文獻(xiàn)[19]則提出Hermite多項(xiàng)式展開的方法進(jìn)行求解。本文采用Hermite多項(xiàng)式方法進(jìn)行求解，該方法具有較高的精度，同時(shí)避免了無窮積分的計(jì)算。

對(duì)CK,t和Ck,T分別進(jìn)行Cholesky分解得到下三角矩陣Bt和Dk：

(11)

根據(jù)Nataf逆變換，可將獨(dú)立標(biāo)準(zhǔn)正態(tài)分布矩陣UK×T轉(zhuǎn)化相關(guān)標(biāo)準(zhǔn)正態(tài)分布矩陣YK×T，如下：

(12)

式中：U、Z和Y分別為獨(dú)立、空間相關(guān)和時(shí)空相關(guān)的標(biāo)準(zhǔn)正態(tài)分布隨機(jī)變量組成的矩陣。

最后，根據(jù)等概率逆變換原理，得到具有時(shí)空相關(guān)性的風(fēng)速相對(duì)誤差樣本

(13)

2.2 LHS方法介紹

LHS是一種分層采樣的方法，通過產(chǎn)生更加均勻的樣本來提高計(jì)算效率。LHS方法的步驟主要分為采樣和排序兩步。

a)采樣。采樣的目的是產(chǎn)生已知分布的樣本。對(duì)于R個(gè)隨機(jī)變量依次進(jìn)行采樣，形成R×N階采樣矩陣X。

b)排序是對(duì)所得的樣本進(jìn)行重新排序的過程，目的是改變樣本間的相關(guān)性，使之滿足問題需要。首先對(duì)各個(gè)隨機(jī)變量的樣本獨(dú)立地進(jìn)行隨機(jī)排序，此時(shí)各隨機(jī)變量樣本間的相關(guān)系數(shù)接近0，近似認(rèn)為各隨機(jī)變量相互獨(dú)立。然后，利用2.1節(jié)提出的時(shí)空相關(guān)性處理方法，將R組獨(dú)立標(biāo)準(zhǔn)正態(tài)分布隨機(jī)變量樣本轉(zhuǎn)化為滿足一定時(shí)空相關(guān)性的任意分布的樣本。

2.3 DPPF計(jì)算

確定性潮流計(jì)算可以用非線性方程組表示，考慮節(jié)點(diǎn)負(fù)荷和風(fēng)電場(chǎng)出力的不確定性，DPPF可以用非線性隨機(jī)方程組表示。

DPPF計(jì)算流程如圖1所示。

圖1 DPPF的計(jì)算流程Fig.1 Calculation flow chart of DPPF

3 算例分析

3.1 IEEE 14節(jié)點(diǎn)系統(tǒng)

選取江蘇某地某季節(jié)的實(shí)際風(fēng)速作為樣本，采用反向傳播(back propagation, BP)神經(jīng)網(wǎng)絡(luò)法作為日前風(fēng)速預(yù)測(cè)方法[20]，預(yù)測(cè)時(shí)間尺度為日前24 h，預(yù)測(cè)結(jié)果時(shí)間分辨率為1 h，得到該地這一季節(jié)風(fēng)速相對(duì)誤差樣本，利用非參數(shù)核密度估計(jì)得到相對(duì)誤差的概率分布如圖2所示。

圖2 風(fēng)速相對(duì)誤差的概率分布曲線Fig.2 Probability distribution curve of wind speed relative error

從圖2可以看出，利用非參數(shù)核密度估計(jì)可以根據(jù)歷史數(shù)據(jù)得到風(fēng)速相對(duì)誤差的概率分布，從而客觀描述風(fēng)速相對(duì)誤差的分布特性。

以文獻(xiàn)[21]提供的IEEE 14節(jié)點(diǎn)為例，采用基于LHS-MCS計(jì)及預(yù)測(cè)誤差時(shí)空相關(guān)性的DPPF模型進(jìn)行計(jì)算。在節(jié)點(diǎn)10、11、13和14處分別加入額定功率50 MW的風(fēng)電場(chǎng)，對(duì)應(yīng)風(fēng)電場(chǎng)標(biāo)號(hào)分別為1、2、3、4。風(fēng)電機(jī)組均以恒功率因數(shù)0.98運(yùn)行，風(fēng)電機(jī)組的切入風(fēng)速、額定風(fēng)速和切出風(fēng)速分別取為3 m/s、12 m/s和22 m/s，α取為0.08 MW·s/m。

模擬得到風(fēng)電場(chǎng)的實(shí)際出力的散點(diǎn)圖如圖3所示。由圖3可以看出，風(fēng)電場(chǎng)的功率特性采用考慮輸出功率誤差的概率模型，其功率散點(diǎn)按概率分布在基于預(yù)測(cè)值的一個(gè)范圍內(nèi)，更容易匹配到實(shí)際風(fēng)電場(chǎng)的出力情況。

采用BP神經(jīng)網(wǎng)絡(luò)法進(jìn)行日前風(fēng)速預(yù)測(cè)，得到4座風(fēng)電場(chǎng)的預(yù)測(cè)風(fēng)速序列如圖4所示。以圖2所示的概率分布曲線作為相對(duì)誤差的概率分布，考慮風(fēng)速相對(duì)誤差的時(shí)空相關(guān)性。算例中給定負(fù)荷數(shù)據(jù)作為負(fù)荷的期望，負(fù)荷期望的5%作為負(fù)荷的標(biāo)準(zhǔn)差，并考慮負(fù)荷的空間相關(guān)性。

圖3 模擬風(fēng)電場(chǎng)輸出功率Fig.3 Simulation of wind power output

圖4 風(fēng)速預(yù)測(cè)曲線Fig.4 Forecast curves of wind speed

采用LHS-MCS進(jìn)行500次采樣的DPPF計(jì)算，得到各節(jié)點(diǎn)電壓和支路潮流的概率分布曲線和概率密度曲線。所選節(jié)點(diǎn)的位置和節(jié)點(diǎn)的個(gè)數(shù)并不影響仿真結(jié)果，因此隨機(jī)選取節(jié)點(diǎn)9的電壓幅值作為研究對(duì)象。節(jié)點(diǎn)9電壓幅值的概率密度曲線如圖5所示。

圖5 節(jié)點(diǎn)9電壓幅值的概率密度曲線Fig.5 Probability density curves of node 9 voltage magnitudes

從圖5可以看出，基于LHS的DPPF計(jì)算可以準(zhǔn)確得到節(jié)點(diǎn)電壓幅值的概率密度曲線，并且能夠有效處理服從任意分布的風(fēng)速誤差之間的時(shí)空相關(guān)性。

針對(duì)以下4種方案考慮時(shí)空相關(guān)性對(duì)系統(tǒng)運(yùn)行特性的影響：

a)方案1，考慮時(shí)空相關(guān)性；

b)方案2，只考慮空間相關(guān)性；

c)方案3，只考慮時(shí)間相關(guān)性；

d)方案4，不考慮相關(guān)性。

各方案下，節(jié)點(diǎn)9在各時(shí)段的節(jié)點(diǎn)電壓幅值標(biāo)準(zhǔn)差如圖6所示。

圖6 各時(shí)段節(jié)點(diǎn)9的電壓幅值的標(biāo)準(zhǔn)差Fig.6 Standard deviations of node 9 voltage magnitudes in each period

從圖6可以看出，風(fēng)速預(yù)測(cè)誤差的時(shí)間和空間相關(guān)性都會(huì)增大節(jié)點(diǎn)電壓的波動(dòng)性，并且系統(tǒng)運(yùn)行特性的影響受空間相關(guān)性的影響較大。這是因?yàn)楫?dāng)某風(fēng)電場(chǎng)在某時(shí)刻的風(fēng)速誤差較大時(shí)，與之地理位置相鄰的風(fēng)電場(chǎng)在相鄰時(shí)段上也會(huì)以較大的概率出現(xiàn)大的風(fēng)速誤差。風(fēng)速誤差時(shí)空相關(guān)性的存在，導(dǎo)致輸出電壓的波動(dòng)性增加，波動(dòng)范圍增大。

通過預(yù)測(cè)數(shù)據(jù)得到節(jié)點(diǎn)9的電壓幅值預(yù)測(cè)曲線、預(yù)測(cè)日實(shí)際的電壓幅值曲線、通過傳統(tǒng)方法假設(shè)風(fēng)速服從Weibull分布進(jìn)行DPPF得到的電壓幅值置信度為99%的置信區(qū)間，以及本文的DPPF方法得到電壓幅值置信區(qū)間的對(duì)比如圖7所示。

從圖7可以看出，通過計(jì)及預(yù)測(cè)誤差時(shí)空相關(guān)性的DPPF所得結(jié)果與傳統(tǒng)的基于風(fēng)速服從Weibull分布所得的結(jié)果相比，計(jì)算結(jié)果的保守性得到降低。由于傳統(tǒng)方法每個(gè)時(shí)段的風(fēng)速都服從Weibull分布，風(fēng)速會(huì)在一個(gè)較大的區(qū)間波動(dòng)；本文方法的風(fēng)速基于預(yù)測(cè)值以誤差所在區(qū)間變化，通過歷史數(shù)據(jù)擬合得到誤差分布，風(fēng)速只在一個(gè)較小的區(qū)間波動(dòng)，降低了結(jié)果的保守性，更加接近實(shí)際運(yùn)行情況，便于調(diào)度運(yùn)行人員了解系統(tǒng)的真實(shí)情況。

圖7 不同方法得到的節(jié)點(diǎn)9電壓幅值Fig.7 Node 9 voltage magnitudes calculated by different methods

3.2 IEEE 118節(jié)點(diǎn)系統(tǒng)

以文獻(xiàn)[21]提供的IEEE 118節(jié)點(diǎn)為例，采用基于LHS-MCS的DPPF模型進(jìn)行計(jì)算。在節(jié)點(diǎn)20、35、38和51分別加入額定功率50 MW的風(fēng)電場(chǎng)。

所選支路的位置和個(gè)數(shù)并不影響仿真結(jié)果，因此隨機(jī)選取支路19-20基于預(yù)測(cè)值得到的有功功率預(yù)測(cè)曲線、實(shí)際的有功功率曲線，以及基于本文和傳統(tǒng)的DPPF方法得到的有功功率置信區(qū)間對(duì)比如圖8所示。

圖8 不同方法得到的線路19-20有功潮流Fig.8 Active power in line 19-20 calculated by different methods

從圖8可以看出，根據(jù)負(fù)荷和風(fēng)速預(yù)測(cè)值得到的支路潮流與實(shí)際情況存在差異，這種誤差是常規(guī)的基于日前預(yù)測(cè)的潮流計(jì)算所無法避免的。而本文提出的基于預(yù)測(cè)值計(jì)及預(yù)測(cè)誤差時(shí)空相關(guān)性的DPPF方法可以減少這種差異。通過區(qū)間的形式可以反映所有可能的系統(tǒng)運(yùn)行情況，具有良好的工程應(yīng)用性。

4 結(jié)論

本文提出了一種計(jì)及風(fēng)速預(yù)測(cè)誤差時(shí)空相關(guān)性的DPPF算法。計(jì)算分析表明：

a)風(fēng)速預(yù)測(cè)誤差的時(shí)空相關(guān)性均會(huì)增大節(jié)點(diǎn)電壓的波動(dòng)性，而空間相關(guān)性對(duì)系統(tǒng)運(yùn)行特性的影響較大；

b)本文所提的計(jì)及預(yù)測(cè)誤差時(shí)空相關(guān)性的DPPF方法能夠反映系統(tǒng)所有可能的運(yùn)行情況，降低了結(jié)果的保守性，更加符合實(shí)際運(yùn)行情況，具有良好的工程應(yīng)用前景。