馬立元，李永軍，王天輝，李世龍

（軍械工程學(xué)院 導(dǎo)彈工程系，石家莊050003）

作為子空間模態(tài)參數(shù)辨識算法的一個分支，頻域子空間辨識算法出現(xiàn)較晚，最早由T.Mckelvey等人于1994年提出。對于復(fù)雜結(jié)構(gòu)，其輸出點較多，關(guān)心的頻域范圍較廣，因此導(dǎo)致非常大的計算量。很多情況下，實驗室的設(shè)備很難滿足這種大計算量的要求。為解決這一矛盾，筆者提出頻域子空間辨識算法的改進方法，不需要采用濾波手段，能在不減少數(shù)據(jù)量的前提下，減小計算量，得到關(guān)注頻域內(nèi)的模態(tài)參數(shù)。

1 頻域子空間模態(tài)參數(shù)辨識基本原理［1－2］

多自由度振動系統(tǒng)離散時間系統(tǒng)的狀態(tài)方程描述為［3］

式中：ξk＝ej2πk／N為傅里葉算子；U、Y 分別為輸入u、輸出y的傅里葉變換。

則式（2）寫成

由于Y、X和U均為復(fù)數(shù)，則上式可以寫成

Re（·）表示復(fù)數(shù)的實部，Im（·）表示復(fù)數(shù)的虛部。

根據(jù)輸出獲得系統(tǒng)狀態(tài)方程中的矩陣A、C的估計矩陣＾A和＾C，具體算法如下：

構(gòu)建矩陣

式中：A／B表示矩陣A 到矩陣B 的行空間上的投影。對RT22進行奇異值分解

可得到Or的估計值

式中：n為模型的階次，系統(tǒng)的階次常根據(jù)Σ中較大的奇異值個數(shù)確定。

根據(jù)最小二乘法估計矩陣A，C

式中：上標(biāo)“?”表示 Moore-Penrose逆；N0為輸出點個數(shù)。

采用離散系統(tǒng)建模，模態(tài)頻率、阻尼和振型的計算式如下［4］：

式中：λk為矩陣A的特征值；Ψ為矩陣A的特征向量。

2 分段頻域子空間辨識算法

上述的頻域子空間算法，頻域子空間辨識方法中，矩陣Y和U的塊行數(shù)q需要仔細(xì)考慮。為了能進行辨識計算，必須滿足q≥n，n為系統(tǒng)階次。經(jīng)數(shù)值計算經(jīng)驗證明，對于小阻尼系統(tǒng)，q值取的較大為好。一般需要滿足q≥10n［7］。而為了獲得較為準(zhǔn)確的傅里葉變換，又必須有足夠大的數(shù)據(jù)長度。一般情況下，實驗室設(shè)備無法滿足要求。以某一梁結(jié)構(gòu)為例，其輸入點1個，輸出點為5，在一臺式計算機上進行計算，計算機參數(shù)為主頻2.6GHz，內(nèi)存4GB。當(dāng)頻域數(shù)據(jù)點超過2 500時，則產(chǎn)生了內(nèi)存溢出問題。

為此，筆者提出了一種分段計算方法，以減小對計算設(shè)備要求。其方法如下：U（k）∈CN×l，Y（k）∈CN×L.分別表示時域輸入輸出數(shù)據(jù)u，y的離散傅立葉變換，N為采用點數(shù)，l為輸入點數(shù)目，L為輸出點數(shù)目。

將U，Y分成K 段，第i段為

Ni為該段內(nèi)的頻率點數(shù)。

當(dāng)采用離散模型建模時，

表示z域內(nèi)單位圓上的某一點。對于第i段頻域數(shù)據(jù)Ui和Yi，其對應(yīng)的廣義算子為

k0為該段的起始點在全部頻率點中的排序序號。

其值為

對上式按照頻域子空間算法的一般算法進行模態(tài)參數(shù)辨識，則可獲得在第i個頻率段內(nèi)的模態(tài)參數(shù)。將所有頻率段內(nèi)的模態(tài)參數(shù)聚集在一起，則構(gòu)成系統(tǒng)的所有模態(tài)參數(shù)。

分段算法的理論證明如下：

式（5）可以寫為

式中（·）i表示第i段數(shù)據(jù)分解結(jié)果。則

式（7）變?yōu)?/p>

可得到Or的估計值

式中，ni為該頻域段內(nèi)系統(tǒng)階次。

而對系統(tǒng)矩陣的估計為

因此，對頻域分段辨識得到的特征值為系統(tǒng)矩陣特征值的一部分，辨識所得模態(tài)參數(shù)為系統(tǒng)在該頻域段內(nèi)的模態(tài)參數(shù)。但是，需要特別指出的是，在分段辨識模態(tài)參數(shù)時，因其只分析該段內(nèi)的數(shù)據(jù)，系統(tǒng)的階次將發(fā)生變化，需要重新確定系統(tǒng)階次。

3 分段計算辨識矩陣與系統(tǒng)系數(shù)矩陣之間的關(guān)系

因分段計算得到的系統(tǒng)矩陣僅為某一頻段辨識結(jié)果，因此，最終系統(tǒng)矩陣需要進行綜合計算。設(shè)輸入輸出的傅里葉變換數(shù)據(jù)分成K段，第i段數(shù)據(jù)得到的系統(tǒng)矩陣為Ai、Bi、Ci、Di，則

式中（·）（i）表示第i段傅里葉變換數(shù)據(jù)。

由上式可得最終辨識的系統(tǒng)矩陣為對D進行平均的原因是各個頻域數(shù)據(jù)段辨識結(jié)果中的矩陣Di不受頻域段的影響。可以認(rèn)為D≈D1≈…≈Di…≈DK。

4 變截面軸上的應(yīng)用

某裝備的傳動部件受環(huán)境、自身轉(zhuǎn)速等因素影響，經(jīng)常出現(xiàn)故障。為了減少故障率，要使傳動軸的轉(zhuǎn)速避開其固有頻率，因此有必要對傳動軸的模態(tài)參數(shù)進行辨識。

實驗對象是一個橫截面為方形的臺階軸結(jié)構(gòu)，材料是45號鋼，傳感器及激勵點如圖1所示。由于該軸的對稱性，可以只考慮其彎曲模態(tài)，減少傳感器數(shù)量需求。

圖1 變截面軸測點布置

在進行模態(tài)分析實驗時，對實驗對象采用力錘進行激振。采樣頻率為10 000Hz，每個響應(yīng)信號取10 000個采樣點。實驗時利用動態(tài)信號測試儀測得響應(yīng)數(shù)據(jù)，利用Matlab軟件編程實現(xiàn)模態(tài)參數(shù)識別。系統(tǒng)階次根據(jù)SVD分解結(jié)果，確定為4階。

對采集到的數(shù)據(jù)進行傅里葉變換后，頻域點分為3段，分別為：1—2 000點；2 001—4 000點；4 001—5 000點。有限元計算（FEM）能夠使分析者對該框架結(jié)構(gòu)的特征有初步了解，也可以為實驗室的測試工作和參數(shù)辨識計算提供參考，使分析者在實驗之前就能夠大致掌握振型節(jié)點的所在之處，結(jié)構(gòu)的基本振型等情況。對于小型金屬結(jié)構(gòu)，已經(jīng)驗證有限元計算結(jié)果的可靠性［5－6］。表1給出了分段計算辨識結(jié)果與有限元計算結(jié)果的對比。

表1 辨識結(jié)果與有限元結(jié)果比較

圖2—圖5列出了ANSYS計算的振型和分段頻域子空間方法識別的振型對比。

通過表1和圖2—圖5對比發(fā)現(xiàn)，試驗?zāi)B(tài)參數(shù)識別得到的各階模態(tài)參數(shù)與有限元分析得到的分析模態(tài)參數(shù)基本一致，證明了本研究所采用的模態(tài)參數(shù)識別方法基本正確，所得的模態(tài)參數(shù)準(zhǔn)確反映了該部件的動態(tài)特性。

圖2 第1階振型對比

圖3 第2階振型對比

5 結(jié)論

從頻域子空間基本算法出發(fā)，提出了一種采用采用分段頻域子空間的模態(tài)參數(shù)識別方法，其有效性及準(zhǔn)確性被試驗分析所證明。

圖4 第3階振型對比

圖5 第4階振型對比

經(jīng)試驗證明，在頻域數(shù)據(jù)點數(shù)較小的情況下，該算法仍然能夠獲得相對準(zhǔn)確的模態(tài)參數(shù)，從而能夠減小計算量和對計算資源的需求。

［1］ 唐煒.飛機顫振模態(tài)參數(shù)識別方法研究［D］.西安：西北工業(yè)大學(xué)，2006.

［2］ Mckelvey T.Subspace Methods for Frequency Domain Data［A］.In Proceeding of the 2004American Control Conference［C］∥Boston，Massachusetts，USA，2004.

［3］ 常軍.隨機子空間方法在橋梁模態(tài)參數(shù)識別中的應(yīng)用研究［D］.上海：同濟大學(xué)，2006.

［4］ 樊江玲.基于輸出響應(yīng)的模態(tài)參數(shù)辨識［D］.上海：上海交通大學(xué)，2007.

［5］ 薛瑋飛，張智，陳進等.空調(diào)配管模態(tài)仿真分析與實驗研究［J］.機械科學(xué)與技術(shù)，2010，29（6）：822-825.

［6］ 孟杰，陳小安，陳鋒.高速電主軸的試驗?zāi)B(tài)分析［J］.機械設(shè)計，2009，26（6）：70-78.

［7］ Schoukens J，Pintelon R.Estimation of nonparametric noise models［C］∥IEEE Instrumentation and Measurement Technology Conference.Sorrento，Italy，2006.