基于響應傳遞率的非白隨機激勵僅輸出結構模態參數辨識

周思達，劉 莉，楊 武，馬志賽

（北京理工大學宇航學院飛行器動力學與控制教育部重點實驗室，北京 100081）

基于響應傳遞率的非白隨機激勵僅輸出結構模態參數辨識

周思達，劉 莉，楊 武，馬志賽

（北京理工大學宇航學院飛行器動力學與控制教育部重點實驗室，北京 100081）

在對工作狀態的航天器等結構進行僅輸出模態參數辨識時，激勵非白造成了無法剔除諧波引起的虛假模態。為解決此問題，該研究對最新的基于結構響應傳遞率的模態參數辨識方法進行了改進，提出了一種改進的基于多參考點響應傳遞率的僅輸出模態參數辨識方法。首先，推導了多參考點的結構響應傳遞率表達式，建立了其左矩陣分式多項式的參數化模型，進而給出了辨識問題的最小二乘估計，利用正則方程Jacobi矩陣的分塊性質對最小二乘問題矩陣形式完成了縮減，降低了計算量；然后，通過高維伴隨矩陣方法解決了矩陣多項式的特征值求解問題，即多參考多輸出的模態參數求解問題，以及通過矩陣偽逆解決了現有方法中載荷工況數與響應點數和參考點數的約束問題；最后，通過兩個數值算例對所提方法進行了驗證，辨識結果表明：辨識方法能夠很好的辨識出結構的模態頻率、阻尼比和模態振型，且能夠很好的避免激勵中含有的諧波分量對辨識結果造成的影響，解決了傳統僅輸出模態參數辨識中激勵非白對辨識結果造成的影響問題。

模態參數辨識；僅輸出；響應傳遞率；非白隨機激勵

隨著航天科技的快速發展，以及相關領域的迫切需要，對大型空間結構的需求變得愈加迫切。由于航天器對結構重量的限制，其結構往往具有大柔性。這樣的大柔性造成的結構振動，會對整器姿態控制、關鍵儀器性能與力學環境、結構健康等方面產生不可忽視的影響。

通常，對于分布參數的結構振動控制，從整器姿態控制的觀測量中濾除柔性振動量，以及結構健康監測等方面，都涉及到結構的模態參數，甚至基于模態參數。辨識得到的模態為以上應用提供了有效的低階模型。獲取模態參數的途徑有很多，然而不同于理論分析和數值計算，以及地面的模態實驗手段，通過測量在軌運行航天器結構動力學響應，從而通過辨識的手段獲取航天器結構模態參數，可以很好的考慮航天器在軌運行工作狀態中受環境因素的真實影響。

然而，對在軌航天器進行模態參數辨識也面臨巨大的挑戰。在軌航天器承受來自環境的隨機激勵，以及自身機構運動和控制元件的驅動力，這些激勵是無法準確測量的，因此在軌航天器結構模態參數辨識是僅輸出的。這樣針對運行結構的僅輸出的模態分析被稱為工作模態分析，已經在工程［1］、建筑工程［2－3］、航空工程［3］、機械工程［4］等領域得到了一定應用?，F有方法多基于白噪聲假設，即假設作用在結構上的激勵為白噪聲，那么結構響應的功率譜函數就只是模態參數（除模態參與因子）和激勵幅值的函數，而后者由于白噪聲假設，其功率函數為關于頻率的常函數，從而可僅由結構響應辨識出結構的模態參數。

由于航天器工作的特殊環境，以及動量輪等姿控器件、液體晃動等產生的諧波激勵，造成了作用在航天器結構上的激勵為非白噪聲，因此現有基于白噪聲假設的僅輸出模態參數辨識方法無法滿足要求。

旋轉機械的工作模態分析也遇到了相似的非白噪聲激勵問題。過去，研究者利用預置陷波濾波器的方法濾除諧波影響，取得了一定的效果［5］。此類方法需要預先知道激勵中諧波激勵的頻率和水平，但是在很多情況下是無法滿足的。最近，有學者提出了基于響應傳遞率的技術處模態參數辨識方法，以克服工作模態分析中激勵非白對辨識結果的影響。Devriendt等［6］第一次提出了基于響應傳遞率模態參數辨識的概念，之后，他們又在多輸出、辨識算法、求解方法方面做了一些改進［7－9］。

為解決激勵非白對結構模態參數辨識的影響問題，本文對現有基于響應傳遞率的模態參數辨識方法進行改進，提出了一種利用矩陣偽逆的頻域僅輸出模態參數辨識方法。

1 多參考的結構響應傳遞率

由于傳遞率定義為結構不同位置力或響應量頻譜的比值，表達的是結構不同位置力或響應量的傳遞關系，通常用來評價結構隔振的效果［10－11］。然而結構響應傳遞率參數化表達式中并不顯含系統極點信息，因此，長期以來并未如結構系統傳遞函數、功率譜函數等顯含系統極點信息的函數一樣，被廣泛的應用于結構模態參數辨識。

下面通過結構系統的傳遞函數、不同工況下的結構響應來說明結構響應傳遞率的性質和偽傳遞率的定義。

若定義多輸入多輸出（MIMO）的結構系統傳遞函數為H（s）∈CNo，Ni（其中No和Ni分別為輸出和輸入數），那么結構響應為：

式中，X（s）∈CNo×1為響應（輸出），F（s）∈CNi×1為激勵（輸入）。

為方便定義結構響應傳遞率，將響應分為參考響應和非參考響應，即X（s）＝［Xr（s） Xp（s）］T，其中r 和q分別表示參考點集合和非參考點集合。那么，式（1）可寫為：

式中，Xr（s）∈CNr×1，Xq（s）∈C（No－Nr）×1，Nr為參考點數，Hr（s）∈CNr×Ni，Hq（s）∈C（No－Nr）×Ni。

根據上式，文獻［12］定義了多參考點的結構響應傳遞率：

將式（2）代入式（3）得：

若激勵F（s）不為0，則有：

如果Nr＝Ni，即輸入點個數等于參考點個數，式（5）可以寫為：

如果Nr＞Ni，則可以用偽逆來代替矩陣Hr（s）的逆，則式（5）可以寫為：

如果Nr＜Ni，則無法利用如式（6）或式（7）中Hr（s）的逆或偽逆來定義多參考點的結構響應傳遞率。

類似于單參考點的結構響應傳遞率，除在極點附近，不同工況的結構響應傳遞率是不同的，多參考點的結構響應傳遞率也與不同激勵點上作用的激勵相對幅值有關，即與工況有關。

假設在Ni個輸入點上作用種工況，用l指出不同的工況Fl（s）。對于Ni個輸入點，存在Ni組獨立（線性無關）的激勵分布組合出工況Fl（s），如：

式中，Fl，j為第l個工況的第j組獨立的激勵分布，βl，j（s）為第l個工況的第j組獨立的激勵在s域上的相對幅值。將式（8）寫為矩陣形式為：

式中，Fl∈CNi×Ni為激勵分布矩陣，βl（s）∈CNi×1為激勵在s域上的相對幅值。

根據式（2），式（9）可寫為：

多參考點的系統傳遞函數也可寫為公分母形式的分式多項式參數化模型，如H（s）＝N（s）／D（s），因此，式（11）可改寫成：

如式（12）所示，多參考點傳遞函數矩陣的逆與非參考點傳遞函數矩陣相乘導致結構響應傳遞率中的極點信息被消去。因此，要用響應傳遞率進行模態參數辨識，還需要利用關于傳遞函數和傳遞率的更多性質，具體方法將在下文給出。

2 基于多參考偽傳遞率的模態參數辨識

2.1 傳遞率的參數化模型

如式（3），傳遞率與系統傳遞函數的形式一致，表示頻域的兩個向量的傳遞關系，因此可采用矩陣分式模型的形式對多參考點的傳遞率進行參數化建模。又如式（3）定義，傳遞率矩陣T（s）為（No－Nr）×Nr的矩陣，因此采用左矩陣分式多項式［13］對第l個工況下的多參考點的傳遞率建立參數化模型：

式中，

式中，Al（Ω）為（No－Nr）×（No－Nr）矩陣，Bl（Ω）為（No－Nr）×Nr矩陣，nA和nB為對應多項式階數為頻率基函數。為提高參數估計的矩陣數值條件，以z域多項式為頻率基函數，即＝e－jωkiΔT。

2.2 從傳遞率的參數化模型到模態參數

為建立參數化的傳遞率模型與傳遞率性質的關系，即如何利用參數化模型中的待估參數計算模態參數，下面利用了多輸入多輸出傳遞函數矩陣秩為1的性質完成推導。

由于多輸入多輸出傳遞函數矩陣在極點處秩為1，且在極點附近所有列與模態振型成正比［11］，因此，式（18）在極點附近有：

利用Nl種工況，結合式（13）和（19），極點λm和模態振型φm可由以下的矩陣多項式的廣義特征值求解得到。

矩陣多項式的廣義特征值需展成如下形式按照普通矩陣特征值問題進行求解：

式中，n＝max（nA，nB），Γ為方陣，稱為式（20）的高維伴隨矩陣，其中的子矩陣須為方陣。文獻［9］為滿足子矩陣為方陣，給出了約束Nl＝No／（No－Nr）。然而，這個約束在實際應用中是很難滿足的。因此，本文利用矩陣偽逆，實現廣義矩陣的方陣化。即令：

2.3 模型參數的估計

如前述，只要得到不同工況下的待估參數矩陣Al，i和Bl，i，就能夠利用式（22）計算出結構的模態參數。本節將給出待估參數的最小二乘估計方法。

根據傳遞率的定義和傳遞率的參數化模型，此參數估計問題的最小二乘花費函數可表示為：

式中，al，m和bl，m分別為Al（Ω）和Bl（Ω）第m行的行向量，Nf為頻帶內的頻率采樣點個數，最小二乘的方程誤差為：

將式（14）代入式（25），可得：

那么，最小二乘的花費函數也可以寫成相應的矩陣形式：

那么，最小二乘的正則方程為：

為減小最小二乘求解的計算量，將原最小二乘問題轉化為減縮形式。對原正則方程（式（29））進行處理：

由于若傳遞率參數化模型中矩陣Al（Ω）和Bl（Ω）同乘一個常數，而待估的傳遞率T＾l（Ω）并不會變化，因此如式（25）的傳遞率參數化模型是參數冗余的，當它被用來做最小二乘估計時導致了正則方程（式（32））僅有平凡解或無解。那么，為了消除參數化模型中的參數冗余，需對待估參數矩陣Θl中的某一個子塊進行約束。本文將Ψl中Al，nA約束為單位矩陣，那么，約束后的正則方程為：

解方程式（33）得到Ψl，進而通過（31）計算出Ξl，至此待估參數Al，j和Bl，j就通過最小二乘估計得到了。

3 算例一

本節利用一個4自由度的彈簧－阻尼－質量系統對本文研究的方法進行驗證。4自由度的彈簧－阻尼－質量系統如圖1所示。

圖1 四自由度的彈簧－阻尼－質量系統Fig．1 4－DOF system

作用在如圖1所示的激勵為理想白噪聲，即激勵的頻譜是頻率的常函數。工況數Nl＝2，兩種工況的相對幅值如圖2所示。

圖2 兩種工況激勵的相對幅值Fig．2 Relative Amplitude of excitation for the two load cases

圖3 穩定圖（算例一）Fig．3 Stabilization diagram（Example 1）

響應參考點數Nr＝2，分別為響應點1和2。令nA＝nB，并從4取值到20構造穩定圖，如圖3所示（o為不穩定，f、d和v分別為頻率、阻尼比或振型穩定，s為頻率、阻尼比和振型都穩定）。

由如圖3所示的穩定圖可以看出4自由度系統的4階模態（連續的“s”表示），辨識得到的模態頻率和阻尼比與理論值對比如表1所示。

表1 辨識得到的模態頻率和阻尼比與理論值對比（算例一）Tab.1 Comparison between estimated and theoretical modal frequency and damping ratio（Example 1）

辨識得到的模態振型與理論值對比如表2所示。（由于僅輸出模態參數辨識無法得到模態縮放因子，因此為方便對比，辨識值和理論值都進行最大值歸一化）。

表2 辨識得到的模態振型與理論值對比Tab.2 Comparison between estimated and theoretical mode shapes

由表1和表2可以看出，在沒有噪聲干擾的情況下，本文給出的辨識方法能夠十分精確的辨識出結構的模態頻率和阻尼比。

4 算例二

為驗證方法對實際結構模態參數辨識的有效性和正確性，本節給出一個有限元數值算例。

4.1 模型

有限元模型模擬一個在空間飛行的火箭結構，邊界條件為自由－自由。有限元模型如圖4所示，結構采用梁單元模擬，錐形整流罩部分采用變截面梁單元實現。結構總長20 m，直徑1．4 m，分為20個單元，共21個節點。模型阻尼矩陣由Rayleigh阻尼公式得到，其中α＝1．0，β＝0．000 25。

方便問題描述和結果分析，本算例僅研究橫向一個面內的響應和模態特性，即僅考慮如圖4所示的y方向的激勵以及y方向的平動和繞z方向的轉動自由度。

圖4 火箭有限元模型（三維顯示）Fig．4 Finite element model of a rocket（3D display）

4.2 工況與響應仿真分析

假設作用在結構上的有兩種工況：一種是作用在整流罩部分節點上的隨機載荷，模擬氣動力的隨機波動部分；另一種是作用在火箭尾部上的發動機擺動在y方向的隨機分量。通常液體火箭的液體燃料、發動機擺動等導致了周期載荷，為考慮這樣的周期載荷，本算例模擬其極端情況，即考慮兩個單頻率分量的諧波作用在結構上，分別是6 Hz和46．5 Hz，幅值為隨機載荷方差的103（僅為考驗本文的辨識算法，實際情況不會達到此量級）。利用數值積分得到結構在兩種工況下的加速度響應，從而獲取結構響應的頻譜信息X（Ω）。兩個工況下1號點頻率響應如圖5所示。

圖5 兩個工況下1號點頻率響應（上：工況1；下：工況2）Fig．5 Frequency responses of Point 1 （upper：load case 1；lower：load case 2）

如圖5所示，結構的頻率響應中在6 Hz和46．5 Hz處存在由于諧波引起的峰（黑色虛線處）。若采用傳統的基于功率譜函數的僅輸出辨識方法進行模態參數辨識將無可避免的在這兩個峰值處產生虛假模態。

4.3 辨識結果

通過本文給出的基于結構響應傳遞率的模態參數辨識方法，僅利用結構的響應信息進行模態參數辨識，其中，辨識帶寬選為0－120 Hz，21個點的響應皆用于辨識（No＝21），響應參考點數Nr＝2，分別為響應點3和響應點4。令nA＝nB，并從8取值到20構造穩定圖（如圖6所示）。

如圖6所示，穩定圖中的4列連續的穩定極點指出了帶寬內的四階模態。模態頻率和阻尼比辨識結果如表3所示。

圖6 穩定圖（算例二）Fig．6 Stabilization diagram（Example 2）

表3 辨識得到的模態頻率和阻尼比與理論值對比（算例二）Tab.3 Comparison between estimated and theoretical modal frequency and damping ratio（Example 2）

由于僅輸出的模態參數辨識無法得到模態縮放因子，因此為方便對比，辨識值和理論值都對其最大值進行歸一化。模態頻率和阻尼比辨識結果如圖7所示（線為理論值，點為辨識結果）。

圖7 辨識得到的模態振型與理論值對比Fig．7 Comparison between estimated and theoretical mode shapes

由表3和圖7可以看出，在實際激勵的情況下，即使激勵中含有明顯的諧波分量，本文提出的辨識方法能夠很好的辨識出結構的模態頻率、阻尼比和模態振型，且辨識結果不受諧波影響，即辨識結果中不含有諧波造成的虛假模態。

5 結 論

本文推導了多參考點的結構響應傳遞率表達式，建立了其左矩陣分式多項式的參數化模型，進而給出了辨識問題的最小二乘估計格式，利用Jacobi矩陣的分塊性質對最小二乘問題矩陣形式完成了縮減，降低了計算量，通過高維伴隨矩陣方法解決了矩陣多項式的特征值求解問題，即多參考多輸出的模態參數求解問題，以及通過矩陣偽逆解決了現有方法中載荷工況數與響應點數和參考點數的約束問題。最后，本文通過兩個數值算例對所提方法進行了驗證，辨識結果表明：本文給出的辨識方法能夠很好的辨識出結構的模態頻率、阻尼比和模態振型，且能夠很好的避免激勵中含有的諧波分量對辨識結果造成的影響，解決了傳統僅輸出模態參數辨識中激勵非白對辨識結果造成的影響問題。

［1］李火坤，練繼建．高拱壩泄流激勵下基于頻域法的工作模態參數識別［J］．振動與沖擊，2008，27（7）：149－53．

LI Huo-kun，LIAN Ji-jian．Operational modal parameter identification of high arch dam under flood discharge excitation based on frequency domain decomposition method［J］．Journal of Vibration and Shock，2008，27（7）：149－53．

［2］Peeters B，De Roeck G．Reference-based stochastic subspace identification for output-only modal analysis［J］．Mechanical Systems and Signal Processing，1999，13（6）：855－878．

［3］Pickrel C R，White P J．Flight flutter testing of transport aircraft：in-flight modal analysis［C］／／proceedings of the IMAC 21，the International Modal Analysis Conference，Kissimmee（FL），USA，F，2003．

［4］呂中亮，楊昌棋，安培文，等．多點激勵模態參數識別方法研究進展［J］．振動與沖擊，2011，30（1）：197－203．

Lü Zhong-liang，YANG Chang-qi，AN Pei-wen，et al．Progress on modal parameter identification with multipleexcitation［J］．Journal of Vibration and Shock，2011，30 （1）：197－203．

［5］Hermans L，Van Der Auweraer H，Hatami A，et al．InflightModal testing and analysis of ahelicopter［C］／／proceedings of the Proceedings of the 17th International Modal Analysis Conference（IMAC17），Kissimmee，FL，USA，F，1999．

［6］Devriendt C，Guillaume P．The use of transmissibility measurements in output-only modal analysis［J］．Mechanical Systems and Signal Processing，2007，21（7）：2689－2696．

［7］Devriendt C，Guillaume P．Identification of modal parameters from transmissibility measurements［J］．Journal of Sound and Vibration，2008，314（12）：343－356．

［8］Devriendt C，De Sitter G，Vanlanduit S，et al．Operational modal analysis in the presence of harmonic excitations by the use of transmissibility measurements［J］．Mechanical Systems and Signal Processing，2009，23（3）：621－35．

［9］Weijtjens W，Sitter G D，Devriendt C，et al．Operational modalanalysisbasedonmultivariabletransmissibility functions：revisited［C］／／proceedings of the IMAC-XXXI，International Modal Analysis Conference，Garden Grove，California，United States，F，2013．

［10］吳川，劉學文，黃醒春．短型浮置板軌道系統隔振性能研究［J］．振動與沖擊，2008，27（8）：74－76．

WU Chuan，LIU Xue-wen，HUANG Xing-chun．Study on effectiveness of short floating slab track system［J］．Journal of Vibration and Shock，2008，27（8）：74－76．［11］涂奉臣，張洋，陳照波，等．采用磁流變阻尼器的整星隔振平臺試驗研究［J］．振動與沖擊，2012，31（21）：33 －36．

TU Feng-chen，ZHANG Yang，CHEN Zhao-bo，et al．Tests for a whole-spacecraft vibration isolation platform using a magnetorheological damper［J］．Journal of Vibration and Shock，2012，31（21）：33－36．

［12］Maia N M M，Silva J M M，Ribeiro A M R．The transmissibility concept in multi-degree-of-freedom systems ［J］．Mechanical Systems and Signal Processing，2001，15 （1）：129－137．

［13］Allemang R J，Brown D L．A unified matrix polynomial approach to modal identification［J］．Journal of Sound and Vibration，1998，211（3）：301－322．

Output-only structural modal parameter estimation under non-white excitations based on response transmissibility

ZHOU Si-da，LIU Li，YANG Wu，MA Zhi-sai
（School of Aerospace Engineering，Beijing Institute of Technology，Beijing 100081，China）

In cases of output-only modal parameter estimation for operational engineering structures，such as aerospace structures，the non-white excitations lead to non-physical modes being not able to be identified．For these cases，the recently developed response transmissibility-based methods of modal parameter estimation were studied and an improved poly-reference response transmissibility-based output-only approach of modal parameter estimation was proposed．Firstly，the parametric model of the poly-reference response transmissibility was introduced，it was presented by the left matrix fraction polynomials，the formulations of the least square estimation for unknown parameters were proposed，and the normal equation of the least square estimation was reduced by using the block property of the corresponding Jacobian matrix，it largely decreases the computational costs．Secondly，based on the high-dimensional accompany matrix method，the eigenvalue problem of the matrix polynomials was solved，it led to the solution to the modal parameters．In addition，the original constraints between the number of loading conditions and the number of references were removed by using the pseudo inverse of matrix．Finally，two numerical examples validated the proposed method．The results illustrated that the proposed method can be used to estimate the modal parameters very well and has the advantage of avoiding the influence of harmonics，so the problem due to non-white excitations is solved．

modal parameter estimation；output-only；response transmissibility；non-white noise excitation

N945．14；TH113．1；O324；TN911

A

10．13465／j．cnki．jvs．2014．23．009

北京理工大學基礎研究基金（20120142009）

2013－09－12 修改稿收到日期：2013－12－12

周思達男，博士，講師，1984年7月生

劉莉女，教授，博士生導師，1964年生