基于ICEEMDAN的連續梁橋車致振動信號的HHT分析

邢世玲 ，呂雙雙 ，朱利明 ，張 佳

（南京工業大學交通運輸工程學院，江蘇 南京 210009）

為了保證橋梁的服役安全，橋梁使用一定年限后，必須對橋梁進行定期檢測來監測橋梁結構的健康狀況并進行養護管理.損傷識別與定位是結構損傷檢測系統設計的關鍵.我國有大量待檢測的常規橋梁，而這種橋梁往往原始資料缺失，無法獲取橋梁結構的初始狀態.所以，不基于模型的損傷識別方法非常有必要和意義.無模型識別一般可分為時域識別、頻域識別和時頻分析.傅里葉譜分析是常用的頻域分析方法，時頻分析方法包括小波變換、小波包變換和Hilbert-Huang變換（Hilbert-Huang transform，HHT）.近十幾年的基于移動車輛荷載作用的橋梁損傷識別研究中，采用最多的時頻分析方法是小波變換[1-4]，采用HHT法的還較少[5-8].HHT是Huang等[9]1998年提出的，由經驗模式分解（empirical mode decomposition，EMD）和 Hilbert譜分析（Hilbert spectral analysis，HSA）構成，其中 EMD是核心，HHT法是一種比傅里葉變換及小波變換等更具適應性的時頻局部化分析方法.因為HHT具有完備的自適應性，其EMD過程無需提前設定基函數，克服了小波分析方法依賴主觀經驗的缺點.但EMD方法由于自身計算理論的缺陷，在分解過程中容易出現模態混疊現象，影響分解效果，限制了HHT方法在橋梁工程中的應用.所以，EMD自被提出以來一直在改進，經歷了EMD[9]、集合經驗模式分解（ensemble empirical mode decomposition，EEMD）[10]、完備集合經驗模式分解（complementary ensemble empirical mode decomposition，CEEMD）[11]、帶有自適應噪聲的完備集合經驗模式分解（complete ensemble empirical mode decomposition with adaptive noise，CEEMDAN）[12]和改進的帶有自適應噪聲的完備集合經驗模式分解（improved complete ensemble empirical mode decomposition with adaptive noise，ICEEMDAN）[13]的發展過程.

相較于 EMD、EEMD、CEEMD和 CEEMDAN，ICEEMDAN方法對解析信號分解效果最理想[13]，其在橋梁工程領域的應用還很少.下面將其應用于橋梁振動信號的處理中，對移動車輛荷載作用下的連續梁橋的數值模擬信號和實測信號進行HHT分析，提取模擬信號的時頻特征，并以此為基礎對損傷橋的實測信號的頻譜特征進行分析和討論.當前基于移動荷載的損傷識別一般采用簡單的移動常力或移動質量來模擬移動車輛荷載，為了更加接近真實情況，在計算移動車輛荷載作用下的橋梁結構動響應信號時，將考慮隨機橋面不平度的影響.

1 模擬和實測振動信號的獲取

1.1 基本信息

1.1.1 橋梁概況

選一鋼筋混凝土曲線連續梁橋為對象，其基本信息如下：1）跨徑布置為（18.0+23.5+18.0）m（以橋面中心線弧長為參考），曲線半徑為75 m，跨中橫截面形式及細部尺寸見圖1，采用C40混凝土；2）連續梁橋兩端采用雙支座布置，中間均采用單支座且無偏心布置.該橋在外觀檢查中發現橋身存在多處裂縫，跨中段存在多條底板延伸至腹板的橫向裂縫，部分裂縫寬度超過0.1 mm.

圖1 跨中橫斷面（單位：cm）Fig.1 Cross section of mid-span (unit：cm)

1.1.2 數值模擬信號計算方法

本文將車輛和橋梁分離成兩個相對獨立的子系統，首先以車輛振動系統為對象，以橋面不平度的功率譜密度（power spectral density，PSD）作為輸入，推導求得輪胎沖擊荷載的PSD，利用此PSD采用諧波疊加法直接構建車輛與橋面之間的接觸力數據序列；然后以橋梁振動系統為對象建立橋梁有限元模型，將移動車輛視為運動的點源荷載作用于橋梁系統，采用直接積分法求解得到橋梁動響應數據，即模擬信號數據序列.

計算時車輛采用“剛體-彈簧-阻尼”構成的系統進行模擬，車輛計算模型見圖2，圖中：Z為車體沉浮運動位移；V為車輛行駛速度；mv為車體質量；mt1、mt2為簧下質量；ks1、ks2為懸掛剛度；cs1、cs2為懸掛阻尼；kt1、kt2分別為前、后軸輪胎剛度；θ為車體俯仰（或稱點頭）運動位移；Zt1、Zt2分別為前、后軸簧下質量沉浮運動位移；Lc為前后軸之間的距離；c1、c2為系數；c1Lc、c2Lc分別為前、后軸距車體質心的距離.各參數取值見表1（表中Jθ為車體的俯仰運動轉動慣量），其余計算參數取值為：橋面不平度系數為16 × 10?6；空間頻率為0.011～2.830 m?1；空間頻率劃分為211份；標準空間頻率 0.1 m?1.作用于橋面的車輛輪胎動載的具體計算方法見文獻[14].V=10 km/h時，計算得到的一組作用于橋面的輪胎動載數據序列如圖3所示.

圖2 兩軸車輛計算模型Fig.2 Two-axis vehicle model

圖3 作用于橋面的移動荷載Fig.3 Moving load data sequence acting on deck

表1 車輛模型的計算參數Tab.1 Calculation parameters of vehicle model

建立橋梁有限元模型，以移動的輪胎動載作為激勵源，計算時采用Newmak-β直接積分法，結構阻尼比取0.05，具體方法和過程見文獻[14].橋梁的損傷用剛度折減的方法來模擬，即模型中通過減小損傷單元的彈性模量值來實現.

1.1.3 實橋跑車試驗設計

在橋面無障礙的情況下組織跑車試驗，試驗工況為一輛 30 t的汽車分別以 10、20、30 km/h 的速度勻速駛過大橋.測點布置和編號見圖4，共布置5個測點（按行車方向在圖中依次編號1、2、3、4和5），分別于邊跨的跨中，中跨的1/4、1/2、3/4位置布置豎向拾振器，橫橋向布置在防護欄內側；圖中：①～④為沿行車方向橋墩編號.

圖4 跑車試驗測點布置及編號（單位：cm）Fig.4 Measuring points layout and number in running tests (unit：cm)

1.2 模擬信號和實測信號

1.2.1 模擬信號

在第二跨跨中5 m梁長范圍內設置損傷（彈性模量取值減少40%），速度為10 km/h時計算得到的測點1的加速度模擬信號如圖5所示.圖中0時刻對應車輛前輪駛入橋面時間，在不同速度下，車輛前、后輪到達特征截面位置的時間匯總于表2.

表2 移動車輛前、后輪到達特征截面時間Tab.2 Time when front and rear wheels of moving vehicle reach the characteristic sections s

圖5 加速度模擬信號Fig.5 Acceleration simulated signal

對比表2和圖5可看出，模擬加速度信號的振幅較大值出現在車輛后輪作用于橋跨中間梁段，振幅最小段對應于車輛后輪移動到橋墩支點附近.

1.2.2 實測信號

實際橋梁在跑車試驗中行駛速度為30 km/h工況下測點1處豎向拾振器采集的實測加速度記錄如圖6所示.

圖6 加速度實測信號Fig.6 Acceleration measured signal

2 基于ICEEMDAN的信號分解

對加速度模擬信號和加速度實測信號進行ICEEMDAN分解，可得到多個內模式函數（intrinsic mode function，IMF）.關于 ICEEMDAN 算法理論和編程見文獻[15].速度為10 km/h時非損傷橋梁測點1處模擬信號的ICEEMDAN的分解結果見圖7，分解出來的系列IMF分量按頻率從高到低依次表示為imf1，imf2，···，imf10（余項）.

圖7 ICEEMDAN 分解Fig.7 ICEEMDAN decomposition

損傷橋模擬信號與實際橋梁實測信號的ICEEMDAN分解過程類似，且分解出的IMF數量與實測信號相當.

3 振型振動分量識別

由結構動力學理論知識可知，在不考慮初位移和初速度的情況下（零初始條件），結構強迫振動的解主要由按結構固有頻率振動的自由振動部分和按激勵力頻率振動的純強迫振動部分組成.這里純強迫振動部分激勵力頻率主要和車輛自振頻率相關.所以，圖7中的IMF分量系列里應包含了振型振動分量和純強迫振動分量，其中與振型振動對應的IMF分量的頻帶應該在結構自振頻率附近，與純強迫振動對應的IMF分量的頻帶應在車輛自振頻率附近.下面首先解決與振型振動對應的IMF分量的識別問題.

3.1 無損傷橋梁的模擬信號

通過有限元模型的特征值分析，計算得到橋梁前三階豎彎自振頻率分別為5.815、8.910、10.303 Hz.若把車輛前、后軸視為兩個獨立的振動系統，前、后軸豎向振動頻率分別約為3.900、3.200 Hz[14].

計算各IMF分量的Hilbert邊際譜，圖8為無損傷狀態下速度為10 km/h工況下測點1的模擬信號各IMF分量的Hilbert邊際譜.

圖8 無損傷狀態各 IMF 分量的 Hilbert邊際譜Fig.8 Hilbert marginal spectrum of IMF components in undamaged state

從圖8可看出：1）原始信號中有3個顯著分量，分別為imf2、imf3和imf4；2）imf2 的峰值頻率為9.063 Hz，與有限元計算的二階自振頻率8.910 Hz接近；imf3的峰值頻率為5.563 Hz，與有限元計算的結構一階自振頻率5.815 Hz接近；imf4對應的峰值頻率為2.813 Hz，與車輛的自振頻率3.200 Hz接近；3）圖中imf1分量雖然不顯著，但峰值頻率與有限元計算的第三階豎彎頻率10.303 Hz接近.

結合理論知識，可判斷出：imf2分量對應于自由振動部分的二階豎彎振型振動；imf3分量對應自由振動部分一階豎彎振型振動；imf4對應純強迫振動部分；imf1對應自由振動部分三階豎彎振型振動.

3.2 有損傷橋梁的模擬信號

速度為10 km/h工況下，有損傷狀態下測點1的模擬信號的各IMF分量的Hilbert邊際譜如圖9所示.

圖9 有損傷狀態各 IMF 分量的 Hilbert邊際譜Fig.9 Hilbert marginal spectrum of IMF components with damage

由圖9得到與圖8類似的結論：imf2、imf3、imf4是原始信號中的顯著分量，分別對應結構二階振型振動、一階振型振動和純強迫振動部分；imf1對應三階豎彎振型振動，但不顯著.從圖9中還發現顯著分量中，與二階豎彎振型振動和一階豎彎振型振動對應的imf2和imf3分量的邊際譜圖均出現了明顯的雙峰現象，且二階振型振動分量imf2的邊際譜明顯變寬，可利用此特征進一步研究結構損傷識別、定位和定量的方法.

3.3 橋梁的實測信號

V=10 km/h 時，測點 1 實測信號的各 IMF 分量的 Hilbert邊際譜見圖10（a）.可以看出，信號分解出現了模態混疊，影響了主成分的識別效果.這里選取該測點30 km/h時實測加速度信號來進行主成分分量識別，其各IMF分量的Hilbert邊際譜見圖10（b）.

由圖10（b）可知：1）原始信號有 imf1、imf2 和imf3 3個主要分量，分別對應于結構二階豎彎振型振動、一階豎彎振型振動和純強迫振動部分；2）實測信號中沒有發現與三階豎彎振型對應的IMF分量.

由上可見：實測信號和模擬信號總體識別效果是相當的，主要差別體現在實測信號中與強迫振動對應的分量imf3不是很顯著，沒有模擬信號那么好識別；實測信號中與三階豎彎振型振動對應的IMF分量沒有出現，這說明模擬信號中不顯著的振型振動分量試驗時可能量測不到.

與模擬信號的Hilbert邊際譜相比，圖10（b）中實測信號的IMF分量的Hilbert邊際譜圖的平滑性較差，主要是由于該實測信號采樣頻率取值不高導致的，試驗時采樣頻率取值僅為50 Hz，模擬信號的采樣頻率為100 Hz.

圖10 各 IMF 分量的 Hilbert邊際譜Fig.10 Hilbert marginal spectrum of IMF components

4 一階振型振動分量的Hilbert瞬時頻率譜特征分析

4.1 損傷橋梁模擬信號

識別出與振型振動分量對應的IMF分量后，選擇與一階豎彎振型振動對應的IMF分量進行Hilbert變換得到其 Hilbert瞬時頻率譜.圖11是V=10 km/h時的Hilbert瞬時頻率譜（藍色曲線），為了方便分析，將測點加速度原始信號波形也在圖中示出（紅色曲線）.圖11中瞬時頻率上下波動的均值與該分量邊際譜中的峰值處的頻率基本是對應的.可以看出：1）各測點一階振型振動分量的瞬時頻率在差不多相同時刻出現了明顯的頻率降低現象；在降低區域瞬時頻率波動也相對穩定，即出現了局部移頻；2）各跨中測點均表現出了相同的特征，說明基于瞬時頻率的損傷識別，不受測點位置影響.

圖11 模擬信號及其一階振型振動分量的Hilbert瞬時頻率譜Fig.11 Simulated acceleration signal and Hilbert instantaneous frequency spectrum of first-order mode vibration components

圖11能很好地說明瞬時頻率能發現結構的局部剛度變化，而通過傅里葉變換得到的頻率只能反應結構的整體剛度情況.因此可以利用瞬時頻率的此特征進行結構損傷定位.

4.2 橋梁實測信號

圖12 是V=10 km/h 時的實測信號通過分析，得到的各跨中測點與一階振型振動對應的IMF分量的Hilbert瞬時頻率譜（藍色曲線），測點加速度原始信號的波形（紅色曲線）.

圖12 實測信號及其一階振型振動分量的Hilbert瞬時頻率譜Fig.12 Measured acceleration signal and Hilbert instantaneous frequency spectrum of first-order mode vibration components

由圖12可看出：各測點均在差不多相同時間附近出現了局部瞬時頻率降低現象，根據模擬信號瞬時頻率曲線呈現的損傷特征不難類推出這種規律性的頻率降低是由損傷引起的，并通過頻率首次出現降低的時刻可以推斷損傷的位置，從而完成損傷的定位.

4.3 討 論

圖13是取實際橋梁脈動試驗的試驗數據為對象得到的各測點一階振型分量的Hilbert瞬時頻率譜.

圖13 環境激勵下一階振型分量的Hilbert瞬時頻率譜Fig.13 Hilbert instantaneous frequency spectrum of firstorder mode component under environmental excitation

從圖13可看出：環境激勵下，各測點一階瞬時頻率曲線沒有在相同的位置出現規律性的瞬時頻率局部降低現象，這表示損傷特征是需要移動車輛激發才能表現出來.所以，采用本文的方法對連續梁橋進行損傷識別時移動荷載作用是必要條件.

原則上也可以選擇二階或更高階豎彎振型分量的瞬時頻率作為損傷識別的指標.但從前面的分析可看出：模擬信號和實測信號雖然在識別振型分量上總體相當，但實測信號由于受環境的干擾較大，其振型分量的識別條件比模擬信號嚴格（譬如，第一跨跨中測點1的實測信號在速度為10 km/h工況下通過Hilbert邊際譜只能發現一階振型分量頻率，根本識別不出高階振型）；另外在模擬信號中能識別出來的第三階豎彎振型分量，實測信號中根本就沒發現；而且二階以上振型分量的識別還會受到測點位置的限制，譬如第二跨跨中測點正好位于二階振型的零點，利用此測點的數據就很難識別出二階振型振動分量.雖然，高階振型分量在識別時會有諸多影響，但無論是模擬信號還是實測信號，一階振型的識別還是比較穩定的，基本所有測點都能識別.

5 結 論

1）通過Hilbert邊際譜能較好地識別出與低階豎彎振型對應的IMF分量，這種識別是信號經過ICEEMDAN分解后的各IMF分量之間的Hilbert邊際譜幅值顯著性比較來實現的，而不需要進行模擬信號和實測信號的之間對比，避免了模型誤差.

2）模擬信號中通過Hilbert邊際譜能夠識別出來的振型振動分量比實測信號多，其中模擬信號中表現不顯著的高階豎彎振型分量，在實測信號中沒有出現.

3）與一階振型振動對應的IMF分量的瞬時頻率對局部損傷敏感，且其Hilbert瞬時頻譜橫坐標為時間，很容易轉換為位置函數，可作為特征參數用于橋梁結構的損傷識別和損傷定位.

4）以移動荷載作用下橋梁的動響應信號為對象，采用ICEEMDAN方法，通過Hilbert邊際譜和Hilbert瞬時頻率譜進行損傷識別的基本過程為：采用ICEEMDAN對原始信號進行分解，得到一系列IMF分量；通過各IMF分量的Hilbert邊際譜幅值的比較識別與一階振型振動對應的IMF分量；以一階振型振動分量為對象，通過其Hilbert邊際譜和瞬時頻率譜進行損傷識別、定位和定量.

5）依賴移動荷載在橋梁上的遍體歷經性而完成的損傷識別結果受測點位置的影響較小，說明在實際應用中利用移動荷載進行損傷識別，不一定要將拾振器放在損傷點上或附近才能進行損傷識別.