基于平方和方法的電力系統(tǒng)穩(wěn)定域估計

摘" 要： 目前對電力系統(tǒng)穩(wěn)定性分析主要使用小信號分析方法，但當大信號擾動發(fā)生時，電力系統(tǒng)內(nèi)在的非線性無法忽略，小信號分析方法的有效性值得商榷。在大信號擾動時，構造合理的暫態(tài)能量函數(shù)并確定系統(tǒng)臨界穩(wěn)定對應的函數(shù)值（臨界能量），從而確定穩(wěn)定域并非易事。為此，提出一種基于多項式表達的平方和方法來估計電力系統(tǒng)的穩(wěn)定域。該方法將電力系統(tǒng)通過多項式近似的方式表達，利用平方和方法將李雅譜諾夫函數(shù)的構造問題轉化為線性矩陣不等式的可行解問題，通過數(shù)值計算得到構造的李雅譜諾夫函數(shù)，再通過優(yōu)化計算得到穩(wěn)定域的邊界。最后采用數(shù)值算例驗證了所提估計方法的可行性和有效性。

關鍵詞： 電力系統(tǒng)； 穩(wěn)定域； 平方和方法； 多項式近似； 李雅譜諾夫函數(shù)； 線性矩陣不等式； 算例驗證

中圖分類號： TN911.1?34； TM712" " " " " " " " " 文獻標識碼： A" " " " " " " " " " 文章編號： 1004?373X（2025）04?0097?05

Power system ROA estimation based on SOS

ZHANG Weiwei， LI Zhaoming， SHI Hongtao， GAO Feng， ZHANG Bai

（School of Electrical and Information Engineering， North Minzu University， Yinchuan 750021， China）

Abstract： The existing research on stability analysis of power system mostly focuses on small?signal analysis methods. However， when large?signal disturbances occur， the intrinsic nonlinearity of the power system becomes inevitable， so that the small?signal analysis method may be no longer valid. How to construct a reasonable transient energy function and determine the function value （critical energy） corresponding to the critical stability of the system to determine the region?of?attraction （ROA） is not easy when large?signal disturbances occur. The sum of squares （SOS） method based on polynomial expression is proposed to estimate ROA of power system. In this method， the power system can be expressed by means of polynomial approximation. The construction problem of Lyapunov function is transformed into the feasible solution problem of linear matrix inequalities by means of SOS method. The constructed Lyapunov function is obtained by means of numerical calculation， and the boundary of ROA is obtained by means of the optimization calculation. The numerical examples verify the feasibility and effectiveness of the proposed method.

Keywords： power system; region?of?attraction; sum of squares method; polynomial approximation; Lyapunov function; linear matrix inequality; example verification

0" 引" 言

伴隨新能源的迅猛發(fā)展，高比例可再生能源經(jīng)電力電子接口匯集并網(wǎng)，新能源的波動性與電力電子的混雜控制改變了電力系統(tǒng)的動態(tài)特性，經(jīng)典穩(wěn)定性定義和分析方法亟需擴展[1?2]。

目前針對電力系統(tǒng)穩(wěn)定的研究主要借鑒經(jīng)典的暫態(tài)穩(wěn)定性分析方法，即等面積法和能量函數(shù)法[3?4]。但是，上述兩種方法均忽略了電力系統(tǒng)數(shù)學模型中存在的非線性項，導致穩(wěn)定性判定結果既可能保守又可能冒進。同時，時域仿真法存在計算速度慢、不能給出穩(wěn)定裕度的缺點。

李雅普諾夫直接法在定量分析非線性系統(tǒng)的大范圍穩(wěn)定性方面具有優(yōu)勢，通過構造暫態(tài)能量函數(shù)，從系統(tǒng)能量的角度出發(fā)，可以在不計算整個系統(tǒng)運動軌跡的前提下進行穩(wěn)定分析和判斷，得到了長足的發(fā)展[5?9]。首次積分法[2，7]、Zubov方法[8]、模態(tài)解耦方法[9]等常用來構造能量函數(shù)。但是對于電力系統(tǒng)而言，構造合理的暫態(tài)能量函數(shù)并確定系統(tǒng)臨界穩(wěn)定對應的函數(shù)值（臨界能量），從而確定穩(wěn)定域并非易事。一方面，由于電力系統(tǒng)的復雜性，建立的數(shù)學模型通常都簡化為二階經(jīng)典模型，采用候選的二次李雅普諾夫函數(shù)，對穩(wěn)定域的估計通常偏于保守；當電力系統(tǒng)模型階次較高或者非線性較復雜時，很難找到合適的候選李雅普諾夫函數(shù)；另一方面，電力系統(tǒng)的平衡點隨著潮流遷移，通過坐標平移將平衡點平移到零平衡點很難進行分析，通常假設其平衡點變化很小，只在特定的平衡點的鄰域內(nèi)分析臨界能量，這會導致臨界能量偏小時得到的結果偏于保守，臨界能量偏大時得到的結果偏于樂觀。

為更準確地估計系統(tǒng)穩(wěn)定域，本文在對電力系統(tǒng)非線性模型進行多項式表達的基礎上，基于平方和方法構造高次的候選李雅普諾夫函數(shù)，通過數(shù)值計算方法嘗試構造暫態(tài)能量函數(shù)，將尋找臨界能量問題轉化為能量函數(shù)的最值問題，進而求得系統(tǒng)的穩(wěn)定域，并對單機無窮大系統(tǒng)進行仿真驗證，同時與首次積分法和Zubov方法進行比較。仿真結果表明，通過所提方法求得的穩(wěn)定域可逼近真實的穩(wěn)定域，實現(xiàn)了對穩(wěn)定域的近似估計，表明了該方法的可行性。

1" 理論基礎

對于自治非線性系統(tǒng)：

[xt=fxt] （1）

式中：[xt∈Rn]是狀態(tài)向量；[f：Rn→Rn]是局部Lipschitz映射，且[f0=0]。

設[?ζ，t]是系統(tǒng)在初始條件[x0=ζ]下[t]時刻的解，則對應系統(tǒng)平衡點[x=0]的吸引域為[ROA0：=ζ∈Rn：limt→∞?ζ，t=0]，如果系統(tǒng)的解[?ζ，t∈M，?t≥0]，且[ζ∈M]，則稱集合[M]為系統(tǒng)的不變集[10]。

引理[11]設[γgt;0]并存在連續(xù)可微函數(shù)[V：Rn→R]滿足：

[Λ：=x∈Rn：V（x）lt;γ有界] （2）

[V（x）gt;0，V（0）=0，" for x≠0] （3）

[Λ＼{0}?x∈Rn：?V（x）f（x）lt;0] （4）

則對于所有的[ζ∈Λ]，系統(tǒng)的解存在且滿足[?ζ，t∈Λ]，[?t≥0]，且[limt→∞?ζ，t=0]，即[Λ]是系統(tǒng)ROA的一個不變子集。

對一般系統(tǒng)而言，精確計算ROA甚至對ROA進行估計是不容易做到的。因此，本文基于引理，求解系統(tǒng)ROA內(nèi)的最大不變集來作為對ROA的最優(yōu)估計，也即穩(wěn)定域的最大估計。

針對任意自治非線性（非多項式）系統(tǒng)（見公式（1）），在滿足一定精度條件下，可重塑為多項式系統(tǒng)[12]。

對于[n]元多項式[p（x）]，如果能分解為多項式[f1x，f2x，…，fmx]的平方和，即：

[p（x）=i=1mf2ix] （5）

則稱[p（x）]為平方和多項式，用[Σn]表示所有平方和多項式的集合。顯然，如果多項式能夠寫成平方和的形式，則該多項式全局半正定。因此可以通過驗證多項式是否存在平方和分解代替半正定條件[13]，而驗證多項式是否存在平方和分解在計算上更容易處理。

將常用的二次型多項式[xTQx]推廣到高階多項式的形式，表示為：

[p（x）=zT（x）Qz（x）] （6）

式中：[Q]是對稱半正定矩陣（Gram矩陣）；[z（x）]表示由次數(shù)不大于[deg（p（x））2]的單項式構成的向量。平方和分解等價于存在[z（x）]和使得[Q]對于式（6）成立[14]。

因此，可以先利用高階多項式構造候選的李雅普諾夫函數(shù)，即可以通過平方和分解加入引理中的條件（見式（4）），構造滿足引理中條件[V（x）gt;0]的李雅普諾夫函數(shù)。在得到李雅普諾夫函數(shù)的基礎上可以通過計算其最大值確定穩(wěn)定邊界[γ]，從而得到穩(wěn)定域的最大估計。

[n]元多項式的非負性可以通過判斷是否存在平方和分解來得出，只是充分條件，不是必要條件。因此，通過平方和方法判定多項式正定具有一定保守性。

2" 算法分析

2.1" 電力系統(tǒng)的多項式近似

一般非線性系統(tǒng)的多項式近似表達，在足夠高的近似階數(shù)下，多項式近似系統(tǒng)與原系統(tǒng)平衡點可任意接近，且其不穩(wěn)定平衡點類型可保持不變[15]。因此，可以用多項式近似系統(tǒng)來逼近原系統(tǒng)的穩(wěn)定域邊界。

本文采用泰勒級數(shù)展開的方法將電力系統(tǒng)中的非線性部分在工作點附近展開，并取一定階次的近似，忽略較高階次的等價無窮小。具體階次的取舍與候選李雅普諾夫函數(shù)的階次相關，這在后面的算例分析中將進一步討論。

2.2" 李雅普諾夫函數(shù)的構造

平方和分解可以代替半正定條件，為滿足引理中[?V（x）f（x）lt;0]的條件，由Positivstellensatz定理[12，16]和S?procedure方法[17]引入給定的正定多項式[li（x）gt;0]（可以取[li（x）=10-6xTx]），[i=1，2]，使得候選的李雅普諾夫函數(shù)滿足：

[V（x）-l1（x）∈Σn] （7）

[-?V（x）f（x）+l2（x）∈Σn] （8）

從上面的分析中可以知道SOS 問題的求解是先要將多項式[p（x）]轉換成式[p（x）]中單項式和Gram矩陣相乘的形式，即向量[z（x）]和正定矩陣[Q]的求解問題，這一步稱為 SOS 分解（SOS decomposition， SOSP）[18]。SOS分解可轉化為半正定規(guī)劃問題（SDP），在確定[z（x）]組成單項式的基礎上，即SOS分解的LMI可行性問題。為求得[V（x）]的可行解，做如下處理。

直接指定[V（x）]的階次為偶數(shù)次，組成[z（x）]的單項式最高階次數(shù)[deg（V（x））2]，例如，當系統(tǒng)變量為2維時，指定[deg（V（x1，x2））=4]時，[z（x）= x1 x2 x21 x22x1x2]。

確定[z（x）]后，求解矩陣[Q]的過程等價于下面LMI問題：

[Q（μ）=Q0+iμiQi≥0] " " （9）

式中[μi]表示[Q]的自由變量，不管它們?nèi)『沃?，都得到相同形式的多項式[p（x）]。

該LMI的求解過程可通過SEDUMI、SDPT3等工具實現(xiàn)。LMI 有可行解時，將[Q]的可行解代入式（6），即可得到李雅普諾夫函數(shù)。

1） [z（x）=x]時，[V（x）]退化為二次型的李雅普諾夫函數(shù)，二次型李雅普諾夫函數(shù)是多項式李雅普諾夫函數(shù)階次為2的情形，多項式李雅普諾夫函數(shù)可以看作二次型李雅普諾夫函數(shù)的推廣。

2） [z（x）]的選取不是唯一的，[Q]也不是唯一的。這為李雅普諾夫函數(shù)的構造提供了更多的自由度。系統(tǒng)無法構造二次型李雅普諾夫函數(shù)時，可以嘗試構造多項式李雅普諾夫函數(shù)。

3） 多項式李雅普諾夫函數(shù)的階次可以指定為任意的偶數(shù)，這為SOSP的求解提供了新的自由度，當李雅普諾夫函數(shù)的階次較低，LMI無法求得可行解時，可以嘗試高階李雅普諾夫函數(shù)。

2.3" ROA的估計

在得到李雅普諾夫函數(shù)的基礎上，可以通過求解如下優(yōu)化問題計算其最大值確定穩(wěn)定邊界[γ]，從而得到ROA的最大估計。

[max γs.t." s（x）∈Σns（x）（V（x）-γ）∈Σn] （10）

式中[s（x）]為引入的決策變量。

若考慮[V（x）]是待定的多項式，將待求的不變子集[ΛV，γ]也加入到優(yōu)化設計過程中，可在[ΛV，γ]內(nèi)定義一個區(qū)域[Λh，β]。

[Λh，β?ΛV，γΛh，β：=x∈Rn|h（x）≤β] （11）

式中[h（x）]為已知的二次型多項式，[h（x）=xTRx]。

式（11）的作用是用已知的二次型來定義一個內(nèi)接于[ΛV，γ]的區(qū)域[Λh，β]，求解[β]的最大值以得到盡可能大的ROA的不變子集。由S?procedure方法轉化為優(yōu)化問題：

[max βs.t." s1（x）∈Σns2（x）∈ΣnV（x）-l1（x）∈Σn-（β-h（x））s1（x）+（V（x）-γ）∈Σn-?V（x）f（x）+l2（x）+s2（x）（V（x）-γ）∈Σn] （12）

式中：[s1（x）]、[s2（x）]為引入的決策變量。

至此，利用多項式近似系統(tǒng)逼近原系統(tǒng)的穩(wěn)定域，基于平方和方法構造多項式李雅普諾夫函數(shù)，并通過優(yōu)化問題計算其最大值確定穩(wěn)定邊界，可實現(xiàn)電力系統(tǒng)ROA的最大估計。

3" 算例分析

考慮單機無窮大系統(tǒng)，其數(shù)學描述如下：

[δ=ωω=Pmaxωs2H（sinδs-sin（δ+δs））-D2Hω] （13）

式中：[δ]和[ω]分別表示發(fā)電機轉子功角和轉速；[D]和[H]分別為阻尼系數(shù)和轉動慣量常數(shù)；[ωs]和[δs]為同步轉速和穩(wěn)定功角；[Pmax]為最大傳輸功率。相關參數(shù)取值[19]如下：[δs=15°]，[ωs=2π×60 rad/s]，[D=1]，[H=3 s]，[Pmax=1.7 p.u.]。

分別取系統(tǒng)在平衡點（SEP）處3～9階的多項式近似，通過第2節(jié)的方法可求得不同階次多項式近似下的李雅普諾夫函數(shù)。為了和文中的結果比較，這里給出系統(tǒng)三階多項式近似下求得的李雅普諾夫函數(shù)（指定李雅普諾夫函數(shù)的階次為6），其余多項式近似系統(tǒng)的李雅普諾夫函數(shù)類似，限于篇幅，這里不一一列舉。

[V（6）（δ，ω）=4.4δ6-0.380 6δ5ω+4.655δ5+" " " " " " " " " " " "0.087 92δ4ω2-0.534 7δ4ω-108.0δ4-" " " " " " " " " " " "0.008 743δ3ω3+0.154 5δ3ω2+2.238δ3ω-" " " " " " " " " " " "76.13δ3+0.000 501δ2ω4-0.015 2δ2ω3-" " " " " " " " " " " "1.065δ2ω2+1.029δ2ω+849.3δ2-" " " " " " " " " " " "5.519×10-5δω5+0.000 810 9δω4+" " " " " " " " " " " "0.022 36δω3-0.197 6δω2+0.986 1δω-" " " " " " " " " " " "0.154 5δ+1.877×10-6ω6-8.567×10-5ω5-" " " " " " " " " " " "0.004 863ω4+0.015 32ω3+8.144ω2-" " " " " " " " " " " "6.523×10-5ω] （14）

為比較多項式近似系統(tǒng)的穩(wěn)定域和原系統(tǒng)穩(wěn)定域的關系，通過逐點掃描得到原系統(tǒng)的穩(wěn)定域，圖1和圖2分別給出了3階和9階多項式近似系統(tǒng)和原系統(tǒng)的穩(wěn)定域的比較圖。從圖可以看出，高階多項式近似系統(tǒng)可以更好地逼近原系統(tǒng)的穩(wěn)定域，求得的穩(wěn)定域可以作為原系統(tǒng)穩(wěn)定域的近似估計。圖3為3～9階多項式近似系統(tǒng)的穩(wěn)定域。從圖中可以看出，隨著多項式近似系統(tǒng)階次的增加，多項式近似系統(tǒng)的穩(wěn)定域有增大的趨勢。

與文中采用首次積分法和Zubov方法給出的李雅普諾夫函數(shù)求得的穩(wěn)定域比較，首次積分法構造的李雅普諾夫函數(shù)為：

[Vint（δ，ω）=0.5ω2+51.59δ2-4.608δ3-4.299δ4] （15）

Zubov方法構造的李雅普諾夫函數(shù)為：

[VZubov（δ，ω）=6.291×10-4ω2+9.692×10-6ωδ+" " " " " " " " " " " " " 0.064 91δ2+8.386×10-9ω3+" " " " " " " " " " " " " 4.193×10-9ω2δ+1.299×10-6ωδ2-" " " " " " " " " " " " " 5.797×10-3δ3-1.771×10-7ω4+" " " " " " " " " " " " " 7.75×10-9ω3δ-3.654×10-5ω2δ2+" " " " " " " " " " " " " 1.16×10-6ωδ3-7.294×10-3δ4+" " " " " " " " " " " " " 3.172×10-11ω5+2.811×10-11ω4δ+" " " " " " " " " " " " " 8.192×10-9ω3δ2+3.269×10-6ω2δ3+" " " " " " " " " " " " " 4.281×10-7ωδ4+3.369×10-4δ5] （16）

利用不穩(wěn)定平衡點計算穩(wěn)定域邊界的方法[19]，分別求得[Vint（δ，ω）=101.3]，[VZubov（δ，ω）=0.105 2]。采用這兩種方法分別得到的穩(wěn)定域和3階、9階多項式近似系統(tǒng)得到的穩(wěn)定域邊界比較，結果如圖4所示。

從圖4中可以看出，Zubov方法得到的穩(wěn)定域最保守，首次積分法得到的穩(wěn)定域和3階多項式近似得到的穩(wěn)定域大小基本一致，但都比9階多項式近似系統(tǒng)的穩(wěn)定域小。因此，本文采用高階多項式近似系統(tǒng)得到的穩(wěn)定域更寬松，但與原系統(tǒng)的穩(wěn)定域相比也偏于保守，這是由于李雅普諾夫直接法分析系統(tǒng)穩(wěn)定的條件是充分條件，通過平方和方法判定多項式正定也是充分條件，得到的結果通常偏于保守。此外，在[Λδ，ω=δ∈[2，4]，ω∈[-14.8，0]]的區(qū)域內(nèi)，也存在較大的偏差。這是由于原系統(tǒng)真實的穩(wěn)定域是一個非凸集合，采用SOS方法求解得到的穩(wěn)定域估計是凸集所致。

4" 結" 論

本文針對電力系統(tǒng)的穩(wěn)定域估計，通過對原系統(tǒng)的多項式近似，提出了一種基于平方和優(yōu)化的穩(wěn)定域估計方法。該方法將李雅普諾夫函數(shù)的構造問題轉化為LMI的可行解問題，通過數(shù)值計算得到構造的李雅普諾夫函數(shù)，通過優(yōu)化計算得到穩(wěn)定域的邊界值，相比首次積分法和Zubov方法，得到的穩(wěn)定域估計能更好的逼近原系統(tǒng)的穩(wěn)定域。但是，相比原系統(tǒng)真實的穩(wěn)定域，該方法同樣存在一定的保守性，如何加入寬松變量改變方法的保守性，還需要進一步的研究。此外，針對建立的高階電力系統(tǒng)模型，在不簡化到低階模型的前提下實現(xiàn)穩(wěn)定域的估計，也是下一步需要研究的內(nèi)容。

注：本文通訊作者為李兆明。

參考文獻

[1] 姜齊榮，趙崇濱.并網(wǎng)逆變器的電磁暫態(tài)同步穩(wěn)定問題[J].清華大學學報（自然科學版），2021，61（5）：415?428.

[2] GRIGSBY Leonard. Power system stability and control [M]. Boston： CRC Press， 2012.

[3] 張宇，張琛，蔡旭，等.并網(wǎng)變換器的暫態(tài)同步穩(wěn)定性分析：穩(wěn)定域估計與鎮(zhèn)定控制[J].中國電機工程學報，2022，42（21）：7871?7884.

[4] KABALAN M， SINGH P， NIEBUR D. Large signal Lyapunov?based stability studies in microgrids： a review [J]. IEEE transactions on smart grid， 2017， 8（5）： 2287?2295.

[5] ZHANG W， DU Z， HAN D. Transient stability domain estimation of AC/DC systems considering HVDC switching characteristics based on the polynomial Lyapunov function method [J]. International journal of electrical power amp; energy systems， 2021， 129（4）： 106875.

[6] AWREJCEWICZ J， BILICHENKO D， CHEIB A K， et al. Estimating the region of attraction based on a polynomial Lyapunov function [J]. Applied mathematical modelling， 2021， 90： 1143?1152.

[7] PAI M A. Energy function analysis for power system stability [D]. Boston， MA， USA： Kluwer Academic Publishers， 1989.

[8] CHIANG H D， WU F F. Foundations of the potential energy boundary surface method for power system transient stability analysis [J]. Circuits amp; systems IEEE transactions on， 1988， 35（6）： 712?728.

[9] WANG B， SUN K， XU X. Nonlinear modal decoupling based power system transient stability analysis [J]. IEEE transactions on power systems， 2019， 34（6）： 4889?4899.

[10] KHALIL H K. Nonlinear systems. 3rd [M]. New Jersey， USA： Prentice?Hall， 1996.

[11] TAN W， PACKARD A. Stability region analysis using polynomial and composite polynomial Lyapunov functions and sum?of?squares programming [J]. IEEE transactions on automatic control， 2008， 53（2）： 565?571.

[12] PAPACHRISTODOULOU A， PRAJNA S. On the construction of Lyapunov functions using the sum of squares decomposition [C]// Proceedings of the 41st IEEE Conference on， Decision and Control. Las Vegas， NV， USA： IEEE， 2002： 3482?3487.

[13] SHOR N Z. Class of global minimum bounds of polynomial functions [J]. Cybernetics， 1987， 26（3）： 731?734.

[14] PRAJNAL S， PAPACHRISTODOULOUL A， WU F. Nonlinear control synthesis by sum of squares optimization： a Lyapunov?based approach [C]// Proceeding of the Asian Control Conference. Melbourne， VIC， Australia： IEEE， 2004： 157?165.

[15] 孫玉嬌，劉鋒，梅生偉.非線性系統(tǒng)的多項式近似表示及電力系統(tǒng)應用（Ⅰ）：理論篇[J].電機與控制學報，2010，14（8）：19?23.

[16] PARRILO P A. Structured semidefinite programs and semialgebraic geometry methods in robustness and optimization [D]. Pasadena： California Institute of Technology， 2000.

[17] BOYD S， EL GHAOUI L， FERON E， et al. Linear matrix inequalities in system and control theory [M]. Philadelphia， PA： SIAM， 1994.

[18] PRAJNA S， SEILER P， PARRILO P A. SOSTOOLS： sum of squares optimization toolbox for MATLAB [EB/OL]. [2023?04?17]. https：//optimization?online.org/2002/05/483/.

[19] WANG B， KAI S. Formulation and characterization of power system electromechanical oscillations [J]. IEEE transactions on power systems， 2016， 31（6）： 5082?5093.

作者簡介：張巍?。?983—），男，寧夏固原人，博士研究生，高級工程師，研究方向為電力系統(tǒng)穩(wěn)定性分析和仿真。

李兆明（2000—），男，回族，寧夏固原人，碩士研究生，主要研究方向為新型電力系統(tǒng)的穩(wěn)定與控制。

師洪濤（1984—），男，河北保定人，博士研究生，副教授，主要研究方向為分布式發(fā)電系統(tǒng)控制與并網(wǎng)、風電功率預測等新能源并網(wǎng)。

高" 峰（1983—），男，陜西府谷人，博士研究生，高級工程師，主要研究領域為電力系統(tǒng)仿真與分析。

張" 白（1981—），男，回族，內(nèi)蒙古烏拉特前旗人，博士研究生，副教授，主要研究方向為精密測量技術與儀器、人工智能算法與裝備等。