基于動力學系統的狀態估計算法

趙天騏，江曉東，李 鵬，白 浩，于 力

（1.天津大學電氣自動化與信息工程學院，天津 300072；2.南方電網科學研究院有限責任公司，廣州 510080）

電力系統狀態估計是能量管理系統EMS（ener?gy management system）的主要功能之一，在監測和控制電力系統的可靠運行方面發揮著重要作用。狀態估計器通過電壓、電流、線路功率、注入功率等測量來估計整個系統的狀態，即每個母線的電壓幅值和相角。這些測量由數據采集與監視控制系統的遠動終端提供。狀態估計結果為其他EMS功能，如安全分析、最優潮流、電壓穩定性分析等提供基礎。因此狀態估計的收斂性是電網在線實時監測和控制的關鍵。

基于高斯-牛頓迭代法的加權最小二乘WLS（weighted least squares）算法是最常用的狀態估計算法。然而牛頓法對初始值比較敏感，初始點選擇不當會造成算法發散，尤其在以下的情況難以獲得狀態估計解：

（1）如果測量系統中出現嚴重的測量錯誤或電網數據中有拓撲錯誤，基于牛頓迭代法的WLS估計器可能會發生不收斂；

（2）如果系統中含有不可觀測的狀態變量，這時增益矩陣會發生奇異，傳統基于牛頓迭代的WLS算法將無法進行計算。

從數值的角度來看，當存在很大的注入測量或者采用很大的權值時，狀態估計問題可能會出現病態。正交變換法[1]可以在一定程度上解決病態問題。Peters-Wilkinson方法[2]在求解速度和數值穩定性之間進行了很好的平衡。Hachtel增廣矩陣法[3]和Cholesky分解法[4]也在克服數值病態問題上有重要貢獻。文獻[5-6]指出正交分現震蕩發散而能給出狀態估計的解，提出了基于信賴域（trust region）的狀態估計算法來克服這個問題。

本文提出基于動力學系統的狀態估計算法，基本原理是將狀態估計問題轉化為求解非線性動力系統的穩定平衡點SEP（stable equilibrium point），通過構造商梯度系統QGS（quotient gradient system），使其退化的穩定平衡流形對應于WLS的解，從非線性動力學系統理論出發，證明了這個動力學系統是漸進穩定的，因此通過對QGS進行軌跡追蹤，可以得到狀態估計的解。由于所提出的方法不需要對增益矩陣進行求逆，因此能夠克服傳統WLS算法的局限性。

1 電力系統狀態估計問題

WLS算法是電網狀態估計最常用的方法，其實現簡單并且計算速度快。該方法旨在找到測量值與其對應測量方程之間的最佳擬合，從而確定出系統的狀態變量。由于潮流方程是非線性的，因此狀態估計問題也是非線性最小二乘問題。定義以下的列向量。

對于一個有N個節點的電力網絡可表示為

式中：z為測量向量；x為系統的狀態變量；h(x)為測量方程；hi(x)為關于測量i的非線性測量方程；e為測量誤差列向量。

最小二乘估計是目標函數J(x)的最小化問題，可表示為

式中，Rii為測量誤差方差矩陣的對角元素。

以上最小二乘狀態估計問題通常采用高斯-牛頓法求解一階最優條件，可表示為

迭代方程由g(x)的一階泰勒展開并忽略高階項導出為

式中：Δx為迭代步長；xk為第k次迭代后的狀態變量。

增益矩陣G可以表示為

式中，R為測量誤差方差矩陣。

高斯-牛頓法的收斂取決于增益矩陣G是否可逆，也就是說必須有足夠的測量分布在整個網絡上，從而可以推斷出2N-1個狀態變量的信息。當有足夠的測量用于所有狀態變量時，稱該系統為可觀測的。一個可觀測系統的狀態可以用測量來唯一地確定，然而測量并非總是在整個網絡中均勻分布。系統的一些子區域可能會因為測量不足或測量丟失而導致失去可觀測性。在這種情況下，傳統的高斯-牛頓法就不能給出狀態估計的解，而對于不可觀測的系統，狀態估計的解不唯一。

2 基于動力學系統的狀態估計算法

基于動力學系統的方法基本思想是將狀態估計問題轉化為適當構造的動力系統，狀態估計的解可以通過求解對應動力系統的穩定平衡點來得到。因此，動力系統的構造和性質是該方法的關鍵。首先介紹非線性動力系統的一些基本概念。

2.1 非線性動力學系統概述

非線性自治動力學系統可以表示為

式中：F:?n→?m為從?n空間到?m空間的映射；為狀態變量x的一階導數。定義動力學系統的軌跡為t=0時刻起從x0出發的曲線，記作

下面介紹平衡流形和穩定平衡流形的定義，以及相關的定理。

定義1 平衡流形：對于動力學系統式（6），稱Σ={F (x)≡0|x∈?n}為該動力學系統的一個平衡流形。即平衡流形Σ上的任意一點x都滿足F(x)=0。

定義2 穩定平衡流形：若對于任意ε＞0，存在δ=δ(ε)＞0 ，使對于任意的 x∈Bδ(Σ )，有 Φ(t,x)∈Bε(Σ )，其中 Bδ(Σ)={x ∈ ?n|‖x-y‖＜δ,?y∈Σ }，則動力學系統式（6）的平衡流形Σ是漸近穩定的，稱作穩定平衡流形SEM（stable equilibrium mani?fold）。反之，平衡流形Σ是不穩定的，稱作不穩定平衡流形UEM（unstable equilibrium manifold）。

定理1 設動力學系統的平衡流形Σ和任意一點x∈Σ，DF(x)為F在x處的雅可比矩陣。若DF(x)在正規空間Nx(Σ)上的特征向量（所有不在切空間上的特征向量）所對應的特征值的實部都為負值，則x是穩定的，Σ為穩定平衡流形。若DF(x)在正規空間Nx(Σ)上的特征向量所對應的特征值至少有一個含有正實部，則x是不穩定的，Σ是不穩定平衡流形。如果有k個含有正實部的特征值，稱Σ為k型不穩定平衡流形（type-k UEM）。

2.2 理論基礎

考慮如下QGS[7]：

式中，DH(x)為 H(x)的雅可比矩陣，H(x)=[z -h(x)]。

定理2 動力學系統式（7）的一個穩定平衡流形 Σs是函數E(x)=‖H (x)‖22的一個局部最小點。

證明：

設Σs是動力學系統式（7）的穩定平衡流形，因此存在一個 ε＞0使得對任意的y∈Bε(Σs)，其軌跡Φ(t,y)→Σs。因為函數 E(x)=‖H(x)‖22 沿著軌跡的導數 E(x)=-‖DH (x)TH(x)‖2≤0是非正的，也就是說對于任意的y∈Bε(Σs)，E(x)是一個非增函數，即Σs是下面最小化問題的一個局部最優解：

定義3 退化穩定平衡流形DSEM（degenerat?ed stable equilibrium manifold）：動力學系統式（7）的一 個 穩 定 平 衡 流 形 Σs，若 H(Σs)≠0 且DH(Σs)TH(Σs)=0，則稱 Σs為退化穩定平衡流形。

由于測量誤差的存在以及測量在局部具有一定的冗余度，根據定理2和定義3，動力學系統將只存在退化的穩定平衡流形，并且該退化的穩定平衡流形是函數E(x)=‖H (x)‖22的一個非零局部最小點。這里2E(x)=Hs(x)THs(x)即為最小二乘準則的目標函數。若將Hs(x)中的約束乘以權值其中測量i的標準差σi反映了對應測量的精確度，則有E(x)=Hs(x)TWHs(x)，這時WLS的解便為QGS的退化穩定平衡流形。因此，可以通過對QGS進行軌跡追蹤，得到其退化穩定平衡流形上一點，從而得到WLS模型的狀態估計解。

下面討論DSEM的維數與狀態估計中不可觀測狀態變量維數之間的關系。設集合ΣDSEM為QGS的退化穩定平衡流形，如果ΣDSEM中的每個點對應的H(x)值不同，那么集合ΣDSEM在?n中的Lebesgue測度為0。定理的證明見文獻[8]。

當狀態估計問題中系統的網絡結構和測量能夠保證系統的可觀測性時，定理3中的假設是滿足的。因此對應QGS的退化穩定平衡流形，即WLS的解是孤立的點。然而，當狀態估計問題中出現不可觀測變量時，由于不可觀狀態量的取值對WLS目標函數值沒有影響，因此定理3中的假設不能滿足，此時給出下面的定理。

定理4 設系統中含有k個不可觀狀態變量，那么QGS動力學系統的退化穩定平衡流形ΣDSEM的維數為k。

證明：

設x為退化穩定平衡流形ΣDSEM上一點，因此x=QH(x)=-DH(x)TH(x)=0。設系統中有k維不可觀測狀態變量，由于不可觀狀態量的取值對QH(x)的值沒有影響，則存在映射Ρ:U→?n，其中U是?k中的一個開集，滿足QH[ ]Ρ(s)≡0，s∈U。因此ΣDSEM={Ρ (s)|s∈U }為一個k維平衡流形。證畢。

定理4說明當狀態問題中含有k維不可觀測變量時，對應QGS的退化穩定平衡流形也是k維的，在k維退化穩定平衡流形上的任意一點都可以使WLS目標函數最小，因此WLS不具有唯一解。這時QGS收斂到退化穩定平衡流形上的哪一點將和初始值有關。

2.3 算法流程

根據上面闡述的理論基礎，本節給出基于動力學系統算法求解WLS模型的求解流程，當傳統基于高斯牛頓的WLS算法不收斂時，啟動基于QGS的狀態估計算法，通過對QGS的軌跡進行追蹤，得到其退化穩定平衡流形，從而得到WLS估計模型的狀態估計解。

步驟1 讀入電力系統網絡參數和測量數據；

步驟2 生成狀態估計測量方程；零注入平衡方程；

步驟3 采用高斯牛頓法求解WLS狀態估計問題，如果收斂，輸出狀態估計結果；否則進入步驟4；

步驟4 并構造對應的QGS如式（7）所示；

步驟5 對QGS的軌跡進行跟蹤，采用數值積分方法，求解到QGS的退化穩定平衡流形，收斂判據為‖QH(x)‖＜ε

求解過程中使用的數值積分方法可以采用隱式歐拉法、龍格庫塔法等常用的數值方法進行軌跡的追蹤。仿真軟件Matlab中也提供了多種積分器可供選擇，本文的測試算例均采用隱式剛性積分器ode15s。

3 算例分析

3.1 IEEE 14節點測試系統少量測算例

本節給出在IEEE14節點測試系統上進行的仿真結果，圖1給出了IEEE14節點測試系統的接線圖，以及本算例使用的測量類型和位置。該測試系統狀態變量維數為n=2×14-1=27，測量個數為m=36，零注入約束3個，分別是7號母線的有功和無功零注入約束，8號母線的有功零注入約束。所有的測量數據均由已知的潮流解推得，并疊加標準差σ=0.02的隨機誤差。

圖1 IEEE 14節點系統單線圖及其測量配置情況Fig.1 Single-line diagram of IEEE 14-bus test system and its measurement arrangement

該系統存在不可觀測狀態變量Va10（10號母線的電壓相角），這時WLS算法中的增益矩陣奇異，因此傳統基于高斯-牛頓法的WLS算法無法計算。按照如上介紹的基于動力學系統的狀態估計算法流程，可以求解到QGS的穩定平衡流形，從而得到可觀測變量的狀態估計解。

初始點設置為平啟動（電壓幅值為1.0，相角為0°），得到測量殘差、狀態變量估計誤差和收斂過程如圖2所示。圖2（a）為QGS的收斂過程，橫坐標是指積分運算中的時間參數，不是實際計算時間。縱坐 標 為QGS方 程 ‖QH(x)‖，收 斂 判 據 為‖QH(x)‖ ＜10-6。圖2（b）為測量殘差，圖2（c）為估計誤差，可以看到除不可觀變量外，其他狀態變量均有較高的估計精度，估計誤差均值為D=0.004 175。本算例中Va10的實際值為-0.264 rad，不可觀測變量的估計值和初始點選擇有關。

圖2 在測量含有標準差為0.02的隨機誤差時的收斂過程和測量殘差以及估計誤差Fig.2 Convergence process,measurement residuals and estimation errors when measuring random errors with standard deviation of 0.02

根據定理4，該算例含有1維不可觀測變量，則QGS的退化穩定平衡流形也為1維，下面在Va9×Va10×Vm8的3維空間中展示QGS的退化穩定平衡流形以及不同初始點的收斂軌跡。

如圖3所示，在Va9×Va10×Vm8的3維空間中隨機采樣20個不同的初始點，通過對QGS進行軌跡追蹤，均收斂到QGS的退化穩定平衡流形上，圖3中與Y軸Va10平行的一維曲線。不同的初始值會收斂到退化穩定平衡流形上的不同點。由于不可觀變量Va10的取值對WLS目標函數沒有影響，因此一維退化穩定平衡流形上的任意一點都可以使WLS目標函數最小。圖3中曲線的灰度表示WLS目標函數，數值越小，黑色越深。

圖3 IEEE14節點不可觀測測試算例中QGS的退化穩定平衡流形以及不同初始點的收斂軌跡Fig.3 DSEM of QGS for IEEE 14-bus unobservable test system and the convergence trajectory with different initial points

4.2 IEEE 118節點少量測測試算例

下面使用IEEE 118測試系統對基于動力學系統的狀態估計算法進行測試。系統中含有54個發電機，186條線路和91個帶有P-Q負荷的母線。本算例使用完整的電氣量測量一共有1 098個，包括118個電壓幅值測量，236個功率注入測量和744個線路首末端功率測量。本節使用的測量均在真實值基礎上疊加了標準差為0.05的隨機誤差。

從全部1 098個量測數據開始，隨機地去掉量測直到只剩1個量測，對基于動力學系統的狀態估計算法進行測試。顯然，當逐個去掉測量時，系統測量的冗余度下降，會有部分測量變為關鍵量測，如果繼續隨機地去掉測量，系統會出現不可觀狀態變量。圖4給出了測量數目逐漸減少時，基于動力學系統算法和傳統高斯-牛頓法的估計精度對比。當測量數目小于460個時，由于增益矩陣的奇異，傳統WLS算法無法獲得狀態估計解，而本文提出的基于動力學系統的算法能夠在系統失去可觀測性時提供可以接受的估計值。當然，隨著測量的減少，不可觀狀態變量的增多，估計誤差也迅速變大。

當測量數目減少為300個時，經可觀測性分析，系統中將有4個不可觀測狀態變量分別為Va112、Va20、Va38、Vm112。這時WLS算法的增益矩陣奇異，條件數為cond(G)=7.143 9×1019，傳統高斯牛頓法不能收斂。圖5給出了此時增益矩陣G的特征值。圖6給出了QGS的收斂過程。

圖4 隨著測量數目減少，基于動力學系統算法和傳統高斯-牛頓法的估計精度對比Fig.4 Comparison of estimation accuracy between the method based on dynamical system and the traditional Gauss-Newton method when the number of measurements is decreasing

圖5 在測量數目減少為300個時，WLS算法中增益矩陣G的特征值Fig.5 Eigenvalue of the gain matrix G in WLS method when the number of measurements is reduced to 300

圖6 在測量數目減少為300個時，基于動力學系統算法的收斂曲線Fig.6 Convergence curve of the method based on dynamical system when the number of measurements is reduced to 300

3.3 IEEE 14節點存在拓撲錯誤時的收斂性測試

采用文獻[5]中的測試算例，假設6-12線路被認為斷開，而實際上這條線路沒有斷開。該算例有42個測量，包括7個電壓幅值測量，11個注入功率測量，24個線路功率測量。

圖7給出了IEEE 14節點測試算例含有一個拓撲錯誤時，高斯-牛頓算法和基于QGS動力學系統算法的收斂曲線對比。上方的橫坐標表示高斯-牛頓法的迭代次數，下方的橫坐標表示對QGS動力學系統的積分時間參數，可以看出高斯-牛頓法呈現出震蕩的現象，未能求解到WLS解。而本文提出的算法可以收斂。

圖7 IEEE 14節點測試算例含有一個拓撲錯誤時，高斯牛頓法和動力學系統算法的收斂曲線對比Fig.7 Comparison of convergence curve between the Gauss-Newton method and the method based on dynamical system for IEEE 14-bus test system with one single topological error

牛頓（或高斯-牛頓）法在初始點接近解時具有優良的收斂速度，但缺點在于對初始點很敏感，因此當狀態估計問題中出現很大的測量殘差時可能造成求解失敗。當測量數據中含有嚴重的壞數據或者拓撲錯誤時，往往需要先求解出一個WLS解，然后對殘差或標準化殘差進行分析，從而對壞數據和錯誤拓撲進行識別。然而從本算例可以看出，在出現嚴重錯誤時，可能造成牛頓法不收斂，從而無法進行下一步的壞數據和拓撲錯誤識別。本文提出的算法在這種情況下能夠保證收斂，為錯誤測量和拓撲識別提供基礎。

4 結語

本文提出了一種新的基于動力學系統的方法來求解狀態估計問題，構造商QGS，其退化的穩定平衡流形對應于WLS的解，并證明QGS是漸進穩定的，因此可以保證算法的收斂性。該方法可以可靠地處理因拓撲錯誤而可能引起的不收斂問題，當系統中存在不可狀態變量時，可以繼續進行求解并得到有效的狀態估計解，同時用數值算例對不可觀算例的WLS解進行了刻畫。所提出的算法在IEEE 14節點和IEEE 118節點系統上進行測試，驗證了其可靠性和有效性。