孫鑫暉， 郝木明， 王淮維

(中國石油大學(xué)(華東)化學(xué)工程學(xué)院，山東 東營 266555)

實驗?zāi)B(tài)分析經(jīng)歷了幾十年的發(fā)展歷程，已經(jīng)從單自由度發(fā)展為多自由度，由單輸入單輸出發(fā)展為多輸入多輸出，由局部估計發(fā)展為整體估計［1]，并且新的識別方法不斷出現(xiàn)。其中PloyMAX方法［2，3]是近幾年發(fā)展起來的頻域識別算法，該方法也稱作多參考點最小二乘復(fù)頻域算法(poly－reference least square complex frequencydomain)［4]。由于比較好的解決了頻域算法中數(shù)值病態(tài)問題，該方法可以處理很寬的頻率范圍和很高的模態(tài)階次。并且能夠產(chǎn)生非常清晰的穩(wěn)定圖，適用于大阻尼與高噪聲數(shù)據(jù)，目前已經(jīng)成為一種工業(yè)標(biāo)準(zhǔn)。

本文通過分析該算法實現(xiàn)過程中的矩陣結(jié)構(gòu)，給出了縮減正則方程的近似計算方法，從而避免了矩陣求逆以及乘法運算，使得計算量大為降低，提高了模態(tài)參數(shù)識別速度。

1 PolyMAX算法描述

1.1 頻響函數(shù)的右矩陣分式模型

對于一個具有No個輸出，Ni個輸入的多輸入多數(shù)出系統(tǒng)(MIMO)，系統(tǒng)的頻響函數(shù) H(ω)∈CNo×Ni可以通過右矩陣分式模型進行描述:

其中矩陣多項式 B(ω)∈CNo×Ni、A(ω)∈CNi×Ni定義為:

其中 Ωj(ω)=exp(－iωTs·j)為基函數(shù)(Ts為采樣周期)，n為模型階次，Aj、Bj(j=0…n)矩陣多項式系數(shù)。

1.2 正則方程的縮減

通過對式(1)右乘A(ω)進行線性化處理，然后利用最小二乘法求解系數(shù)矩陣。但是直接求解存在以下缺點是:

(1)當(dāng)頻響函數(shù)規(guī)模較大、模型階次較高時導(dǎo)致正則方程的維數(shù)太高，直接求解難以實現(xiàn)。

(2)分母系數(shù)矩陣Aj不能獨立于分子系數(shù)矩陣Bj求解，而模態(tài)參數(shù)的求解只需要分母系數(shù)矩陣Aj。

文獻［4] 給出了一種通過縮減正則方程求解分母系數(shù)矩陣的方法，該方法能夠大幅減小正則方程的維數(shù)且無需求解分子系數(shù)矩陣Bj。這里給出縮減正則方程的表達(dá)式，具體推導(dǎo)參見文獻［4] 。

式(4)中Xo、Yo的表達(dá)式如下:

其中?為 kronecker積，Nf為識別頻段內(nèi)包含的譜線數(shù)。

1.3 通過系數(shù)矩陣得到模態(tài)參數(shù)

當(dāng)系數(shù)矩陣Aj確定之后，模態(tài)參數(shù)可以通過式(6)中友矩陣進行特征值分解得到:

其中極點e－λiTs位于Λ的對角線上，頻率與阻尼可以通過式(7)求得，模態(tài)參數(shù)與向量位于V最后Ni行對應(yīng)的列向量。

2 PolyMAX算法的快速實現(xiàn)

PolyMAX算法的識別可以按照以下步驟進行:

① 確定出最高計算階次nmax、最低計算階次nmin。

② 通過式(4)得到階次n對應(yīng)的Ro、So、To。

③ 計算縮減的正則方程M，然后通過最小二乘求解系數(shù)矩陣Aj。根據(jù)式(6)、式(7)得到模態(tài)參數(shù)。

④ 改變系統(tǒng)階次建立穩(wěn)定圖。

在計算 Ro、So、To的過程中，當(dāng) Xo、Yo維數(shù)很高時，一般不直接通過式(4)中矩陣乘法計算，而是采用FFT 運算代替矩陣乘法［4，5]。當(dāng) Ro、So、To計算完成之后，計算量主要集中在正則方程中矩陣M上，因為在建立穩(wěn)定圖的過程中，每一個階次的改變，都需要重新計算矩陣M，通過式(3)可以看出，這需要大量的矩陣求逆以及矩陣乘法運算。對于復(fù)雜結(jié)構(gòu)，輸出數(shù)目N0很多、計算階次n很高，在這種情況下計算量尤其突出。

本文通過分析正則方程中矩陣R0的結(jié)構(gòu)，給出其逆矩陣的近似計算方法。從而避免了矩陣的求逆以及乘法運算。通過將基函數(shù)表達(dá)式 Ω.j(ωk)=exp·(－iωkTs·j)=exp(－i(k－1)/(Nf－1)·j)帶入到X0，然后由式(4)得到Ro的表達(dá)式如下:

由式(8)可以得出:當(dāng) m=n時，［Ro]mn=Nf

當(dāng)m－n為偶數(shù)時:

當(dāng)m－n為奇數(shù)時:

因此Ro具有如下形式:

由于PolyMax為寬頻帶模態(tài)參數(shù)識別算法，識別范圍內(nèi)的譜線數(shù)Nf?1，因此:

其中In+1為n+1維的單位矩陣。則M可以寫成:

3 仿真算例

圖1 七自由度系統(tǒng)Fig.1 Seven DOF system

圖2 穩(wěn)定圖Fig.2 Stabilization diagram

圖3 渦輪盤測點布置(箭頭為參考點)Fig.3 Measurement setup of gas turbine disk

表1 常規(guī)實現(xiàn)與快速實現(xiàn)的識別結(jié)果比較Tab.1 Identification results between normal and fast implementation

4 實測算例

實測算例來自于某航空發(fā)動機渦輪盤模態(tài)實驗，目的是驗證當(dāng)數(shù)據(jù)規(guī)模比較大時，快速實現(xiàn)方法所能帶來的效率提高。實驗通過錘擊法進行測試，測點布置如圖3所示。其中FRF中包含168個輸出，3個輸入。分析頻率范圍0－8 000 Hz，譜線數(shù)Nf=1 600。在分析范圍內(nèi)大約有50多個結(jié)構(gòu)模態(tài)，根據(jù)式(14)中模型階次n與結(jié)構(gòu)模態(tài)數(shù)目Nm之間的關(guān)系，可以得出所需的模型階次n大約為100階。

圖4為在不同階次下的計算時間，每一次計算中階次范圍都為20(運行環(huán)境為Pentium M 、主頻1.6 GHz、1 G內(nèi)存、Matlab6.5編程)。從圖中可以看出在計算階次不高的情況下，二者的計算時間相差不大。隨著計算階次的提高，快速實現(xiàn)方法所用時間明顯少于常規(guī)實現(xiàn)方法，并且階次越高，差距越明顯。這說明快速實現(xiàn)方法適合于復(fù)雜結(jié)構(gòu)的大規(guī)模運算。圖5為在最高計算階次100下得到的穩(wěn)定圖。在此情況下，常規(guī)實現(xiàn)所用時間為123.6 s，而快速實現(xiàn)方法所用時間為42.7s。圖6為識別的部分振型(前9階)。

圖4 計算時間的比較Fig.4 Comparison of computation time

圖5 穩(wěn)定圖Fig.5 Stabilization diagram

圖6 識別的振型Fig.6 Identified modeshapes

5 結(jié)論

本文給出了PolyMAX算法的一種快速實現(xiàn)方法，該方法通過避免縮減正則方程中的求逆運算與矩陣乘法，減少了運算量。實測算例表明:在數(shù)據(jù)量比較大的情況下，與常規(guī)實現(xiàn)相比較，采用快速實現(xiàn)后參數(shù)識別時間顯著減少，證明了該方法的可行性。并且輸出數(shù)目越多、系統(tǒng)階次越高的情況效果越明顯，因此該方法適合于復(fù)雜結(jié)構(gòu)的模態(tài)參數(shù)識別。

［1] 傅志方，華宏星.模態(tài)分析理論與應(yīng)用［M] .上海:上海交通大學(xué)出版社，2000.

［2] Peeters B，Van Der Auweraer H，Guillaume P，et al.The PolyMAX frequency-domain method:a new standard for modal parameter estimation ［J] .Shock and Vibration，2004，11:395－409.

［3] Peeters B，Guillaume P.Automotive and aerospace applictions of the LMS PolyMAX modal parameter estimation method［C] //Proc of the 22th International Modal Analysis Conference，Dearborn，USA，January，2004.

［4] Guillaume P，Verboven P，Vanlanduit S，et al.A polyreference implementation of the least-squares complex frequency domain estimator［C] // Proc of the 21th International Modal Analysis Conference.Kissimmee，USA，F(xiàn)ebruary，2003.

［5] Schoukens J，Rolain Y，Gustafsson F，et al.Fast calculation of least-squares for system identification［C] //Proc on the 37th IEEE Conference on Decision ＆ Control，F(xiàn)lorida，USA，December，1998，3408 －3410.

［6] Li Z F，Hua H X，Song H W，et al.Extracting modal parameters from structures undergoing ambient excitation excitation［J] .Journal of Shanghai Jiaotong University，2001，2:117－122.

［7] Shih C Y，Tsuei Y G，Allemang R J.Complex mode indication function and its applications to spatial domain parameter estimation［J] .Mechanical System and Signal Processing，1988，2(4):367 －377.