李曉，李滿禮，胡云龍，劉國(guó)靜，胡曉燕，倪明

（1．南瑞集團(tuán)有限公司（國(guó)網(wǎng)電力科學(xué)研究院有限公司），江蘇省 南京市 211106；2．智能電網(wǎng)保護(hù)和運(yùn)行控制國(guó)家重點(diǎn)實(shí)驗(yàn)室，江蘇省 南京市 211106；3．國(guó)網(wǎng)江蘇省電力有限公司經(jīng)濟(jì)技術(shù)研究院，江蘇省 南京市 210008）

0 引言

隨著氣候變化加劇，傳統(tǒng)化石能源供應(yīng)日益緊張，伴隨而來的能源效率低、供需矛盾突出、結(jié)構(gòu)不合理等問題迫切需要能源朝著可再生、綜合互補(bǔ)利用和低碳化的方向發(fā)展[1]。構(gòu)建清潔可持續(xù)的新型電力系統(tǒng)、保障能源安全已成為各個(gè)國(guó)家的重要戰(zhàn)略目標(biāo)，中國(guó)也提出分別于2030年前和2060年前實(shí)現(xiàn)碳達(dá)峰和碳中和的“雙碳”目標(biāo)，推動(dòng)能源低碳轉(zhuǎn)型[2]。傳統(tǒng)的單一供能網(wǎng)絡(luò)缺少不同種類能源之間的相互協(xié)調(diào)，能源利用率較低，而能源之間必須相互協(xié)調(diào)才能夠最優(yōu)化資源配置，達(dá)到低碳和可持續(xù)的目的。因此，能夠滿足多種能源互補(bǔ)運(yùn)行的綜合能源系統(tǒng)（integrated energy system, IES）應(yīng)運(yùn)而生[3]。其中，利用化石燃料、可再生能源、電能等方式同時(shí)發(fā)電和供熱的熱電聯(lián)供網(wǎng)絡(luò)在冷、熱、電等多種能源形式聯(lián)合供應(yīng)的多能流網(wǎng)絡(luò)中發(fā)展最為迅速[4]。

在電力系統(tǒng)中，狀態(tài)估計(jì)（state estimation, SE）利用能量管理系統(tǒng)（energy management system, EMS）收集到的冗余信息數(shù)據(jù)，按照一定準(zhǔn)則對(duì)原始數(shù)據(jù)進(jìn)行處理，估算出更加精確的狀態(tài)量以及系統(tǒng)所處的運(yùn)行狀態(tài)，是電力系統(tǒng)運(yùn)行與控制的重要基礎(chǔ)。該領(lǐng)域研究已較為成熟，在抗差估計(jì)、計(jì)及虛假數(shù)據(jù)注入攻擊、計(jì)算效率等方面均有涉及[5-8]。隨著“雙碳”背景下綜合能源系統(tǒng)的快速建設(shè)，跨區(qū)級(jí)和園區(qū)級(jí)綜合能源系統(tǒng)的規(guī)模越來越大，采集和傳輸?shù)臄?shù)據(jù)種類、數(shù)據(jù)量也愈加龐大和復(fù)雜。為更準(zhǔn)確地掌握綜合能源系統(tǒng)的運(yùn)行狀態(tài)，需要對(duì)綜合能源系統(tǒng)進(jìn)行狀態(tài)估計(jì)，得到精度更高的全局一致基態(tài)潮流解，從而為后續(xù)的系統(tǒng)在線安全評(píng)估、優(yōu)化調(diào)度提供堅(jiān)實(shí)基礎(chǔ)[9]。目前，對(duì)綜合能源系統(tǒng)狀態(tài)估計(jì)的研究還處于起步階段。文獻(xiàn)[10]對(duì)比分析了電力系統(tǒng)SE和IES-SE的異同，對(duì)已有的靜態(tài)及動(dòng)態(tài)IES-SE模型和方法進(jìn)行了剖析，指出了加快研究面向多能流的IES-SE的迫切性。文獻(xiàn)[4]和文獻(xiàn)[11]利用加權(quán)最小二乘法（weighted least squares,WLS）分別對(duì)電熱耦合系統(tǒng)和電氣互聯(lián)系統(tǒng)進(jìn)行了狀態(tài)估計(jì)，并驗(yàn)證了聯(lián)合估計(jì)相比獨(dú)立系統(tǒng)估計(jì)可降低整體估計(jì)誤差，具有重要的理論意義和應(yīng)用價(jià)值。但是WLS沒有抗差性，當(dāng)量測(cè)數(shù)據(jù)含有粗差時(shí)，精度會(huì)受到很大影響。文獻(xiàn)[12]提出了一種用于電-氣-熱IES的三階段分布式狀態(tài)估計(jì)方法，該方法能夠獲得具有足夠精度和可靠收斂性的全局IES估計(jì)結(jié)果，顯著提高了估計(jì)效率，但該文也未考慮不良數(shù)據(jù)的影響。

基于WLS的狀態(tài)估計(jì)是誤差數(shù)據(jù)服從高斯分布前提下的最優(yōu)估計(jì)算法，具有簡(jiǎn)單、快速等優(yōu)點(diǎn)，但實(shí)際測(cè)量系統(tǒng)中誤差數(shù)據(jù)的分布不一定遵循高斯分布，會(huì)對(duì)狀態(tài)估計(jì)結(jié)果產(chǎn)生負(fù)面影響。為克服WLS估計(jì)的非抗差性，業(yè)內(nèi)學(xué)者、專家相繼提出了多種抗差估計(jì)方法。常見方法有：極大似然性估計(jì)（maximum likelihood estimation, ME）、順序統(tǒng)計(jì)量線性組合估計(jì)（linear combination estimation of order statistics, LE）、非參數(shù)秩檢驗(yàn)估計(jì)（nonparametric rank test estimation,RE）等[13]。其中，ME估計(jì)是電力系統(tǒng)中應(yīng)用最為廣泛的抗差估計(jì)方法，如加權(quán)最小絕對(duì)值法（weighted least absolute value, WLAV）[14-15]、指數(shù)估計(jì)[16]均屬于ME估計(jì)。由于LE和RE估計(jì)方法不容易推導(dǎo)線性觀測(cè)方程，在電力系統(tǒng)狀態(tài)估計(jì)中的應(yīng)用相對(duì)較少。此外，當(dāng)前電熱綜合能源系統(tǒng)狀態(tài)估計(jì)領(lǐng)域也實(shí)現(xiàn)了一些抗差估計(jì)算法的應(yīng)用，如WLAV、指數(shù)估計(jì)、雙線性抗差估計(jì)等。總的來說，WLAV方法應(yīng)用較多，但其抗差效果還具有提升空間；指數(shù)型抗差估計(jì)方法抗差效果較好，但計(jì)算求解效率偏低；雙線性抗差估計(jì)在將非線性模型轉(zhuǎn)換為線性模型的過程中并不一定能保證求解精度。而屬于ME估計(jì)的Huber抗差估計(jì)方法因其算法特點(diǎn)，可以實(shí)現(xiàn)抗差精度與求解效率的良好折中，對(duì)IES-SE問題具有較強(qiáng)的適用性。

綜上所述，為提高IES狀態(tài)估計(jì)抑制粗差的能力，本文以電熱綜合能源系統(tǒng)為研究對(duì)象，提出一種基于Huber估計(jì)的綜合能源系統(tǒng)抗差估計(jì)方法。首先建立電熱綜合能源系統(tǒng)量測(cè)模型，其次將Huber抗差估計(jì)方法引入到電-熱聯(lián)合狀態(tài)估計(jì)模型中，最后，通過與WLS方法在不同比例的不良數(shù)據(jù)條件下進(jìn)行對(duì)比，驗(yàn)證所提Huber估計(jì)方法的抗差性和有效性。

1 Huber抗差估計(jì)數(shù)學(xué)原理

1.1 最小二乘估計(jì)

設(shè)觀測(cè)向量為z= [z1, z2, ···,zn]T，狀態(tài)向量為x= [x1, x2, ···,xn]T，誤差向量為v=[v1,v2, ···,vn]T，則線性觀測(cè)方程可記為

式中：h表示非線性函數(shù)關(guān)系。

作為當(dāng)前求解狀態(tài)估計(jì)問題應(yīng)用最為廣泛的WLS法，根據(jù)狀態(tài)變量與量測(cè)矢量之間的關(guān)系式，其本質(zhì)為求解如下數(shù)學(xué)問題[17]：

式中：w為權(quán)重矩陣，一般取量測(cè)誤差向量的方差矩陣的逆，即

式中：R為量測(cè)誤差向量的方差矩陣，形式為對(duì)角矩陣，且取值如下：

在大多數(shù)情況下，正態(tài)分布接近于真實(shí)的觀察世界。當(dāng)量測(cè)數(shù)據(jù)服從正態(tài)分布時(shí)，WLS估計(jì)靈活，模型簡(jiǎn)單且易于計(jì)算，因此有著廣泛的應(yīng)用[19]。然而，當(dāng)量測(cè)量存在粗差時(shí)，將對(duì)估計(jì)結(jié)果產(chǎn)生很大影響，使估計(jì)結(jié)果精度降低。而Huber抗差估計(jì)對(duì)量測(cè)量進(jìn)行等價(jià)權(quán)處理，當(dāng)量測(cè)數(shù)據(jù)有粗差時(shí)，可以得到比WLS更精確的估計(jì)結(jié)果[20]。

1.2 Huber抗差估計(jì)

WLS估計(jì)精度會(huì)受到殘差過大的量測(cè)值的影響，為克服WLS估計(jì)的非抗差性，本文重點(diǎn)研究了屬于ME估計(jì)類的Huber抗差估計(jì)方法。

ME準(zhǔn)則可由如下目標(biāo)函數(shù)表示：

令ρ=-lnf，則ME的目標(biāo)函數(shù)可轉(zhuǎn)換為

令φ=ρ′，由上式求極值可得到：

式中：ai為h函數(shù)的雅可比矩陣第i行行向量。

令φi=ψivi，式（7）可記為

令P=ψw，w為原量測(cè)量的權(quán)重。式（8）的矩陣形式為

式（9）的解如下：

Huber估計(jì)能結(jié)合均值無偏最小方差估計(jì)的靈敏度和絕對(duì)值無偏估計(jì)的穩(wěn)健性，其ρ函數(shù)為

則φ函數(shù)、ψ函數(shù)如下：

式中：當(dāng)粗差比例在1%～10%時(shí)，k值在1～2之間[21]。由式 （8）可知，當(dāng)k=∞時(shí)，即ψ=1，Huber估計(jì)等價(jià)于WLS估計(jì)。

由以上推導(dǎo)過程可知，Huber抗差估計(jì)是通過改變殘差較大的量測(cè)量權(quán)重實(shí)現(xiàn)的。當(dāng)測(cè)量值對(duì)應(yīng)的殘差超過極限閾值時(shí)，通過權(quán)重處理降低測(cè)量值對(duì)應(yīng)的權(quán)重值，當(dāng)不超過閾值時(shí)，權(quán)重值保持不變。因此，在每次迭代中，使殘差保持在一定范圍內(nèi)，從而降低較大殘差對(duì)估計(jì)結(jié)果的影響，保證抗差估計(jì)抵御粗差的能力。

2 電熱綜合能源系統(tǒng)量測(cè)與狀態(tài)估計(jì)模型

2.1 電力系統(tǒng)量測(cè)模型

在電力系統(tǒng)中，需要估計(jì)的狀態(tài)量為電壓幅值（V）和相角（θ），記電力系統(tǒng)狀態(tài)量為x V

e=[θ;]。電力系統(tǒng)的量測(cè)模型如下：

式中：上標(biāo)MEA代表量測(cè)值；vVi為Vi的量測(cè)誤差；Pi、Qi分別為節(jié)點(diǎn)i注入有功、無功功率；Gij、Bij分別為節(jié)點(diǎn)i和節(jié)點(diǎn)j之間的互電導(dǎo)和互電納；Pij、Qij分別為支路ij上節(jié)點(diǎn)i流向節(jié)點(diǎn)j的有功、無功功率；gij、bij分別為支路ij的電導(dǎo)與電納；yc為對(duì)地電納；Vi為節(jié)點(diǎn)i的電壓幅值；θij為支路ij兩端的相角差，且θij= θi- θj。

2.2 熱力系統(tǒng)量測(cè)模型

供熱網(wǎng)絡(luò)的分析包括水力模型與熱力模型兩方面，水力模型量測(cè)量包括節(jié)點(diǎn)壓強(qiáng)p、支路流量mij、節(jié)點(diǎn)注入流量mq，熱力模型中的量測(cè)量包括節(jié)點(diǎn)供熱溫度Ts、節(jié)點(diǎn)回?zé)釡囟萒r、節(jié)點(diǎn)熱功率?[22]，狀態(tài)量一般選取節(jié)點(diǎn)壓強(qiáng)、節(jié)點(diǎn)供熱和回?zé)釡囟龋薀崃ο到y(tǒng)的狀態(tài)量可記為xh=[p; Ts; Tr]。

熱力系統(tǒng)量測(cè)方程可表示為式（15）的形式：

式中：pi為熱網(wǎng)節(jié)點(diǎn)i的壓強(qiáng)；vpi為pi的量測(cè)誤差；mij為支路ij的水流量；mqi為節(jié)點(diǎn)i的流入水流量之和；Kij為與管道直徑、材料以及管道內(nèi)液態(tài)種類有關(guān)的系數(shù)，具體計(jì)算方法可參考文獻(xiàn)[23]；?i是節(jié)點(diǎn)i的熱功率；Cp為水的比熱容系數(shù)；Tsi、Tri分別代表節(jié)點(diǎn)i的供熱溫度和回?zé)釡囟取?/p>

2.3 電熱耦合元件量測(cè)模型

電熱綜合能源系統(tǒng)中常見的耦合元件有熱電聯(lián)產(chǎn)（combined heat and power, CHP）機(jī)組、電鍋爐（electric boilers, EB）、循環(huán)泵等。CHP機(jī)組是熱電聯(lián)產(chǎn)系統(tǒng)中最重要的耦合元件，既可以產(chǎn)生電能，也可以產(chǎn)生熱能。EB也是常見的能量轉(zhuǎn)換耦合單元，依靠電氣元件來產(chǎn)生熱能，同時(shí)不產(chǎn)生燃燒類的化學(xué)反應(yīng)，無黑煙、二氧化硫、二氧化碳等廢棄物排放。CHP機(jī)組、EB能量轉(zhuǎn)換具有如下的數(shù)學(xué)關(guān)系：

式中：?、P分別表示熱功率和電功率；c是熱電生產(chǎn)比；η是電加熱的轉(zhuǎn)換效率。

2.4 電熱綜合能源系統(tǒng)狀態(tài)估計(jì)模型

在電熱耦合網(wǎng)絡(luò)狀態(tài)估計(jì)問題中，狀態(tài)量和量測(cè)量可分別表示為

式中：xe、xh分別是電網(wǎng)、熱網(wǎng)狀態(tài)量；ze、zh分別是電網(wǎng)、熱網(wǎng)量測(cè)量。其中，ze=[ Vi; Pi; Qi; Pij; Qij]，zh=[ pi; mij; mqi;?i]。若采用加權(quán)最小二乘法對(duì)電熱綜合能源系統(tǒng)進(jìn)行狀態(tài)估計(jì)，可表示為如下數(shù)學(xué)問題：

式中：he、hh分別為電網(wǎng)、熱網(wǎng)的非線性函數(shù)；ze、zh分別為電網(wǎng)、熱網(wǎng)的量測(cè)量（詳見式（14）和式（15））；ωe、ωh分別為電網(wǎng)、熱網(wǎng)量測(cè)量權(quán)重；等式約束c(x)為耦合元件約束（詳見式（16）和式（17））。

Huber估計(jì)中通過引用ψ函數(shù)，降低殘差過大的量測(cè)量權(quán)重，從而達(dá)到抗差效果，可轉(zhuǎn)換為如下數(shù)學(xué)問題：

式中：ψ函數(shù)具體計(jì)算方法參見式（13）。

2.5 冗余度分析

量測(cè)冗余度可定義為量測(cè)量個(gè)數(shù)n與待估計(jì)狀態(tài)量個(gè)數(shù)m的比值n/m，狀態(tài)估計(jì)就是利用較大的冗余度來提升數(shù)據(jù)精度。在電熱綜合能源系統(tǒng)中，由于耦合元件增大了系統(tǒng)的冗余度，使得整個(gè)電熱綜合能源系統(tǒng)的狀態(tài)估計(jì)更加準(zhǔn)確。

對(duì)于電熱綜合能源系統(tǒng)，其潮流方程可寫為

式中：ge為電網(wǎng)真值；gh為熱網(wǎng)真值；fe為電網(wǎng)潮流方程；fh為熱網(wǎng)潮流方程；f為耦合元件約束條件。

在電熱綜合能源系統(tǒng)中，若已知部分狀態(tài)量，可通過潮流方程或耦合元件約束求得另一系統(tǒng)狀態(tài)量，從而在量測(cè)量數(shù)目不變的條件下，降低狀態(tài)量個(gè)數(shù)，增大系統(tǒng)冗余度。例如，電網(wǎng)節(jié)點(diǎn)和熱網(wǎng)節(jié)點(diǎn)由于CHP機(jī)組而產(chǎn)生熱電耦合，若對(duì)電熱綜合能源系統(tǒng)進(jìn)行聯(lián)合狀態(tài)估計(jì)，此時(shí)若電網(wǎng)節(jié)點(diǎn)的狀態(tài)量xe為已知，通過電網(wǎng)潮流方程可計(jì)算電熱耦合節(jié)點(diǎn)相應(yīng)電功率PCHP。由于式（16）的存在，與之相對(duì)應(yīng)的熱網(wǎng)節(jié)點(diǎn)的熱功率?CHP即為已知。此時(shí)，若已知Pi，Tsi，Tri三者中任何兩個(gè)量，即可確定第三個(gè)量，也就是狀態(tài)量個(gè)數(shù)減1。根據(jù)冗余度的定義，此時(shí)量測(cè)量不變的前提下，狀態(tài)量數(shù)目減少，此時(shí)冗余度變大。隨著耦合元件增多，冗余度也隨之變大，也越有利于狀態(tài)估計(jì)。

3 算例分析

本文電熱綜合能源系統(tǒng)算例由IEEE 14節(jié)點(diǎn)電網(wǎng)和2個(gè)獨(dú)立供熱網(wǎng)絡(luò)構(gòu)成，如圖1和圖2所示。熱網(wǎng)包括2個(gè)CHP機(jī)組以及14個(gè)熱負(fù)荷節(jié)點(diǎn)。熱源節(jié)點(diǎn)CHP1和CHP2分別對(duì)應(yīng)電網(wǎng)中的2號(hào)和3號(hào)節(jié)點(diǎn)。在熱力網(wǎng)絡(luò)2號(hào)和11號(hào)節(jié)點(diǎn)安裝了EB，分別對(duì)應(yīng)電網(wǎng)節(jié)點(diǎn)6和節(jié)點(diǎn)8。該系統(tǒng)詳細(xì)參數(shù)可參見文獻(xiàn)[24]。

為驗(yàn)證本文所提方法的有效性，采用MATLAB R2016b平臺(tái)進(jìn)行仿真，仿真測(cè)試的計(jì)算機(jī)配置為Intel Core i5-8250 U @1.80 GHz CPU，內(nèi)存為8 GB，操作系統(tǒng)為Win10，仿真用時(shí)0.87 s。

量測(cè)數(shù)據(jù)以潮流計(jì)算結(jié)果為真值的基礎(chǔ)上疊加2%的誤差構(gòu)成。在驗(yàn)證抗差性能的過程中，人為隨機(jī)添加壞數(shù)據(jù)，本文設(shè)置壞數(shù)據(jù)的方式有置零、取反、減小一半、增大一倍4種，如表1所示。

表1 壞數(shù)據(jù)設(shè)置方式Table 1 Setting method of bad data

本文使用平均誤差EM和最大誤差EMA評(píng)價(jià)抗差性能，定義如下：

式中：xi是狀態(tài)量真值；是狀態(tài)量的估計(jì)值。

電力側(cè)電壓幅值與相角誤差結(jié)果如圖3和圖4所示，熱力側(cè)節(jié)點(diǎn)壓強(qiáng)與供熱、回?zé)釡囟日`差結(jié)果如圖5、圖6和圖7所示。圖示結(jié)果為電網(wǎng)、熱網(wǎng)所有節(jié)點(diǎn)上的誤差平均值；對(duì)比的WLS估計(jì)算法未考慮不良數(shù)據(jù)辨識(shí)環(huán)節(jié)。

由圖3可知，Huber估計(jì)相比于WLS估計(jì)，對(duì)電壓幅值和相角具有很好的抗差效果，能夠有效抑制粗誤差對(duì)電壓幅值和相角平均誤差的影響。從圖4中可以看出，Huber估計(jì)在降低電壓幅值與相角最大誤差方面，性能不如降低平均誤差優(yōu)異，但也能有效降低最大誤差。

在圖5中，2種估計(jì)方法的壓強(qiáng)平均誤差和最大誤差增長(zhǎng)趨勢(shì)大致相同，Huber估計(jì)隨著壞數(shù)據(jù)比例的上升，降低誤差的效果也越來越明顯。從圖6、圖7可以看出，熱網(wǎng)供熱、回?zé)釡囟鹊恼`差相差不大，平均誤差和最大誤差隨壞數(shù)據(jù)比例的增長(zhǎng)趨勢(shì)也基本一致，抗差效果相似，Huber估計(jì)的平均誤差和最大誤差相較于WLS均顯著降低。

從圖3—圖7可以看出，在沒有壞數(shù)據(jù)的情況下，Huber估計(jì)與WLS估計(jì)的估計(jì)誤差相差不大。當(dāng)壞數(shù)據(jù)比例較低時(shí)，抗差效果并不明顯。隨著壞數(shù)據(jù)比例逐漸增大，2種方法的平均誤差與最大誤差均隨之增大，但Huber估計(jì)能有效降低狀態(tài)量的平均誤差和最大誤差。在最大壞數(shù)據(jù)比例下，電壓幅值、相角的平均誤差分別下降49.7%、43.4%，壓強(qiáng)平均誤差下降40.9%。圖6中供熱溫度、回?zé)釡囟绕骄`差分別下降52.7%、48.9%，圖7中供熱溫度、回?zé)釡囟茸畲笳`差分別下降50.2%、49.8%，均可獲得良好的抗差效果。因此，Huber估計(jì)方法總體上具有良好的抗差性，能有效抵御壞數(shù)據(jù)對(duì)IES狀態(tài)估計(jì)結(jié)果的不良影響。

4 結(jié)論

本文提出一種考慮量測(cè)誤差的電熱綜合能源系統(tǒng)抗差狀態(tài)估計(jì)方法，以抵御量測(cè)數(shù)據(jù)中的壞數(shù)據(jù)對(duì)狀態(tài)估計(jì)結(jié)果精度的影響。算例仿真結(jié)果表明，本文采用的Huber估計(jì)方法相比于WLS估計(jì)，可在壞數(shù)據(jù)比例增大時(shí)顯著降低狀態(tài)估計(jì)的平均誤差與最大誤差，且隨著壞數(shù)據(jù)占比提高，抗差效果越突出。綜上，本文提出的抗差估計(jì)方法可以有效降低壞數(shù)據(jù)對(duì)系統(tǒng)狀態(tài)估計(jì)精度的影響，提升“雙碳”背景下綜合能源系統(tǒng)狀態(tài)感知的能力，保障系統(tǒng)的安全運(yùn)行。

本文只考慮了包含2種能量形式的電熱綜合能源系統(tǒng)，后續(xù)研究中還可以考慮更多能量形式的耦合系統(tǒng)。另外，針對(duì)不良數(shù)據(jù)的相關(guān)性強(qiáng)弱差異進(jìn)行算法抗差性能分析，以及熱網(wǎng)的動(dòng)態(tài)時(shí)變特性也將在未來工作中進(jìn)一步研究。