趙洪山,薛 寧,時(shí) 寧
(華北電力大學(xué) 新能源電力系統(tǒng)國(guó)家重點(diǎn)實(shí)驗(yàn)室,河北 保定 071003)
目前,在電力系統(tǒng)計(jì)算、仿真、分析與控制過程中,經(jīng)常采用非線性高階動(dòng)態(tài)模型來(lái)提高描述系統(tǒng)動(dòng)態(tài)行為的能力,但同時(shí)帶來(lái)了系統(tǒng)模型動(dòng)態(tài)方程維數(shù)高、模型計(jì)算的復(fù)雜度(包括計(jì)算時(shí)間、速度以及存儲(chǔ)空間等)增加等問題,給工程研究人員帶來(lái)了巨大的挑戰(zhàn)。尤其對(duì)于多機(jī)電力系統(tǒng),隨著系統(tǒng)規(guī)模不斷增大,其動(dòng)態(tài)模型的階數(shù)也隨之增加,這使得通過計(jì)算系統(tǒng)模型來(lái)快速完成各種功能變得非常困難[1-3]。
模型降階是解決上述非線性電力系統(tǒng)模型高維動(dòng)態(tài)行為計(jì)算復(fù)雜、難以分析等問題的有效方法之一。學(xué)者們也已進(jìn)行了大量研究,如奇異值攝動(dòng)法[4-7]、解耦算法[8-10]、同調(diào)等值理論[11-14]、模態(tài)分析方法[15-16]、基于相關(guān)性的模型降階方法[17]以及基于Krylov子空間的模型降階方法[18]。但上述系統(tǒng)模型降階多用于分析研究,較少關(guān)注系統(tǒng)的輸入控制,或?qū)⑤斎肟刂谱骱愣ㄌ幚?;而?yīng)用于控制設(shè)計(jì)的系統(tǒng)模型降階,則需要考慮降階模型是否能夠保持原系統(tǒng)的動(dòng)態(tài)行為和輸入/輸出特性的問題。
在眾多模型降階方法中,Gramian方法已在線性電力系統(tǒng)廣泛應(yīng)用。文獻(xiàn)[19]為了簡(jiǎn)化區(qū)域電網(wǎng)振蕩阻尼控制器的設(shè)計(jì),利用Gramian方法研究了電力系統(tǒng)模型降階,并驗(yàn)證基于降階系統(tǒng)設(shè)計(jì)的控制器沒有降低原系統(tǒng)閉環(huán)特性的結(jié)論;文獻(xiàn)[20]也在抑制低頻振蕩問題上利用Gramian降階方法研究了更大規(guī)模電力系統(tǒng)的模型降階問題。但對(duì)于一般非線性動(dòng)態(tài)電力系統(tǒng)模型降階,由于缺少有效的理論方法,其模型降階仍是研究的難點(diǎn)。本文利用經(jīng)驗(yàn)Gramian方法和平衡實(shí)現(xiàn)理論,探討非線性電力系統(tǒng)動(dòng)態(tài)模型的降階過程,并以某實(shí)際多機(jī)非線性電力系統(tǒng)為例,驗(yàn)證經(jīng)驗(yàn)Gramian降階方法對(duì)非線性電力系統(tǒng)模型降階的有效性。
具有輸入和輸出的非線性系統(tǒng)模型一般可以表示為:
其中,f(x,u)描述非線性系統(tǒng)的動(dòng)態(tài)行為;x(t)?Rn為非線性系統(tǒng)的狀態(tài)變量;u(t)?Rp為非線性系統(tǒng)的控制輸入;y(t)?Rq為非線性系統(tǒng)的輸出。
模型降階的基本思想就是對(duì)形如式(1)的非線性動(dòng)態(tài)系統(tǒng),通過已知標(biāo)準(zhǔn)列正交轉(zhuǎn)換矩陣,將原系統(tǒng)變換到一個(gè)低維空間Rk(k≤n)上,并可用低維空間中的降階系統(tǒng)來(lái)表示原系統(tǒng),從而降低了原系統(tǒng)理論分析難度和數(shù)據(jù)運(yùn)算量。假設(shè)已知轉(zhuǎn)換矩陣P,由原始非線性系統(tǒng)(1)可以得到降階系統(tǒng):
目前,模型降階中最重要和最基本的方法就是Krylov子空間類方法,經(jīng)驗(yàn)Gramian方法也基于Krylov子空間原理,即通過構(gòu)建轉(zhuǎn)換矩陣P來(lái)進(jìn)行模型降階。其中,矩陣P是由Galerkin投影獲得的。
Galerkin投影:一旦得到子空間的特征函數(shù),就可以利用Galerkin投影將系統(tǒng)(1)投影到該子空間。Galerkin投影矩陣P通??杀硎緸椋?/p>
其中,Ik×k為 k×k 階的正交矩陣;0k×(n-k)為 k×(n-k)階的全零矩陣。 因此,非線性降階系統(tǒng)(2)中,有 class="picture_character" src="images/2bdd96436f03e63f1ee327b703d8a6f2.png" />(t)=Px(t)。
對(duì)于一般非線性系統(tǒng)(1),其輸入與輸出特性通常由系統(tǒng)的可控Gramian矩陣Wc和可觀Gramian矩陣Wo來(lái)描述,可控Gramian矩陣描述控制輸入變化對(duì)系統(tǒng)狀態(tài)的影響,而可觀Gramian矩陣描述輸出對(duì)系統(tǒng)狀態(tài)變化的反應(yīng)。
對(duì)于線性系統(tǒng),其模型降階可通過求解Lyapunov方程獲得可控Gramian矩陣和可觀Gramian矩陣,然后通過分析Gramian矩陣的特征值,來(lái)確定線性系統(tǒng)模型的降階階數(shù)[21]。而對(duì)于非線性系統(tǒng),由于無(wú)法獲得Lyapunov方程,目前仍沒有直接計(jì)算可控Gramian矩陣和可觀Gramian矩陣的解析方法。Sirovich在文獻(xiàn)[22]中提出:通過仿真或?qū)嶒?yàn)獲得系統(tǒng)狀態(tài)變量和輸出變量樣本,再根據(jù)采樣的經(jīng)驗(yàn)數(shù)據(jù)樣本,計(jì)算得到非線性系統(tǒng)模型的經(jīng)驗(yàn)可控Gramian矩陣和經(jīng)驗(yàn)可觀Gramian矩陣。
經(jīng)驗(yàn)Gramian方法得到的系統(tǒng)數(shù)據(jù)樣本應(yīng)滿足在其穩(wěn)定域內(nèi)這一條件,并且要盡可能反映各種可能的擾動(dòng)行為。為此,定義如下集合:
其中,Tn為一組r個(gè)正交的n×n階矩陣,r為激勵(lì)或者擾動(dòng)方向的矩陣數(shù)量;M為一組s個(gè)正常數(shù),s為在每個(gè)方向不同擾動(dòng)大小的數(shù)量;Ep為Rp中標(biāo)準(zhǔn)單位矢量。
定義1設(shè)Tn、Ep和M為以上所給定集合,對(duì)于非線性系統(tǒng)(1),經(jīng)驗(yàn)可控Gramian矩陣Wc定義如下:
對(duì)于具有輸入和輸出的可控非線性系統(tǒng),只研究輸入與狀態(tài)的行為是不全面的,同時(shí)也應(yīng)考慮輸出反映系統(tǒng)狀態(tài)的行為。下面定義給出了經(jīng)驗(yàn)可觀Gramian矩陣的構(gòu)建方式。
定義2設(shè)Tn、Ep和M為以上所給定集合,對(duì)于非線性系統(tǒng)(1),經(jīng)驗(yàn)可觀Gramian矩陣Wo定義如下:
對(duì)經(jīng)驗(yàn)可控和可觀Gramian矩陣分別進(jìn)行奇異值分解,其非零特征值所對(duì)應(yīng)的特征向量就可生成一個(gè)新子空間,然后可以將系統(tǒng)(1)投影到該子空間,得到降階系統(tǒng)(2)。由此,實(shí)現(xiàn)了反映系統(tǒng)可觀性與可控性非線性系統(tǒng)的模型降階。
由上文可知,一旦得到非線性系統(tǒng)(1)的可控與可觀Gramian矩陣,便可對(duì)其進(jìn)行奇異值分解,并利用可控或可觀Gramian矩陣的奇異值大小來(lái)確定Galerkin投影矩陣P,進(jìn)而得到降階模型的階數(shù)。但對(duì)于一般非線性系統(tǒng),直接得到的可控或可觀Gramian矩陣往往是不相同的,如果以可控Gramian矩陣的奇異值進(jìn)行降階,降階系統(tǒng)的可控性獲得較好的性能,但系統(tǒng)的可觀性稍差;反之,以可觀Gramian矩陣的奇異值進(jìn)行降階,降階系統(tǒng)可得到較好的可觀性,可控性能稍差。
為了使降階系統(tǒng)的可控性和可觀性都能獲得較好的性能,利用平衡實(shí)現(xiàn)理論對(duì)非線性系統(tǒng)(1)進(jìn)行處理,將經(jīng)驗(yàn)可控與可觀Gramian矩陣變換成平衡的可控和可觀Gramian矩陣,即將系統(tǒng)(1)變換成平衡系統(tǒng)。所謂平衡系統(tǒng),是指系統(tǒng)的可控Gramian矩陣與可觀Gramian矩陣是相等的,并且是對(duì)角的。下面給出了利用經(jīng)驗(yàn)可控Gramian矩陣Wc和可觀Gramian矩陣Wo計(jì)算系統(tǒng)(1)的平衡系統(tǒng)的步驟。
a.對(duì)經(jīng)驗(yàn)Wc和Wo,進(jìn)行Cholesky因式分解,得到矩陣X和Y。
b.構(gòu)建積矩陣YTX,并對(duì)積矩陣進(jìn)行奇異值分解,得到對(duì)角矩陣∑、正交矩陣U和V。
其中,∑ ?Rn×n為 Hankel奇異值矩陣。
c.根據(jù)式(8)可獲得平衡轉(zhuǎn)換矩陣T。
至此,利用該平衡變換矩陣T就可獲得系統(tǒng)(1)的平衡系統(tǒng):
平衡系統(tǒng)與其原系統(tǒng)的動(dòng)態(tài)行為是完全相同的,只是平衡系統(tǒng)的狀態(tài)變量變?yōu)橐唤M新的狀態(tài)變量。從線性系統(tǒng)理論而言,平衡系統(tǒng)任一狀態(tài)都是原系統(tǒng)所有狀態(tài)的線性組合,即。 因此,對(duì)非線性系統(tǒng)(1)進(jìn)行模型降階,若可控和可觀Gramian矩陣不等,就可先把系統(tǒng)轉(zhuǎn)換成平衡系統(tǒng),再對(duì)平衡系統(tǒng)進(jìn)行模型降階,該過程稱為平衡實(shí)現(xiàn)。
針對(duì)非線性系統(tǒng)對(duì)應(yīng)的平衡系統(tǒng),計(jì)算其經(jīng)驗(yàn)Gramian矩陣,并對(duì)經(jīng)驗(yàn) Gramian矩陣和進(jìn)行奇異值分解。根據(jù)Hankel矩陣∑奇異值σ的大小,當(dāng)滿足σk?σk+1或者當(dāng)前k個(gè)奇異值能量所占所有奇異值能量絕大部分時(shí),就認(rèn)為k后面的奇異值所對(duì)應(yīng)的狀態(tài)變量對(duì)系統(tǒng)的輸入影響和輸出影響較小,將其略去以實(shí)現(xiàn)模型降階的目的。因此,可以將分為主要和非主要兩部分,如下:
假設(shè)多機(jī)電力系統(tǒng)有Z臺(tái)發(fā)電機(jī),同步發(fā)電機(jī)采用計(jì)及轉(zhuǎn)子超瞬變過程且轉(zhuǎn)子q軸要考慮阻尼繞組的6階動(dòng)態(tài)模型。第i臺(tái)發(fā)電機(jī)模型描述如下:
其中,i=1,2,…,Z;δi為第 i臺(tái)發(fā)電機(jī)轉(zhuǎn)子角;ωi為第i臺(tái)發(fā)電機(jī)角頻率;E′di、E′qi、E″di、E″qi分別為第i臺(tái)發(fā)電機(jī)d軸和q軸暫態(tài)電勢(shì)和次暫態(tài)電勢(shì);ωs為發(fā)電機(jī)額定角頻率;ωref為發(fā)電機(jī)參考角頻率;Pmi為第i臺(tái)發(fā)電機(jī)的原動(dòng)機(jī)輸出機(jī)械功率;Pei為第i臺(tái)發(fā)電機(jī)的電磁功率;Hi為第i臺(tái)發(fā)電機(jī)機(jī)組慣性時(shí)間常數(shù);T′d0i、T′q0i、T″d0i、T″q0i分別為第i臺(tái)發(fā)電機(jī)d軸和q軸開路暫態(tài)、次暫態(tài)時(shí)間常數(shù);Efi為第i臺(tái)發(fā)電機(jī)機(jī)組的勵(lì)磁電壓;Di為第i臺(tái)發(fā)電機(jī)定常阻尼系數(shù);rai、xdi、xqi、x′di、x′qi、x″di、x″qi分別為第i臺(tái)發(fā)電機(jī)定子電阻和d 軸和 q 軸暫同步、瞬變、超瞬變電抗;Idi、Iqi、Udi、Uqi分別為第i臺(tái)發(fā)電機(jī)d軸和q軸暫態(tài)電流和機(jī)端電壓,其關(guān)系為Udi=E″di-raiIdi-x″qiIqi,Uqi=E″qi+x″diIdi-raiIqi。
發(fā)電機(jī)與網(wǎng)絡(luò)的接口方程為:
其中,Ixy=CSIdq;Uxy=CUdq;C 為 dq-xy坐標(biāo)系的變化矩陣;S為網(wǎng)絡(luò)容量和發(fā)電機(jī)容量之間的變化矩陣;Ys為電力系統(tǒng)節(jié)點(diǎn)導(dǎo)納矩陣。
將xy坐標(biāo)下的Ixy和Uxy轉(zhuǎn)換到dq坐標(biāo)系下,代入式(12),消去Idi和Iqi,得到如下形式的具有輸入和輸出多機(jī)非線性電力系統(tǒng)模型:
其中,u(t)=[Pm1… PmZEf1… EfZ]T?Rp;x(t)=[δ1…δZω1…ωZE′q1…E′qZE′d1…E′dZE″q1…E″qZE″d1…E″dZ]T?Rn;f(x,u)描述電力系統(tǒng)動(dòng)態(tài)行為;y(x)為系統(tǒng)的輸出函數(shù),由關(guān)于狀態(tài)x(t)的非線性函數(shù)來(lái)表示。
其中,Ui為第 i臺(tái)發(fā)電機(jī)機(jī)端電壓;y(t)?Rq。
對(duì)于非線性電力系統(tǒng)(13),其經(jīng)驗(yàn)可控Gramian矩陣和可觀Gramian矩陣的計(jì)算,由于涉及積分計(jì)算,計(jì)算量非常龐大。假設(shè)系統(tǒng)的樣本是在離散時(shí)間t1、…、tm得到的,那么狀態(tài)或輸出的相關(guān)矩陣可以寫成如下離散形式:
其中,xk(t)為系統(tǒng)的樣本值;xss為系統(tǒng)狀態(tài)的穩(wěn)態(tài)值。
下面分別給出了經(jīng)驗(yàn)可控和可觀Gramian矩陣的離散定義。
設(shè)Tn、Ep和M為前文所給定的集合,p為系統(tǒng)的控制輸入個(gè)數(shù),則系統(tǒng)經(jīng)驗(yàn)可控Gramian矩陣Wc定義為:
設(shè)Tn、Ep和M為前文所給定的集合,則系統(tǒng)經(jīng)驗(yàn)可觀Gramian矩陣Wo定義為:
利用經(jīng)驗(yàn)樣本數(shù)據(jù),通過式(15)和(16)便可計(jì)算出非線性電力系統(tǒng)(13)的經(jīng)驗(yàn)可控Gramian矩陣Wc和經(jīng)驗(yàn)可觀Gramian矩陣Wo,其計(jì)算精度取決于樣本個(gè)數(shù),樣本數(shù)量越大,精度越高;此外,還受外加激勵(lì)大小的選擇cj、外加激勵(lì)方向的選擇Tlei以及外加激勵(lì)的個(gè)數(shù)p等因素的影響。
對(duì)于非線性電力系統(tǒng)動(dòng)態(tài)模型(13),利用經(jīng)驗(yàn)Gramian進(jìn)行平衡降階方法的步驟如下。
a.根據(jù)定義(15)和(16),計(jì)算非線性多機(jī)電力系統(tǒng)(13)的經(jīng)驗(yàn)可控Gramian矩陣Wc和經(jīng)驗(yàn)可觀Gramian 矩陣 Wo。
b.利用式(6)—(8)計(jì)算平衡變換矩陣 T,得到非線性多機(jī)電力系統(tǒng)(13)的平衡系統(tǒng)模型。
c.利用平衡變換矩陣T,計(jì)算平衡系統(tǒng)的經(jīng)驗(yàn)可控Gramian矩陣和經(jīng)驗(yàn)可觀Gramian矩陣:
e.確定非線性多機(jī)電力系統(tǒng)降階模型的階數(shù)。如果滿足
則平衡系統(tǒng)的降階模型的階數(shù)為r。其中,ε為接近于1的值。利用階數(shù)r構(gòu)造Galerkin投影矩陣P。
f.利用平衡變換矩陣T和Galerkin投影矩陣P,計(jì)算得到形如式(10)的電力系統(tǒng)平衡降階模型。
對(duì)某實(shí)際20機(jī)電力系統(tǒng)模型進(jìn)行降階仿真分析,系統(tǒng)接線圖見圖1。
首先,形成多機(jī)電力系統(tǒng)分析模型,每臺(tái)同步發(fā)電機(jī)采用6階動(dòng)態(tài)模型,負(fù)荷為恒阻抗模型。由于發(fā)電機(jī)模型采用dq坐標(biāo),網(wǎng)絡(luò)模型采用xy坐標(biāo),因此,需要xy-dq坐標(biāo)變換實(shí)現(xiàn)發(fā)電機(jī)與網(wǎng)絡(luò)的互聯(lián),其過程見式(11)—(13),最后形成一個(gè) 120階的形如式(13)的非線性電力系統(tǒng)動(dòng)態(tài)模型。
圖1 某實(shí)際20機(jī)系統(tǒng)接線圖Fig.1 Connection diagram of an 20-generator power system
圖1 某實(shí)際20機(jī)系統(tǒng)接線圖Fig.1 Connection diagram of an 20-generator power system
其次,利用經(jīng)驗(yàn)Gramian平衡降階算法,計(jì)算上述非線性電力系統(tǒng)的平衡變換矩陣T,得到該系統(tǒng)的平衡系統(tǒng)。對(duì)平衡系統(tǒng)的經(jīng)驗(yàn)Gramian矩陣進(jìn)行奇異值分解,利用奇異值分解得到Hankel奇異值矩陣∑,奇異值大小分布如圖2所示。
圖2 Hankel奇異值的分布圖Fig.2 Distribution of Hankel singular values
最后,利用Hankel奇異值大小,根據(jù)模型降階階數(shù)判定條件確定降階系統(tǒng)的維數(shù),得到Galerkin投影矩陣P。由此,利用平衡變換矩陣T和投影矩陣P,得到20機(jī)非線性電力系統(tǒng)的降階模型。
從Hankel奇異值角度分析,當(dāng)ε取0.99,可以將原120階非線性電力系統(tǒng)投影到50維的子空間中;當(dāng)ε取0.95,可將原系統(tǒng)降為32階的非線性動(dòng)態(tài)系統(tǒng)。理論上,ε取值越小,降階系統(tǒng)的階數(shù)就越小,但投影到子空間后的動(dòng)態(tài)行為失真就越大。因此,需要合理選擇降階模型的階數(shù)。
對(duì)于該20機(jī)實(shí)際電力系統(tǒng),每臺(tái)發(fā)電機(jī)均考慮調(diào)速控制和勵(lì)磁控制功能,這些控制器的參數(shù)在各降階系統(tǒng)中都保持不變。在考慮三相短路故障(故障點(diǎn)在線路33-34之間50%處)情況下,對(duì)各降階模型的暫態(tài)動(dòng)態(tài)變化過程進(jìn)行仿真分析。為了節(jié)省篇幅,以發(fā)電機(jī)7為例,給出了其相應(yīng)的仿真曲線如圖3—5所示。圖中,功角δ、角速度ω、輸出電壓Ut均為標(biāo)幺值。
圖3 降階系統(tǒng)與原系統(tǒng)功角仿真曲線Fig.3 Simulative angle curves of reduced system and original system
圖4 降階系統(tǒng)與原系統(tǒng)角速度仿真曲線Fig.4 Simulative angular speed curves of reduced system and original system
圖5 降階系統(tǒng)與原系統(tǒng)輸出電壓仿真曲線Fig.5 Simulative output voltage curves of reduced system and original system
從仿真曲線圖3—5中可以看出:對(duì)于50階的降階系統(tǒng),除功角幅值誤差稍大些外,角速率和電壓幅值與原系統(tǒng)輸出相比誤差都很小,而且3個(gè)輸出響應(yīng)的頻率特性(即暫態(tài)過程)與原系統(tǒng)的幾乎相同;當(dāng)系統(tǒng)降到49階時(shí),功角、角速度和電壓幅值3個(gè)輸出響應(yīng)的幅值變化不大,但較50階的降階系統(tǒng),其與原系統(tǒng)的誤差要比50階系統(tǒng)大,并有逐步增大的趨勢(shì)。尤其,當(dāng)系統(tǒng)降到47階或低于47階時(shí),降階系統(tǒng)輸出響應(yīng)的幅值誤差增大,觀察功角和角速度曲線出現(xiàn)頻率偏移,即產(chǎn)生了頻率失真。從Hankel奇異值角度分析,投影到50維子空間的能量約為原系統(tǒng)能量的99.45%,其能保留原系統(tǒng)的大部分能量,即保留了原系統(tǒng)的各種特性和動(dòng)態(tài)行為。
通過對(duì)所有節(jié)點(diǎn)仿真曲線進(jìn)行統(tǒng)計(jì)分析可知:當(dāng)ε選擇大于0.99時(shí),降階模型的誤差都比較小,能夠?qū)⒃到y(tǒng)的輸入和輸出動(dòng)態(tài)行為保留下來(lái)。為了使降階系統(tǒng)得到較為精確的幅頻特性,通常選擇ε為0.99。本例中,降階系統(tǒng)模型的階數(shù)可降到50階,其所有節(jié)點(diǎn)狀態(tài)和輸出與原系統(tǒng)相比,誤差都非常小。
本文提出了一種針對(duì)多機(jī)非線性電力系統(tǒng)模型的經(jīng)驗(yàn)Gramian平衡降階方法,并通過某實(shí)際20機(jī)非線性多機(jī)電力系統(tǒng)模型進(jìn)行仿真驗(yàn)證。仿真結(jié)果表明:降階系統(tǒng)能夠很好地保留原非線性電力系統(tǒng)輸入和輸出的動(dòng)態(tài)行為以及原非線性電力系統(tǒng)的穩(wěn)態(tài)值。對(duì)于該20機(jī)非線性電力系統(tǒng),可以將120階非線性動(dòng)態(tài)模型降到50階,仿真結(jié)果證實(shí)了所提方法在非線性電力系統(tǒng)模型降階中的有效性。