K-power雙線性系統(tǒng)基于Laguerre函數(shù)的保結(jié)構(gòu)模型降階方法

金輝 ,肖志華 ,祁振中

(1.長(zhǎng)江大學(xué)信息與數(shù)學(xué)學(xué)院,湖北 荊州 434023;2.西北大學(xué)數(shù)學(xué)學(xué)院,陜西 西安 710127)

1.引言

在工程應(yīng)用領(lǐng)域中,許多物理系統(tǒng)或現(xiàn)象都可以用數(shù)學(xué)模型來(lái)描述,然而大多數(shù)的數(shù)學(xué)模型規(guī)模都極其復(fù)雜,使得對(duì)這些系統(tǒng)的直接仿真模擬和理論分析變得十分困難.因此,有效減小系統(tǒng)規(guī)模和縮短計(jì)算仿真時(shí)間的需求日益上升.模型降階正是基于此思想,將一個(gè)較大規(guī)模的復(fù)雜系統(tǒng)轉(zhuǎn)化為一個(gè)近似的較小系統(tǒng),同時(shí)保留了原始系統(tǒng)的一些重要性能,例如穩(wěn)定性、結(jié)構(gòu)性和無(wú)源性等[1].模型降階方法已經(jīng)普遍用于超大規(guī)模集成電路模擬、計(jì)算電磁學(xué)、微機(jī)電系統(tǒng)和流體力學(xué)等領(lǐng)域[2].

模型降階自提出以來(lái),已經(jīng)形成了一系列成熟的降階理論和方法,例如Krylov子空間方法、平衡截?cái)喾椒ā⒈菊髡环纸夥椒ê驼欢囗?xiàng)式方法等.其中平衡截?cái)嗍悄Ｐ徒惦A最著名的方法之一[3],其基本思想是計(jì)算出控制系統(tǒng)的可控Gram矩陣和可觀Gram矩陣,進(jìn)一步構(gòu)造合適的映射子空間獲得降階模型.該方法可以構(gòu)造出保持原始系統(tǒng)穩(wěn)定性的降階模型,并且給出近似先驗(yàn)誤差.然而其需求解兩個(gè)大規(guī)模Lyapunov方程,該過(guò)程計(jì)算量大且計(jì)算時(shí)間長(zhǎng),這導(dǎo)致了其應(yīng)用的局限性.在這種情況下,許多學(xué)者作出相應(yīng)改進(jìn),提出了近似平衡截?cái)喾椒╗4].該方法以一種數(shù)值高效的方式獲得近似平衡系統(tǒng),被廣泛應(yīng)用于線性系統(tǒng)及非線性系統(tǒng)的模型降階中.

對(duì)于線性系統(tǒng),相應(yīng)的模型降階方法已經(jīng)日趨完善.而實(shí)際的物理系統(tǒng)本質(zhì)是非線性的,許多物理系統(tǒng)都可以通過(guò)雙線性系統(tǒng)模型來(lái)描述.雙線性系統(tǒng)是一類特殊的非線性系統(tǒng),是連接線性系統(tǒng)和非線性系統(tǒng)之間的橋梁,該系統(tǒng)經(jīng)常出現(xiàn)在非線性電路、熱過(guò)程和生態(tài)系統(tǒng)等科學(xué)領(lǐng)域[5].關(guān)于雙線性系統(tǒng)的模型降階,已經(jīng)發(fā)展了多種方法,例如平衡相關(guān)方法、Krylov子空間方法和插值法等[6-8].值得注意的是,Shaker和Tahavori在文[9]中定義了雙線性系統(tǒng)的交叉Gram矩陣.該交叉Gram矩陣同時(shí)反映了系統(tǒng)的可控性和可觀性,并且當(dāng)雙線性系統(tǒng)有界輸入有界輸出(BIBO)穩(wěn)定時(shí),其為廣義Sylvester方程的解.作為雙線性系統(tǒng)的一個(gè)特殊子類,K-power雙線性系統(tǒng)的應(yīng)用也很廣泛,涉及液壓驅(qū)動(dòng)建模、多項(xiàng)式系統(tǒng)、非線性系統(tǒng)辨識(shí)等領(lǐng)域[10-11].關(guān)于K-power雙線性系統(tǒng)的模型降階方法,在文[10]中,Al-Baiyat和Bettayeb首先證明了K-power雙線性系統(tǒng)具有塊對(duì)角可控Gram矩陣和可觀Gram矩陣,并提出了一種簡(jiǎn)化求解Lyapunov方程計(jì)算的保結(jié)構(gòu)平衡截?cái)喾椒?XIE和Syrmos[11]證明了K-power雙系統(tǒng)的雙線性廣義奇異攝動(dòng)近似(GSPA)可以轉(zhuǎn)化為σ-交互系統(tǒng)的直接截?cái)?基于此,提出了一種基于雙線性σ-交互系統(tǒng)分析的方法.WANG和JIANG[12]還成功地將矩匹配方法和最優(yōu)H2方法應(yīng)用于K-power雙線性系統(tǒng)的模型降階,進(jìn)一步得到保持原始系統(tǒng)結(jié)構(gòu)特征的降階模型.

本文基于Laguerre函數(shù)和雙線性系統(tǒng)的交叉Gram矩陣,提出了一種K-power雙線性系統(tǒng)的保結(jié)構(gòu)模型降階方法.該方法首先將其等價(jià)雙線性系統(tǒng)的交叉Gram矩陣用Laguerre函數(shù)進(jìn)行展開(kāi),并結(jié)合其正交性,進(jìn)一步給出交叉Gram矩陣的低秩因子,最后通過(guò)構(gòu)造K-power雙線性系統(tǒng)各子系統(tǒng)的投影變換得到降階模型.相比于經(jīng)典的平衡截?cái)喾椒?此方法不需要直接求解Sylvester方程來(lái)獲得交叉Gram矩陣,具有較高的計(jì)算效率和靈活度,且具有一定的自適應(yīng)性.

本文首先介紹Laguerre函數(shù)的相關(guān)重要性質(zhì),接著詳細(xì)給出保結(jié)構(gòu)的K-power雙線性系統(tǒng)降階過(guò)程,并給出了相應(yīng)的降階算法.最后,通過(guò)數(shù)值實(shí)驗(yàn)有效驗(yàn)證了算法的有效性.

2.矩陣指數(shù)函數(shù)的Laguerre函數(shù)展開(kāi)

矩陣指數(shù)函數(shù)eAt又稱狀態(tài)轉(zhuǎn)移矩陣,在現(xiàn)代控制理論中,可以將其用于狀態(tài)方程的求解.本節(jié)的目標(biāo)是利用Laguerre函數(shù)[13]對(duì)矩陣指數(shù)函數(shù)進(jìn)行展開(kāi).

類似地,對(duì)于含有一個(gè)穩(wěn)定矩陣A(若Cn×n的所有特征值嚴(yán)格位于復(fù)平面的左半平面,則稱A為穩(wěn)定矩陣) 的矩陣指數(shù)函數(shù)eAt,可作下列Laguerre展開(kāi)[16]

其中Laguerre系數(shù)矩陣Ak滿足

I為單位矩陣.

3.K-power雙線性系統(tǒng)的保結(jié)構(gòu)模型降階算法

考慮如下K-power雙線性系統(tǒng)[1,17]：

雙線性系統(tǒng)(3.3)的交叉Gram矩陣[9]為

其滿足如下廣義的Sylvester方程:

i1,2,···,l,j1,2,···,N ?1.

其中{α1,α2,···,αN}為CF-ADI參數(shù),且Re(αj)<0.假設(shè)所有CF-ADI參數(shù)為?α,我們有

因此,除了符號(hào)的正負(fù)性,低秩因子F與CF-ADI迭代得到的低秩因子Zn相同,即,當(dāng)所有參數(shù)αj與α相同時(shí),上述方法與ADI方法具有等價(jià)性。

由式(3.5)得到交叉Gram矩陣Rl的低秩因子,可以將其應(yīng)用到等價(jià)的雙線性系統(tǒng)來(lái)構(gòu)造投影變換矩陣,進(jìn)而生成降階模型.但該方法沒(méi)有保留K-power雙線性系統(tǒng)自身的結(jié)構(gòu),而在工程設(shè)計(jì)和控制應(yīng)用中,通常要求保留原系統(tǒng)的結(jié)構(gòu).K-power雙線性系統(tǒng)包含k個(gè)子系統(tǒng),且這些子系統(tǒng)依次耦合,其具有保留該系統(tǒng)耦合結(jié)構(gòu)的必要性.為此,通過(guò)構(gòu)造該系統(tǒng)各子系統(tǒng)對(duì)應(yīng)的投影矩陣的思想[19],首先對(duì)F和G作如下分塊:

可得系統(tǒng)(3.1)的降階系統(tǒng)如下

上述降階過(guò)程可由如下算法1描述：

算法1K-power雙線性系統(tǒng)基于Laguerre函數(shù)的保結(jié)構(gòu)模型降階算法

其中σ(A)為A的特征值.對(duì)于更一般的情況,如何選取最優(yōu)的Laguerre參數(shù)α需要進(jìn)一步研究.數(shù)值結(jié)果表明,該最優(yōu)問(wèn)題對(duì)我們算法中參數(shù)α的選取具有指導(dǎo)作用.

降階系統(tǒng)(3.6)可以寫為如下等價(jià)的雙線性系統(tǒng)形式

定理3.1[20]對(duì)于雙線性系統(tǒng)Σ,如果A是穩(wěn)定的,并且可控Gram矩陣P和可觀Gram矩陣Q均存在,則雙線性系統(tǒng)Σ的H2范數(shù)表示為

對(duì)于精確求解雙線性系統(tǒng)可控Gram矩陣P和可觀Gram矩陣Q的平衡截?cái)喾椒╗6]和某些H2優(yōu)化模型降階方法[17],可以得到誤差系統(tǒng)的Hankel范數(shù)或H2范數(shù).一般而言,對(duì)于近似平衡截?cái)喾椒?即,通過(guò)近似求解或分解系統(tǒng)可控Gram矩陣P和可觀Gram矩陣Q,很難得到與精確平衡截?cái)喾椒愃频恼`差估計(jì)式.盡管如此,值得注意的是,XIANG在文[21]中證明了矩陣指數(shù)函數(shù)正交多項(xiàng)式近似的指數(shù)收斂性,受此啟發(fā),我們可以將該結(jié)論推廣應(yīng)用到(雙)線性系統(tǒng)Gram矩陣的低秩分解的收斂性分析,并由此估計(jì)降階系統(tǒng)與原系統(tǒng)的誤差,該問(wèn)題有待進(jìn)一步的研究.

4.數(shù)值實(shí)驗(yàn)

下面通過(guò)一個(gè)數(shù)值算例來(lái)驗(yàn)證上述保結(jié)構(gòu)模型降階算法的有效性.該數(shù)值實(shí)驗(yàn)在Intel(R)Core (TM) i5-10210U CPU (1.60 GHZ) PC上進(jìn)行,內(nèi)存為16.00GB.同時(shí),選用ode15s函數(shù)在MATLAB R2018b中求解系統(tǒng).

考慮由兩個(gè)子系統(tǒng)組成的單輸入單輸出(SISO)K-power雙線性系統(tǒng)[12,17,19],并選取每個(gè)子系統(tǒng)的階數(shù)均為3000.該系統(tǒng)的系數(shù)矩陣如下:

其中B1,R3000.通過(guò)本文給出的算法1(Algorithm 1)對(duì)該系統(tǒng)進(jìn)行降階,該系統(tǒng)優(yōu)化問(wèn)題(3.7)的解為α ≈10.數(shù)值結(jié)果表明,該α值是有效的.取Laguerre函數(shù)展開(kāi)項(xiàng)數(shù)N4,誤差限tol10-10,由此自適應(yīng)得到的兩個(gè)降階子系統(tǒng)的階數(shù)均為4.同時(shí),分別以文[22]中的Krylov子空間方法(Krylov),文[19]中的Laguerre正交多項(xiàng)式方法(Laguerre)和文[10]中的平衡截?cái)喾椒?BT)對(duì)該系統(tǒng)進(jìn)行降階.在Krylov子空間方法中,取q19匹配該系統(tǒng)第一個(gè)子系統(tǒng)的前9階矩,同時(shí)選取p22,q24匹配第二個(gè)子系統(tǒng)的前4階矩[19],由此得到的兩個(gè)降階子系統(tǒng)的階數(shù)為9和8.由Laguerre多項(xiàng)式方法和BT方法分別得到的兩個(gè)降階子系統(tǒng)的階數(shù)也為9和8.

表1分別給出了該系統(tǒng)的輸入變量u1(t)10sin(10t+5)和u2(t)sin(4t)e-0.2t時(shí)各模型降階方法的CPU運(yùn)行時(shí)間(包括降階的時(shí)間和模擬降階系統(tǒng)的時(shí)間)、最大相對(duì)誤差∥y(t)?yr(t)∥2/∥y(t)∥2和加速比,其中yr(t)為降階系統(tǒng)的輸出響應(yīng).

表1 各模型降階方法的CPU運(yùn)行時(shí)間、最大相對(duì)誤差和加速比

圖1和圖2分別給出了輸入變量為u1(t)時(shí)各降階系統(tǒng)和原始系統(tǒng)的瞬態(tài)響應(yīng)及其對(duì)應(yīng)的相對(duì)誤差.

圖1 輸入為u1時(shí)各降階系統(tǒng)的瞬態(tài)響應(yīng)

圖2 輸入為u1(t)時(shí)瞬態(tài)響應(yīng)的相對(duì)誤差

圖3和圖4分別給出了輸入變量為u2(t)時(shí)各降階系統(tǒng)和原始系統(tǒng)的瞬態(tài)響應(yīng)及其對(duì)應(yīng)的相對(duì)誤差.

圖4 輸入為u2(t)時(shí)瞬態(tài)響應(yīng)的相對(duì)誤差

從模擬結(jié)果可以看出,對(duì)于該SISO K-power雙線性系統(tǒng)模型,由算法1得到的降階系統(tǒng)對(duì)原始系統(tǒng)有很好地瞬態(tài)響應(yīng)近似效果.在相同的近似精度下,算法1構(gòu)造的降階系統(tǒng)比BT方法構(gòu)造的降階系統(tǒng)的階數(shù)小.算法1構(gòu)造的降階系統(tǒng)比Krylov子空間方法和Laguerre多項(xiàng)式方法構(gòu)造的降階系統(tǒng)有略好的近似精度,并且算法1得到的降階系統(tǒng)的階數(shù)小于這兩種方法(Krylov,Laguerre)得到的降階系統(tǒng)的階數(shù).此外,算法1構(gòu)造和模擬降階系統(tǒng)的時(shí)間少于其他三種方法(Krylov,Laguerre,BT)所用的時(shí)間.這些結(jié)果表明,本文提出的算法對(duì)于該K-power雙線性系統(tǒng)是有效的.

5.結(jié)論

本文提出了一種K-power雙線性系統(tǒng)基于Laguerre函數(shù)的保結(jié)構(gòu)模型降階方法.該方法首先將K-power雙線性系統(tǒng)轉(zhuǎn)化為等價(jià)的雙線性系統(tǒng),再結(jié)合Laguerre函數(shù)的正交性,將矩陣指數(shù)函數(shù)通過(guò)Laguerre函數(shù)展開(kāi),然后計(jì)算出雙線性系統(tǒng)交叉Gram矩陣的低秩分解因子,從而構(gòu)造K-power雙線性系統(tǒng)各子系統(tǒng)相應(yīng)的投影變換,最后得到保結(jié)構(gòu)的降階模型.通過(guò)所提算法得到的降階系統(tǒng)能很好地近似原始系統(tǒng)的瞬態(tài)相應(yīng),并且該方法計(jì)算高效,具有一定的自適應(yīng)性.最后,數(shù)值實(shí)驗(yàn)有效驗(yàn)證了該算法的有效性.