趙妍,柳旭,孫碩,朱建華,聶永輝
(1.東北電力大學輸變電技術學院,吉林省吉林市 132012;2.國網平湖供電公司,浙江省平湖市 314200;3.國網四平供電公司,吉林省四平市 136000;4.潤電能源科學技術有限公司,鄭州市 450000;5.東北電力大學教務處,吉林省吉林市 132012)
以新能源為主體的新型電力系統呈現高度電力電子化特征。新型電力系統具有發電單元多、單體容量小、交直流相互耦合等特點,這些特點使得傳統電網仿真分析技術的局限性日益凸顯。高比例新能源接入情況下,電力系統運行機理更加復雜,故障不確定性更強,傳統安全防御體系不具備靈活應對不確定性故障的能力,無法保障復雜故障條件下的電網安全運行,連鎖脫網事故多次發生。暫態仿真通過故障還原,提高應對不確定性故障的能力,一定程度上保障了復雜故障條件下的電網安全運行,支撐了構建清潔、低碳、安全、高效的能源體系。
但是,目前尚缺乏認識其機理的有效模型和仿真工具[1-6]。主要原因有:對新能源發電機技術研發持續投入、機組型號眾多,且不斷進行技術更新,通用化建模難度巨大。由于商業保密等原因,廠家提供的多為黑箱模型。風力發電機組的控制新技術快速發展、控制系統控制參數多、測試難度大、模型整定困難。模型的不準確使得電網運行方式更加保守,嚴重制約了新能源的消納。
實際上,若對大電力系統進行全面的電磁暫態仿真,往往耗時很長。而對電力系統進行暫態分析時,一般也只是關注其中的某一部分相關變量的暫態變化過程,因此,可以只對“感興趣”的部分網絡進行詳細的仿真,而對于其余部分網絡(文獻[7]中稱之為“外部系統”)進行建模等效簡化,可以加快暫態仿真的進程[7]。
電力系統外部系統電磁暫態等值方法分為兩類。一類方法是工頻下的ward外部系統等值,這類方法是通過對電力網絡進行星網變換而獲得。對于偽穩態運行的系統具有較高的等值精度。但是當系統發生不對稱故障或者設備飽和時,系統中諧波分量豐富,尤其是需要研究系諧波信號的影響時,這種等值方法的誤差比較大;而且該方法要求接口必須選在電壓或電流數值相對穩定的節點處,否則仿真計算時易出現振蕩[8-12]。另一類是與頻率相關的等值方法(frequency dependent network equivalents,FDNE)[13-19],這類方法是在一個寬頻段范圍內對外部系統進行等值,使得在選定的頻率范圍內,等值系統對研究系統的作用和影響近似或接近于外部系統對研究系統的作用和影響。文獻[13-18]提出了采用基于阻抗(導納)狀態空間表達式的矢量擬合(vector fitting, VF)等值算法。這類方法通過對外部系統進行頻率采樣,利用最小二乘法,求解系數矩陣為稠密矩陣的超定線性方程組,來獲得狀態空間表達式的相應系數矩陣。其等值網絡能夠在較寬頻段內對外部系統均具有較高的近似度,但計算量較大,且所求得的極點可能出現不穩定極點。文獻[19]提出了研究風電并網系統振蕩問題的頻域模式分析方法,使用外部系統阻抗網絡模型描述目標系統,通過求取網絡節點導納矩陣和回路阻抗矩陣的行列式零點,獲得系統的振蕩模式,并根據矩陣的左和右特征向量,計算出各節點和支路對振蕩模式的參與因子、可觀度和可控度指標。這種方法在一個較小頻率范圍內具有較好的等值效果;當頻帶較寬時,等值效果仍不夠理想。
與頻率相關的系統等值方法(例如文獻[13-19])的核心問題是需要獲得準確的狀態空間表達式。確定-隨機子空間(combined deterministic-stochastic identification,CDSI)在寬頻段范圍內是辨識狀態空間表達式的方法。CDSI是隨機子空間(stochastic subspace identification,SSI)算法的提升方法。SSI方法假定白噪聲輸入,僅利用輸出信號,辨識狀態空間表達式,在進行系統參數求解時無法建立確定的激勵-系統-響應模型,存在主導模態提取困難,引入虛假模態[20-22]的問題。CDSI能利用輸入(激勵)和輸出(響應)全部信息進行模態參數識別,建立確定的激勵-系統-響應模型。防止系統特征值解空間的擴大,可以解決密集模態識別困難、虛假模態較多等問題,同時保留了SSI計算量小、極點穩定等優點,提高了狀態空間表達式的辨識精度和等效模型的可靠性。
本文基于外部系統的輸入和輸出數據,利用CDSI對狀態空間表達式的系數矩陣進行辨識,根據零極點構建RLCG等值電路,進行建模等效簡化,實現了新型電力系統外部系統的等值構建。
1.1.1 CDSI原理
CDSI方法是將輸入、輸出數據組成特定的Hankel矩陣,計算Hankel矩陣的行空間投影,對投影進行奇異值分解,從而得到可觀測矩陣和狀態序列的卡爾曼濾波估計[23-25],進而確定系統矩陣A、B、C、D。確定-隨機系統可以表示為圖1。

圖1 確定-隨機系統Fig.1 Deterministic-stochastic identification
圖1中A為狀態矩陣,B為輸入矩陣,C為輸出矩陣,D為傳遞矩陣,u為輸入數據,y為輸出數據,x為狀態向量,w稱為過程噪聲,v稱為量測噪聲。下標k表示序列長度。
假設系統量測到的輸入U、輸出Y的時間間隔為Ts,U、Y分別為長度為k的序列,ui、yi為l個通道在i時刻的觀測向量。如果考慮上述隨機因素的干擾,那么就得到確定-隨機狀態方程:
(1)
設噪聲為白噪聲序列,且協方差矩陣如下:
(2)
式中:E(·)代表數學期望因子;δpq代表kronecker函數;Q、S、R為噪聲序列的協方差矩陣。這樣,CDSI方法的基本目標是:已知輸入數據u0,u1,…,uk和輸出數據y0,y1,…,yk,確定系統合適的階次n后,求得系統矩陣A、B、C、D、Q、R、S。下標p、q表示長度為k的序列中的任意項。
1.1.2 CDSI算法流程
CDSI算法的流程具體如下:
1)通過輸入、輸出數據構建Hankel矩陣[26]。
(3)
式中:U0|2i-1下標的第一個數字代表Hankel 矩陣左上角元素的時間系數,下標第二個數字代表Hankel矩陣左下角元素的時間系數,假設j→∞;將輸入數據Hankel矩陣U0|2i-1分為兩部分,上部分稱為“過去”輸入,用下標“p”表示,即Up;下部分稱為“將來”輸入,用下標“f”表示,即Uf。同理,輸出的Hankel矩陣Y0|2i-1也用同樣的方法定義。
定義輸入和輸出的Hankel矩陣的集合為W0|2i-1,則:
(4)
式中:Wp為“過去”集合;Wf為“將來”集合。
2)對Hankel矩陣進行QR分解:
(5)
式中:R11、R21、R22組成三角矩陣,Q1、Q2組成正交陣。
3)將輸出Yf投影到輸入和輸出Wp空間,由空間投影的性質得出行空間的正交投影Pi:
(6)
4)對投影矩陣Pi進行奇異值分解(singular value decomposition,SVD)確定系統階數:
Pi=USVT
(7)
式中:U和V為正交矩陣;S為降序排列的奇異值組成的對角矩陣。
理論上矩陣S主對角線中非0元素個數就是模型階數n。當受到噪聲污染時,非零奇異值的個數大于模型階數n,但數值都很小。
利用奇異值的相對變化率,即模型階次指標(model order indicator,MOI)作為奇異值突降指標判斷模型階數n。奇異值突降指標記作Moi,i。
σ1≥σ2≥…≥σn
(8)
(9)
式中:σi為降序排列的奇異值。
設初始系統階數為N,由于奇異值是降序排列的,在突降點大的地方,模型階次指標Moi,i也是一個大值,即Moi,i最大值對應的階次為模型的階次n。系統階數N=n+1(n為模型階數,1為噪聲階)。

(10)
(11)
(12)

6)使用卡爾曼濾波狀態序列和輸入輸出代入狀態空間方程求得系統狀態矩陣A、B、C、D和殘差矩陣pw、pv。
(13)
求系統狀態矩陣A的特征值λ(i=1, 2,…,n),進而可獲得頻率f′i。
用yk表示電流I;uk表示電壓U。狀態空間表達式Y(s)為:
(14)
式中:A、B、C、D為已知的系數矩陣;E為單位矩陣;I(s)、U(s)、Y(s)為頻域下輸入、輸出及對應關系,其中s是頻域下的自變量。
將式(14)轉化為頻域表達式,其中m≥n:
(15)
式中:an、bm為降冪排列項系數,下同。
pi表示實數極點;p′k±p″k表示復數極點,可將式(15)寫為:
(16)
裂項得:
(17)
在得到外部系統電流和電壓的傳遞函數后,可用3×3簡單網絡代替外部系統,則該網絡可作為外部系統的等值系統。并聯等值模型如圖2所示。

圖2 并聯等值模型Fig.2 Parallel equivalent model
在圖2中,實數極點對應RL支路:
(18)
(19)
(20)
式中:YRL(s)為RL支路對應關系;R、L、C、G為等效支路元件參數,下同。
復數極點對應RLCG支路:
(21)
(22)
R=[-2p′k+2(c′kp′k+c″kp″k)L]L
(23)
(24)
G=-2(c′kp′k+c″kp″k)CL
(25)
本文的方法包括CDSI辨識態空間表達式和構建RLCG等值電路兩部分。
CDSI部分的流程具體如下:
1)初始化CDSI,構造確定-隨機狀態方程和噪聲協方差矩陣。
2)利用劃定為外部系統的電壓、電流采樣數據構建Hankle矩陣。
3)將卡爾曼濾波狀態序列和輸入輸出代入狀態空間方程,求得系統狀態矩陣A、B、C、D。
RLCG部分的流程具體如下:
1)將狀態空間表達式轉化為頻域形式,求解系統零極點。
2)將實數極點和復數極點分別對應其等效支路,構建外部系統等效電路圖。
本文CDSI等值方法流程圖如圖3所示。

圖3 CDSI等值方法流程圖Fig.3 CDSI equivalent method flowchart
以IEEE 14節點系統電網結構為例進行仿真分析,該拓撲圖包含:14條母線,2臺發電機,3臺同步電動機,11個負荷點,20條輸電線路和40個斷路器[28]。在IEEE 14標準測試電力系統中節點14處,接入永磁直驅風電機組,組成分布式風電場,參數如表1所示[29]。仿真中,在風電場內部設置故障,將風電場視為內部系統,IEEE 14標準測試電力系統視為外部系統。具體接線如圖4所示。

表1 永磁直驅風電機組參數Table 1 PMSG parameters

圖4 接入風電機組的IEEE 14節點系統Fig.4 IEEE 14-node bus system connected to wind power
當風電場穩定運行后,在所連母線15設置持續0.05 s的三相瞬時性短路接地故障,引起分布式風電場功率振蕩。
取節點14切除故障后、恢復穩定運行前三相電壓、電流暫態數據。
構建CDSI狀態方程,通過確定-隨機子空間辨識,SVD分解得模型階次指標Moi,i,如圖5所示。

圖5 模型階次指標MOI曲線Fig.5 The curve of the model order indicator
由圖5確定,模型階數n為6,系統階數N為7。使用卡爾曼濾波狀態序列和輸入輸出代入狀態空間方程式(13)求得系統狀態矩陣A、B、C、D,其中系統矩陣A為:
(26)
根據CDSI辨識結果按照式(14)構建狀態空間表達式,并轉化為頻域形式式(15),求解零級點并裂項式(16)—(17),將1組實數極點項和3組復數極點項按照式(18)—(25)構建RLCG等值電路,結果如圖6所示,實現了外部系統暫態等值構建。

圖6 電力系統外部系統暫態等值結構Fig.6 Transient equivalent structure of the external system of the power system
其中以A-A相為例,其等效電路圖如圖7所示,相關參數見表2。

表2 CDSI算法A-A相RLCG等值電路參數Table 2 A-A phase RLCG equivalent circuit parameters of CDSI

圖7 CDSI算法A-A相RLCG等值電路Fig.7 A-A phase RLCG equivalent circuit of CDSI
1)采用SSI方法[18],假定輸入為白噪聲,輸出為上例中采集的電壓暫態數據,構建Hankle矩陣并進行奇異值分解,求得系統矩陣A和C之后,獲得傳遞函數式(12)。
2)采用VF方法[12-16],將上例中采集的電壓、電流暫態數據組成一個超定線性方程組:
(27)
式中:h和d表示一次項參數和常數。
通過求解線性方程組的最小二乘解可以獲得新的極點,經過6次迭代后結果收斂,獲得誤差范圍內的近似解ai、ci、d和h,獲得傳遞函數式(13)。
將SSI、VF算法所得到的傳遞函數,分別采用本文RLCG等值方法,對外部系統進行等值計算而獲得了相應的3輸入3輸出的等值網絡。其中SSI的A-A相等值電路及相關參數如圖8和表3所示,VF方法的如圖9和表4所示。

表3 SSI算法A-A相RLCG等值電路參數Table 3 A-A phase RLCG equivalent circuit parameters of SSI

表4 VF算法A-A相RLCG等值電路參數Table 4 A-A phase RLCG equivalent circuit parameters of VF

圖8 SSI算法A-A相RLCG等值電路Fig.8 A-A phase RLCG equivalent circuit of SSI

圖9 VF算法A-A相RLCG等值電路Fig.9 A-A phase RLCG equivalent circuit of VF
用等效后的外部系統(見圖6)替換原標準IEEE 14節點系統電網,分別對比CDSI與SSI、VF算法等效前后風力發電機故障后0.1 s電磁轉矩,結果見圖10。

圖10 電磁轉矩對比圖Fig.10 Electromagnetic torque comparison
分別計算電磁轉矩的均方根誤差(root mean square error,RMSE)、平均絕對誤差(mean absolute error,MAE)、平均相對誤差(mean relative error,MRE)以及擬合優度(R2),其結果見圖10和表5。

表5 電磁轉矩誤差Table 5 Electromagnetic torque error
結合圖10和表5可以得到,CDSI等值方所得到的風電機電磁轉矩曲線與未等值時的電磁轉矩曲線幾乎完全重合,誤差參數最優,擬合度最高,可以獲得準確的外部系統的暫態等值建模,優于VF和SSI方法。
VF方法計算量大,在求解暫態信號中出現了不穩定節點。而SSI方法僅利用輸出信號,對主導模態提取時出現了虛假模態,使得最終結果誤差較大。
2015年某月,我國某風電場發生脫網事故,導致300 km外火電機組的轉軸扭振保護(torsional tress relay,TSR)動作[30],引起火電機組停運。
根據錄波數據和事故報告,由于風電場內部發生次同步振蕩并在多級電網中傳播引起電網功率振蕩,火電機組扭振保護裝置動作,造成機組跳閘,最終使電網頻率由50.05 Hz降低至49.91 Hz。為簡化系統,方便分析風電場內部的次同步振蕩形成機理,將故障風電場劃分為內部系統。對其余的外部系統進行等值簡化,其電網結構如圖11所示。

圖11 電網拓撲結構圖Fig.11 Grid topology diagram
已獲得了故障后該地區風電場110 kV母線三相電壓、電流60 s內的錄波結果。限于篇幅,文中僅給出A相的電壓和電流的錄波結果,如圖12、13所示,采樣頻率為12.8 kHz。取波動明顯的25~31 s數據進行分析。

圖12 電壓錄波結果Fig.12 Voltage recording results

圖13 電流錄波結果Fig.13 Current recording results
取25~31 s的電壓和電流數據,利用本文的CDSI辨識等值方法,對量測數據進行參數辨識,獲得狀態空間表達式及傳遞函數,按零極點構建外部系統的RLCG等值電路,得到外部系統等值模型。
重構電流輸出信號,其電流輸出結果和局部放大如圖14、15所示,分別對比CDSI、SSI、VF算法輸出結果如圖16所示,輸出電流的誤差結果如表6所示。

表6 電流輸出的誤差值Table 6 The error value of the current output

圖14 電流輸出結果Fig.14 The current output results

圖15 電流輸出結果局部放大Fig.15 The current output results are locally amplified

圖16 電流輸出結果對比Fig.16 Comparison of current output results
對比本文CDSI外部系統等值結果與未等值前采樣獲得的電流數據,其邊界處三相電流曲線誤差小于±1%,其擬合優度(R2)高于99%。在計算復雜性上,量測數據共計76 800×6個,SSI與本文CDSI方法計算速度相同,平均計算時間80~90 s。VF方法需要先進行等間隔采樣,然后迭代求解超定線性方程組,計算量較大,速度較慢,平均計算時間超過1 000 s。
利用本文方法獲得的外部系統電磁暫態等值模型具有很高的等值精度,解決了新能源機組型號多樣、商業保密等原因造成的建模困難,同時也為后續故障還原、事故分析打下堅實基礎。
本文針對以新能源為主體的新型電力系統,提出了一種基于CDSI的外部系統暫態等值方法。
1)利用本文方法將除故障外的其他部分劃分為外部系統進行等效,可以縮短仿真耗時,加快仿真進程,高效準確地進行故障還原和事故分析。
2)CDSI利用輸入和輸出的全部信息獲得系統狀態空間模型,能夠準確反映外部系統的頻域特性,模型具有很高的等值精度,高于僅利用輸出信息的SSI方法。
將電力系統進行外部系統等效后,下一步重點是解決內部系統模型等值,為故障還原、事故分析和擾動源定位等問題提供依據。