基于Hankel矩陣聯合近似對角化的結構模態識別*

程 琳， 楊 杰， 鄭東健， 任 杰

(1. 西安理工大學西北旱區生態水利工程國家重點實驗室培育基地 西安, 710048)(2. 河海大學水文水資源與水利工程科學國家重點實驗室 南京, 210098)

?

基于Hankel矩陣聯合近似對角化的結構模態識別*

程 琳1,2， 楊 杰1， 鄭東健2， 任 杰1

(1. 西安理工大學西北旱區生態水利工程國家重點實驗室培育基地 西安, 710048)(2. 河海大學水文水資源與水利工程科學國家重點實驗室 南京, 210098)

針對基于二階盲辨識(second order blind identification，簡稱SOBI)的模態參數識別方法存在的不足，提出了一種基于Hankel矩陣聯合近似對角化(Hankel matrix joint approximate diagonalization，簡稱HJAD)技術的結構運行模態分析(operational modal analysis，簡稱OMA)的新方法。該方法通過對隨機子空間類模態識別方法常用的Hankel矩陣進行聯合近似對角化，以分離各階模態響應，進行模態識別。與基于SOBI的模態識別方法相比，在具體實施過程中，僅需要在分析數據中添加與實測振動響應對應的時間延遲的數據，實現難度較小。數值算例和物理模型試驗的分析結果表明，所提出的基于HJAD技術的結構運行模態分析方法，不僅具有魯棒性強和計算效率高的優點，還可以克服傳統的基于SOBI的模態識別方法的模態識別能力受測點數目限制的問題。

運行模態分析；盲源分離；二階盲辨識；Hankel 矩陣聯合近似對角化；振動臺試驗

引 言

根據實測振動數據來識別結構的模態參數并計算模態指標，再應用相關模態指標來評價結構的健康狀態，是基于振動的結構健康診斷方法的基本思路[1]。為了研究大壩等大型工程結構的振動特性，需要對其施加一定的振動激勵。傳統的大型結構原型動力試驗[2]需要施加人工產生的激勵，如起振機和爆破等，不僅耗費巨大，可能對結構產生破壞，而且人工激勵與結構正常運行時可能遇到的振動激勵差異很大。因此，研究各種環境激勵，包括地震、脈動風、脈動水壓力和交通振動等作用下大型工程結構的模態識別方法顯得十分必要。

環境激勵下結構模態參數的識別方法，也被稱為結構的運行模態分析[3-4]，一般僅根據結構的振動響應來識別模態參數。目前常用的運行模態分析方法主要包括頻域類方法、時域類方法、時頻分解算法和SOBI法等。在上述各種方法中：頻域類方法的識別精度有待提高；時域類方法常常要用到矩陣分解技術，計算量很大；時頻分解算法在計算中存在參數選擇困難和邊界效應等問題。模態識別的SOBI方法[5]由于其在模態參數識別方面表現出魯棒性強和計算效率高等優點，受到了廣泛的關注[6-7]。然而，基于SOBI的模態識別方法，只能識別出與測量自由度數相同的模態階數，模態識別能力受到測點數目的限制，并且該方法在系統定階方面也存在困難，限制了其在實際工程中的應用。

為了克服傳統SOBI模態識別方法存在的模態識別能力受測點數目限制的問題，筆者提出了HJAD技術的結構運行模態分析方法，并將該方法與穩態圖相結合，以解決系統定階和剔除虛假模態的問題。通過一個數值算例和一個物理模型試驗，驗證了該方法在識別精度、魯棒性和計算效率等方面的表現。

1 環境激勵下基于SOBI的模態識別

1.1 自然激勵技術

各種環境激勵本身的性質復雜，實際應用中常常將其簡化成隨機的帶限白噪聲。隨機環境激勵下，結構的受迫振動響應也是隨機的，而結構的模態參數識別一般是基于其自由振動響應或脈沖響應的。因此，需要研究從隨機受迫振動響應中提取結構自由響應或脈沖響應的方法。

n自由度結構的振動平衡方程可以表達為

(1)

其中：M，K和C分別為系統的質量、剛度和阻尼矩陣；f(t)為外部激勵；x(t)為位移響應。

環境激勵一般不可直接測量，并且與噪聲相比，二者的幅值常常沒有顯著的差異。因此，環境激勵下結構的振動一般可以表達為以下的隨機離散狀態空間模型[8]

(2)

其中：zk為離散的系統狀態變量；yk為系統觀測量；A為離散系統矩陣；G為觀測矩陣；wk為環境激勵和系統誤差；vk為觀測噪聲。

如果觀測噪聲也是隨機白噪聲，根據式(2)進行推導，可以得到各測點振動觀測間的相關函數矩陣的表達式

(3)

由式(3)可以看出，白噪聲激勵下結構各測點振動觀測間的相關函數與結構的脈沖響應以及自由振動響應類似，僅與結構的初始狀態和系統參數(模態參數)有關。因此，當環境激勵和觀測噪聲接近白噪聲時，可以采用結構各測點響應間的互相關函數代替脈沖響應或自由振動響應來進行模態識別，這就是自然激勵技術(nature excitation technique，簡稱NExT)[9]。

1.2 基于SOBI的模態識別

(4)

如果結構的各階模態響應q(t)相互獨立，結構的自由振動響應可以看作是各階模態響應的線性混合，各階模態響應被看作是虛擬源信號，混合矩陣即為L。這時可以采用盲源分離(blind source separation，簡稱BSS)方法來求解模態識別問題。在BSS的各種實現算法中，基于聯合近似對角化(joint approximate diagonalization，簡稱JAD)技術的SOBI算法[11]在實際應用中表現出魯棒性強和計算效率高等優點。下面基于該方法研究結構模態參數的識別問題。

假定觀測噪聲是高斯白噪聲，這時可以得到各觀測信號間的相關函數矩陣的表達式為

(5)

忽略阻尼影響，各階模態響應相互獨立，Rqq(τi)是一個對角陣。從式(4)可以看出，任意一個時延τi(i= 1, 2,…,S)，理論上應該存在矩陣P=L+使Ryy(τi)對角化，但實際的觀測信號中總存在噪聲，并且由于阻尼的存在，源信號也不能精確滿足相互獨立的條件。這時，對一組不同時延的矩陣進行JAD時，可以轉而求解魯棒性更強的優化問題min:Joff(U)，其中優化的目標函數為

(6)

其中：運算符號off(·) 表示提取矩陣所有非主對角元素；U為需要通過優化得到的對角化矩陣；W為白化矩陣。

白化的目的是使各混合通道的信號經過變換后相互獨立，具有單位方差，以避免優化過程出現零解。一般采用主成分分析(principle component analysis,簡稱PCA)方法來確定白化矩陣W

(7)

其中：矩陣D為以Ryy(0)的特征值為對角元素的對角矩陣；F為其特征向量組成的矩陣。

用某一種數值方法對以上不同時延的協方差矩陣實施JAD，即求解式(6)的優化問題。常用的JAD實現技術包括雅克比方法、交替最小二乘法、并行因子分析、子空間擬合技術和旋轉適量法等。這時，可以得到在JAD意義上的廣義的正交矩陣U(而不是僅使某一個時滯協方差矩陣對角化)。這時，分離矩陣為

P=UHW

(8)

其中：H為共軛轉置算子，對于實數而言即轉置算子T。

源信號，即模態響應的估計為

(9)

混合矩陣為

(10)

其中：W+為W矩陣的廣義逆矩陣。

式(5)是根據結構的自由響應來進行推導的，而實際工程中一般只能得到環境激勵下結構的受迫振動響應。根據NExT可以采用環境激勵下結構受迫振動響應間的相關函數矩陣來代替式(5)所示的結構自由振動間的相關函數矩陣，并進行JAD。分離得到的每一階模態響應代表一個單自由度振動系統，可以采用單模態識別法(如最小二乘擬合法、時間序列法等)來識別對應的模態參數。模態振型可以從混合矩陣L中得到。

2 環境激勵下基于HJAD的模態識別

對于實際工程結構而言，結構的模態階數n是未知的，而且測點數目l常常遠小于結構模態階數n。一些分析結果表明，結構高階模態參數對結構損傷更為敏感。因此，研究采用有限的振動測點，識別出盡可能多階數的結構模態參數，具有十分重要的理論意義和應用價值。采用上述基于SOBI的模態參數識別方法時，最多只能識別與測點數目相同的l階模態。模態識別能力受到測量自由度數l的約束，極大限制了SOBI在實際工程結構運行模態分析中的應用，因此需要研究改進的方法。

定義由實測振動響應及其時延數據形成的向量

(11)

其中：p為時間延遲。

p的設定原則是要使rank(H(τ)) ≥n成立。這時可以采用式(12)計算Hankel矩陣[12]

(12)

對于結構的自由振動，Rqq(τ)=ΣτRqq(0)成立。這時，根據式(5)對式(12)定義的Hankel矩陣進行推導可以得到

(13)

以上Hankel矩陣的分解形式與協方差矩陣的分解形式(5)是一致的。根據式(5)所示的分解形式，協方差函數矩陣可以采用JAD技術來進行聯合近似對角化。對于Hankel矩陣，采用式(13)的分解形式，也同樣可以對一組不同時延的Hankel矩陣，即H(τ1),H(τ2),…,H(τS)采用JAD技術實現近似對角化，并求解模態參數。對一組不同時延的Hankel矩陣進行聯合近似對角化，可以轉化成求解以下的優化問題min:Joff(U)，其中的優化目標函數Joff(U)為

(14)

3 數值驗證

對圖1所示的4自由度集中質量系統的模態參數進行識別。其參數如下：m1=m2=m3= 1 kg；m4= 0.9 kg；k1=k3= 7 kN/m；k2=k4= 8 kN/m；c1=c2= 0.6 (N·s2)/m；c3=c4= 0.55 (N·s2)/m，系統的初始狀態設為零。基礎環境激勵ag(t)采用高斯白噪聲來模擬。結構的振動觀測物理量是絕對加速度，采用加速度響應的計算結果加上一定水平的噪聲來模擬。結構振動響應的計算采用Runge-Kutta 算法來實現，模擬加速度觀測的采樣頻率為100 Hz。結構自振頻率的理論值為5.172，13.524，21.872 Hz和26.019 Hz；阻尼比的理論值分別為0.144%，0.346%，0.510%和0.620%。對于無噪聲的觀測數據，分別采用1，2，3和4個測點來識別結構的模態參數，識別結果見表1。傳統的基于SOBI的模態識別方法僅能識別出與測點數目相同的模態階數，而采用HJAD方法，即使在僅有1個測點的情況下，也可以識別出結構的各階自振頻率和阻尼比。根據1個測點(m4)的振動響應，采用基于HJAD技術的模態識別方法分離得到的結構前2階模態響應及其功率譜如圖2所示。

圖1 4自由度集中質量系統Fig.1 Concentrated mass system of 4 degrees of freedoms

Tab.1 Identification results of nature frequencies and damping ratios

方法4個測點3個測點2個測點1個測點f/Hz阻尼比/%f/Hz阻尼比/%f/Hz阻尼比/%f/Hz阻尼比/%SOBI5．1730．1465．1730．146————13．5220．34813．5180．34013．5690．34113．4860．34321．8490．51321．8320．50221．8390．505——25．9670．598——————HJAD5．1720．1455．1750．1455．1770．1445．1770．14613．5240．34713．5220．34313．5230．34513．5230．34921．8710．51321．8490．50221．8520．51221．8510．51026．0190．61225．9670．62725．9880．61126．1830．625

圖2 模態響應及其功率譜Fig.2 Modal responses and power spectrums

采用測點2，3和4對應的振動響應數據進行結構振型的識別，根據振型的識別值和理論值計算的的模態置信因子(modal assurance criteria，簡稱MAC)見表2。

表2 振型的識別精度

所有振型都采用第4個測點的坐標值進行歸一化。從表2可以看出，HJAD技術可以識別4階振型，而SOBI算法僅能識別3階振型。由于結構的損傷常常在高階模態上表現更為顯著，采用有限的測點，準確地識別結構盡可能多階的模態參數，可以對結構的損傷診斷提供很大的便利。為了驗證筆者提出的模態識別方法的魯棒性，根據4個測點的振動響應在不同信噪比(signal noise ratio，簡稱SNR)水平下，采用特征系統實現算法(eigensystem realization algorithm，簡稱ERA)、隨機子空間方法(stochastic subspace identification，簡稱SSI)和基于HJAD技術的方法進行模態參數識別，自振頻率的識別誤差隨SNR的變化見圖3。

圖3 自振頻率的識別結果對比Fig.3 Comparison of the identification results of nature frequencies

由圖3可以看出，隨著噪聲水平的提高(SNR減小)，基于HJAD技術的模態識別方法的識別結果更穩定。在高噪聲水平下，其模態參數的識別精度更高。

4 試驗驗證

根據某混凝土重力壩的一個擋水壩段，按照幾何相似的原則(比例為1∶160)制作試驗模型，并進行小型振動臺試驗。試驗模型的材料采用混凝土砂漿，質量密度為2 100 kg/m3，水灰比為水∶細砂∶水泥 = 1∶2.9∶2.09。試驗中環境激勵通過振動臺產生的帶限白噪聲來模擬，考慮到模型自振頻率的范圍和采樣的要求，設定帶限白噪聲激勵的頻率范圍為100～2 000 Hz，采樣頻率為4 kHz。模型共布置6個單向的壓電式加速度傳感器，如圖4所示，傳感器等間距布置，間距為10cm。采用DSPACE數據采集系統來測量模型的加速度響應，6個測點實測的加速度響應如圖5所示。

圖4 試驗模型和試驗布置Fig.4 The experimental model and the layout of the experiment

圖5 模型的實測振動響應Fig.5 The measured vibration response of model

采用帶通濾波器來濾除低頻的趨勢分量和無效的頻率分量。頻率的通帶根據有限元計算的結果來確定，然后采用小波閥值去噪法進行去噪處理。根據處理后的振動響應分別采用SOBI方法、SSI方法和HJAD方法進行模態參數識別的結果，以及采用有限單元法(finite element method，簡稱FEM)計算的結果見表3。由表3可以看出，在前7階模態中，SOBI方法僅能識別4階模態，而SSI和HJAD法的識別結果相近，可以識別出6階模態。各種方法都未識別出第4階模態，可能是在環境激勵下，該階模態響應太微弱，無法識別。計算中，SSI識別方法的計算時間大約為HJAD方法的5倍。采用SSI方法、HJAD方法識別的振型與FEM計算結果的對比見圖6。自振頻率和振型識別結果與FEM計算結果有一定的誤差，主要是由有限元模擬誤差(材料和邊界條件等方面)和計算識別誤差等因素造成的。

表3 自振頻率的識別結果

圖6 前2階振型識別結果的對比Fig.6 The comparison of the identification results of the first two modal shapes

5 結束語

基于HJAD技術的模態識別方法既可以看作是在時域模態參數識別方法中，引入計算效率高和魯棒性強的JAD技術，也可看作是將SOBI算法中采用的各測點響應間的協方差函數矩陣，用Hankel矩陣來代替，以實現系統擴階。實例分析表明，基于HJAD技術的結構模態參數識別方法魯棒性強，模態識別結果更穩定，精度更高，克服了傳統的基于SOBI的模態識別方法的模態識別能力受測點數目限制的問題，并且算法的計算效率也很高。這些優點為其在實際工程中的應用提供了良好的前景。

[1] Farrar C, Doebling S, Nix D. Vibration-based structural damage identification [J]. Philosophical Transactions of the Royal Society of London Series A: Mathematical, Physical and Engineering Sciences, 2001, b(359): 131-149.

[2] 張光斗, 張楚漢, 李未顯, 等. 泉水拱壩的振動測量與分析[J]. 中國科學A輯, 1986, 1: 100-112.

Zhang Guangdou, Zhang Cuhan, Li Weixian, et al. Vibration measurement and analysis of Quanshui arch dam [J]. Science China, Series A, 1986, 1: 100-112.(in Chinese)

[3] 張建偉, 李火坤, 練繼建, 等. 基于環境激勵的廠房結構損傷診斷與安全評價[J]. 振動、測試與診斷, 2013, 33(4) : 670-674.

Zhang Jianwei，Li Huokun, Lian Jijian, et al. Damage diagnosis and safety evaluation of underground powerhouse under ambient excitation [J]. Journal of Vibration, Measurement & Diagnosis, 2012, 32(4): 670-674.(in Chinese)

[4] 續秀忠, 華宏星, 陳兆能. 基于環境激勵的模態參數辨識方法綜述[J]. 振動與沖擊，2002, 21(3): 1-5.

Xu Xiuzhong, Hua Hongxing, Chen Zhaoneng. Review of modal identification method based on ambient excitation [J]. Journal of Vibration and Shock, 2002, 21(3): 1-5.(in Chinese)

[5] 曹軍宏, 韋灼彬, 劉樹勇. 改進型盲源分離在結構模態識別中的應用 [J]. 振動、測試與診斷, 2013, 33(4): 689-693.

Cao Junhong, Wei Zhuobin, Liu Shuyong. Application of improved blind source separation in modal parameter identification [J]. Journal of Vibration, Measurement & Diagnosis, 2013, 33(4): 689-693.(in Chinese)

[6] Antoni J, Brun S. Special issue on blind source separation: editorial[J]. Mechanical Systems and Signal Processing, 2005, 19:1163-1165.

[7] McNeill S. Modal identification using blind source separation techniques [D].Houston, Texas:University of Houston, 2007.

[8] Van Overschee P, De Moor B. Subspace identification for linear systems: theory, implementation, applications [M]. Dordrecht, the Netherlands: Kluwer Academic Publishers, 1996: 57-73.

[9] James G H, Carne T G, Lauffer J P. The natural excitation technique for modal parameter extraction from operating structures [J]. Modal Analysis: International Journal of Analytical and Experimental Modal Analysis, 1995, 10(4): 260-277.

[10]McNeill S I. An analytic formulation for blind modal identification [J]. Journal of Vibration and Control, 2012, 18(14): 1-11.

[11]Belouchrani A, Abed-Meraim K, Cardoso J F, et al. A blind source separation technique using second-order statistics [J]. IEEE Transactions on Signal Processing, 1997, 45 (2): 434-444.

[12]Antoni J, Chauhan S. A study and extension of second-order blind source separation to operational modal analysis [J]. Journal of Sound and Vibration, 2012, 18 (14): 2111-2121.

10.16450/j.cnki.issn.1004-6801.2016.06.018

*國家自然科學基金資助項目(51409205)；博士后自然科學基金資助項目(2015M572656XB)；陜西省重點科技創新團隊資助項目(2013KCT-015)；水文水資源與水利工程科學國家重點實驗室開放研究基金資助項目(2014491011)

2014-05-17；

2015-10-29

TV32；TU311

程琳，男，1986年7月生，講師。主要研究方向為水工結構安全監控和無損檢測。曾發表《Two online dam safety monitoring models based on the process of extracting environmental effect》(《Advances in Engineering Software》2013, Vol.57)等論文。 E-mail：wanmu2381@163.com