基于環境激勵的大跨度斜拉橋模態參數和索力識別

葉錫鈞，顏全勝，李 健，王衛鋒，朱添豐，劉明慧

(1.華南理工大學 土木與交通學院，廣州 510641;2.伊利諾伊大學 香檳分校工學院，美國香檳 61801)

基于環境激勵的大跨度斜拉橋模態參數和索力識別

葉錫鈞1，顏全勝1，李 健2，王衛鋒1，朱添豐1，劉明慧1

(1.華南理工大學 土木與交通學院，廣州 510641;2.伊利諾伊大學 香檳分校工學院，美國香檳 61801)

模態參數和索力是評估斜拉橋健康狀態的關鍵參數。建立了崖門大橋的有限元模型并對其進行了基于環境激勵的模態測試和拉索振動測試;提出了基于ERA的多參考點穩定圖算法，設置不同的參考點，利用自然激勵技術結合特征系統實現算法(NExT－ERA)識別模態參數，通過阻尼比、基于輸出矩陣的一致模態指標(CMI_O)和模態置信度(MAC)作為判別標準，識別出崖門大橋的豎向和橫向模態參數，通過和增強頻域分解法(EFDD)識別結果的比較，可知該算法具有良好的識別效果;分析了斜拉索與主梁的共振頻率范圍，通過二次擬合識別較長拉索的低階頻率，根據兩種不同方法的索力識別結果可知，該橋的索力分布比較均勻對稱。

環境激勵;斜拉橋;模態參數識別;多參考點穩定圖;索力識別

對于運營中的大跨度斜拉橋，基于環境激勵的模態測試和頻率法［1］分別是測定其動力特性和索力的最快捷方法，而且幾乎都是唯一選擇。

目前，模態參數識別的方法主要有峰值法(PP)、頻域分解法(FDD)、增強頻域分解法(EFDD)等頻域方法和特征系統實現算法(ERA)、隨機子空間算法(SSI)、時間序列分析法等時域方法。頻域的識別方法簡單、快速而且直觀，但在識別阻尼比時需進行反傅里葉變換，在時域內通過指數衰減法來實現，受截斷誤差的影響，精度不高;而時域的識別方法，其共同的缺點是難以確定結構的系統階次，難于區分真實模態與噪聲模態［2－3］。近年來，穩定圖方法［4－5］的應用越來越廣，被認為是一種比較理想的確定系統階次的方法，可以從諸多模態中判別真假模態，但如何設置判別標準是關鍵步驟。

對于索力識別，近年來國內外學者研究得很多。Zui等［6］提出了利用低階頻率估計索力的實用公式;王衛鋒［7］以二次拋物線代替懸鏈線為基礎，提出了考慮拉索垂度、抗彎剛度、邊界條件和斜度的影響的拉索振動方程，通過迭代求解方程得到索力值。對于各種不同的方法，準確識別拉索的自振頻率是關鍵。

本文建立了崖門大橋的有限元模型并對其進行了基于環境激勵的模態測試和拉索振動測試;提出基于ERA的多參考點的穩定圖算法，識別崖門大橋的模態參數;并利用二次擬合求出長索的低階頻率，根據Zui等［6］的經驗公式和王衛鋒［7］的方法識別索力，進行比較。

1 NExT－ERA

1.1 NExT

NExT(Natural Excitation Technique)方法是由James等［8］于1993年提出，其基本思想是利用白噪聲環境激勵下結構兩點之間響應的互相關函數和脈沖響應函數具有相似表達式的特點，將其代替結構的自由振動響應或脈沖響應函數作為輸入，運用時域模態識別方法進行模態參數識別。NExT法的計算有兩種方法，一種是在時域里采用卷積算法直接計算得到，另一種是先計算出實測信號的互功率譜，再經傅里葉逆變換得到實測互相關函數［9］。本文采用后者，由于其采用了統計平均處理，即使會受到一定的譜泄漏的影響，但所得到的互相關函數的信噪比也有較大幅度的提高。

1.2 ERA

ERA(Eigensystem Realization Algorithm)是多輸入多輸出的時域模態參數識別方法，由 Juang等［10］于1984年提出，其利用實測的脈沖響應函數為基本模型構造廣義Hankel矩陣，運用奇異值分解方法，確定系統的階次和狀態方程中的系統矩陣A、輸入矩陣B和輸出矩陣C，進而求解系統矩陣A的特征值問題，求得極點與留數，從而確定系統的模態參數，具體的算法公式詳見文獻［10］。

2 基于輸出矩陣的一致模態指標

Pappa［11］提出了一致模態指標 CMI，其計算根據EMAC和MPC，公式如下：

其中：i表示第i階振型;j表示第j個輸出點，k表示第k個輸入點;EMACjik表示第i階振型下，對應j－k(輸出－輸入對)的EMAC值;和分別是基于輸出矩陣C和輸入矩陣B的EMAC值。

Yun［12］對CMI進行改進，提出基于輸出矩陣C的模態一致性指標(CMI_O)，并通過數值算例驗證了在基于環境激勵的模態分析過程中，CMI_O比CMI更具可靠性：

式中：i表示第i階振型;j表示第j個輸出點;q表示輸出點的數量。從式(3)可以看出，CMI_O的計算基于，僅取決于輸出響應。

3 基于ERA的多參考點穩定圖算法

對于環境激勵的情況，由于系統的輸入激勵未知，通過NExT/ERA方法識別模態參數時，選擇某輸出點作為參考點并計算其與其他輸出點的互相關函數矩陣作為輸入，參考點非常重要，選擇不當則不能識別出某些較弱的模態。基于Pappa等［13］應用的真實模態判別方法，本文提出一種基于ERA的多參考點穩定圖算法，以基于輸出矩陣的CMI_O代替CMI作為判別準則，本文選取CMI_O閥值為75%。主要包括圖1的三個過程：

圖1 多參考點穩定圖算法的三個過程Fig.1 3 parts of the multiple reference DOFs stabilization diagram algorithm

算法步驟如下：

(1)選定初始參考點r，設定不同的系統階次(如，n+10，n+20，…，max_n)，運用NExT/ERA識別出不同階次下的模態參數;每個階次下識別出的模態參數，若能通過閥值“0% ＜阻尼比＜10%和CMI_O＞75%”，則認為是可信度較高的模態參數，進入下一步判別，共有N個模態參數;

(2)通過設置Δf，將第(1)步中的N個模態歸結為M組(每個組里各個模態相互間的頻率差值小于Δf，可以認為是同一階的頻率)，不符合該要求的模態將被記錄下來加到當參考點為r+1時重復這一步判別;

(3)判斷M組模態中，每一組里模態的數目T，當T＞10時，認為該模態是真實模態(否則將被記錄下來加到當參考點為r+1時重復這一步判別);同時尋找每一組中CMI_O最大值對應的模態作為本階段結果，得到P階模態;

(4)通過步驟(3)可以在穩定圖上清晰看到穩定軸，計算P階模態中任意兩階之間的MAC值，小于70%的則認為都是真實模態;若大于70%，對應CMI_O值較大者則認為是真實模態;得到參考點為r時的模態參數組Set(r)，在穩定圖上以紅點標記;

(5)參考點位置為r+1，依步驟(1)至(4)得到Set(r+1)。比較Set(r)和Set(r+1)同一階模態的CMI_O，若 CMI_O(modei(r))＞ CMI_O(modei(r+1))，則把Set(r+1)中的 modei(r+1)用 Set(r)中的modei(r)代替(CMI_O值較大者認為更準確);若Set(r)中有Set(r+1)沒有出現的模態，則加到Set(r+1)中，得到新的參數組NewSet(r+1);

(6)改變參考點位置為r+2，重復步驟(1)至(4)得到Set(r+2)，再根據步驟(5)比較NewSet(r+1)和Set(r+2)得到NewSet(r+2)，直到最后一個測點作為參考點，得到最終的模態參數。

算法流程圖如圖2所示。

圖2 多參考點穩定圖算法流程圖Fig.2 Flowchart of the multiple reference DOFs stabilization diagram algorithm

4 實例分析

4.1 工程概況

崖門大橋是廣東省西部沿海高速公路建設項目的重點工程之一，大橋為雙塔單索面預應力混凝土斜拉橋，采用塔梁固結柔性墩體系。每側邊跨設置一個輔墩，跨徑組合為(50 m+115 m+338 m+115 m+50 m)，全橋總長668 m。主梁是單箱五室預應力混凝土梁，梁高3.48 m，橋面全寬26.8 m，設雙向四車道。主墩為雙壁柔性墩，橫向寬12 m，雙壁之間中心距6 m，墩高47.6 m。橋面以上塔柱高73.5 m，塔柱斷面為單箱混凝土斷面。全橋共有斜拉索200根，主梁標準索距6 m，如圖3(a)。

利用ANSYS軟件，按照塔梁固結形式，采用Beam4和Link8單元建立大橋的空間桿系有限元模型，斜拉索的初始應變值按照設計索力計算給定，除采用Ernst公式修正斜拉索的彈性模量外，其他參數按照設計參數給定，如圖4。

4.2 模態測試與參數識別

斜拉索的布置為單索面，主梁中間有寬3 m的分隔帶，由于交通限制，只在主梁中間分隔帶布置豎向和橫向測點，見圖3(b)。模態測試的儀器采用中國地震局工程力學研究所的99－1型拾振器、941型8線放大器和杭州憶恒科技的AVANT動態數據采集系統。

4.2.1 豎向模態參數識別

全橋共布置45個豎向測點，以跨中的22、23測點為參考點，分為11組，每組6個測點(第8組除外，表1所示)，設置采樣頻率為160 Hz，每組采樣時間20 mins。

表1 豎向模態測試方案Tab.1 Plan for vertical modal test

從有限元計算結果得知，斜拉橋的前幾階固有頻率都在10 Hz以內，故選取較為平穩的數據段，把采樣頻率降到40 Hz，對數據進行去除趨勢項和平滑去噪處理后采用增強頻域分解法(EFDD)和本文算法進行識別。限于篇幅，本節只列出豎向測試的第4組數據的模態參數識別過程。典型的加速度時程見圖5。

通過EFDD的識別結果可清晰看出，在0～3 Hz范圍內有8階主要的模態，其中第4和第8階模態較弱，其為ERA的時域識別方法提供很好的基礎，見圖6。

圖5 豎向第4組：基準測點22、23的加速度時程Fig.5 Acceleration responses of Channel 22＆23 in group 4

對于NExT－ERA，選擇的參考點不同，識別出的模態數目可能不一樣。以豎向測試的第4組為例，根據本文算法，選定參考點為1通道，通過第(1)～(3)步，由穩定軸觀察到7階模態(圖7，T＞10的共有7階模態)，再根據第(4)步比較7階模態里相互間的MAC值確定是否真實模態，得到模態參數Set(1)(6階，第7階與第6階的MAC值大于70%，排除)，在穩定軸上標記紅點。

同理，改變參考點為2通道，經過(1)～(4)步，得到模態參數Set(2)，第2和4階模態沒被識別出來(圖8);接下來第(5)步，除了比較兩組模態 Set(1)和Set(2)的CMI_O值外，還把Set(1)中的第2和4階模態加到Set(2)中，得到Newset(2)。當參考點為通道3時識別的參數組為Set(3)，比較Newset(2)和Set(3)，得到Newset(3)，直到參考點變為通道6得到第4組的最終模態參數。

同時，對于某些比較弱的模態，往往要通過設置較高的系統階次才能識別出來。如圖7，當通道1為參考點時，直到系統階次為240左右，第1、2階模態才開始出現;而參考點為通道2時，圖8，第1階模態在階次為10時就出現了，但第2階模始終沒被識別出來。

圖6 豎向第4組：EFDD識別結果Fig.6 Singular value by EFDD：group4，vertical

圖7 豎向第4組：參考點為通道1時的穩定圖Fig.7 Stabilization Diagram of vertical group 4：Ref.=Channel 1#

圖8 豎向第4組：參考點為通道2時的穩定圖Fig.8 Stabilization Diagram of vertical group 4：Ref.=Channel 2#

觀察第1～11組的模態參數，在每組里面都出現的模態認為是真實模態。對于EFDD，求11組模態參數的平均值為最終結果;對于本文算法，對比各組里同一階模態的CMI_O值，最大值為最終識別結果，結果列于表2。

從表2看出，本文算法與EFDD識別的結果基本符合，但識別結果與有限元計算結果存在一些差別。

表2 豎向模態參數識別結果Tab.2 Identified vertical modal parameters

按最大位移值歸一化方法，得到橋面豎向振型。由于測試的分組數目較多，識別的振型誤差會有累計，但可以看出，本文算法識別的振型更接近有限元結果。限于篇幅，本文只列出前三階豎向振型，見圖9～圖11。

4.2.2 橫向模態參數識別

橫向測試共布置22個測點，只在圖3(b)中邊跨偶數和中跨基數測點，以19和21測點為參考點，分為7組，每組5個測點(第7組4個測點除外)。

通過EFDD的奇異值分布圖，在0～0.7 Hz范圍能識別出4 階模態(0.332，0.469，0.566，0.664)，圖 12。但通過本文算法最終確定3階真實模態，圖13和表3。識別的前3階橫向振型見圖14～圖16。

表3 橫向模態參數識別結果Tab.3 Identified transverse modal parameters

圖14 第1階橫向振型Fig.14 1st Transverse bending mode

圖15 第2階橫向振型Fig.15 2nd Transverse bending mode

圖16 第3階橫向振型Fig.16 3rd Transverse bending mode

4.3 斜拉索振動測試與索力識別

索力測試采用Imote2無線傳感器節點和SHM_A傳感器板［14］，將其固定在離減震器2 m的位置，測量斜拉索在豎向平面的振動，圖17。設置采樣頻率100 Hz，濾波頻率30 Hz，采樣時間為3 mins，對東塔和西塔中跨各5對不同長度的斜拉索(1、7、13、19、25，每對上、下游各1根，共20根)進行振動信號采集(見圖3(a))。分析采用的降采樣頻率為50 Hz，NFFT=1 024，頻率分辨率為0.005 Hz。

4.3.1 頻率識別

在測試過程中發現，在無風或風荷載很小的時候，有些斜拉索也會產生大振幅的劇烈振動，同時伴有橋面的振動。將各索的自功率譜密度函數依次縮小一定的倍數和主梁的豎向振動頻譜圖畫在同一張圖上(圖19)，可以發現，較長的斜拉索(13、19、25)在0 ～3 Hz范圍內的頻率與橋面的振動頻率相同;可以認定上述拉索的大振幅振動是由橋面振動激發引起的受迫振動，其實測的頻率是受迫振動的頻率，而不是其自振頻率，從頻譜圖上難以觀察出其基頻的位置，只能識別出其第3或4階以后的頻率(圖中箭頭);而在7號索以前的斜拉索，由于其基頻較大(大于2 Hz)，故能從頻譜圖上清晰地識別出來。圖19中，能通過峰值識別出1號索的自振基頻為4.023 Hz，7號索的自振基頻為2.031 Hz。

圖19 東塔上游斜拉索索梁耦合振動的頻譜圖Fig.19 Power spectrum of Deck and Cables：East pylon

對于不能識別前幾階頻率的較長的13、19、25號斜拉索，由于使用頻率間距法確定拉索的低階頻率誤差較大，本文采用二次擬合的方法計算其前幾階的頻率。

圖20(a)是東塔上游25號索的各階頻率擬合結果，從頻譜圖識別第4到20階頻率，以二次曲線進行擬合。實測的和擬合的各階頻率差值見圖20(b)。可見，實測各階頻率的差值最小為0.703 Hz，最大為0.79 Hz，且各階頻率差值的分布無序;而擬合的各階頻率的差值變化是平緩的，從0.731 Hz到0.767 Hz，與文獻［15］的結論較為吻合。

圖20 東塔上游25號索的頻率識別結果Fig.20 Identification result of Cable 25#(upstream)of East Pylon

4.3.2 索力識別

由于安裝減震器后拉索的長度變短，頻率變大，如按原索長計算索力必然造成很大的誤差。故認為2002年竣工驗收［16］時安裝減震器前后的恒載索力不變，按照安裝減震器前后的頻率對計算索長進行修正。利用MATLAB根據Zui等［6］的實用公式編制索力計算程序，以100 m為標準將索初步判定為長索或短索;再以Γ和ξ作為依據進行迭代判斷，該使用7個實用公式中的哪一個;長索用第4階頻率為計算，短索則用第1階或第2階頻率為計算。王衛鋒［7］的索力計算方法采用跟Zui的方法同樣的頻率階數和頻率值。這里應當指出，對于長索且滿足Γ和ξ條件的，Zui的公式只用到第4階頻率，故上一步中擬合的低階頻率誤差對該索力的識別沒有影響。

圖21 東、西塔實測索力結果Fig.21 Result of cable tension estimation of East and West Pylon

從圖21中看出，索力的整體分布比較對稱均勻，兩種方法計算的索力值很接近，誤差在5%之內。除了西塔1號索外，其余斜拉索的實測索力都比設計索力要大，這主要是由于施工過程中為保證主梁線型而作了索力調整，而且經過對比修正后的計算索長也不能保證十分精確。

5 結論

從模態識別結果看出，對于環境激勵的模態分析，使用CMI_O代替CMI作為閥值，通過設置不同的參考點計算互相關函數作為輸入，基于ERA的多參考點穩定圖算法能很好地識別出大橋的低階模態參數，有效地分辨虛假模態和避免真實模態的遺漏，準確地識別出主梁豎向和橫向的前3階振型。

但在分析過程中發現，有三個地方的取值不能是一個定值：第(3)步中判斷模態出現的次數T時(本文取值T=5)，對于某些模態測試環境理想、激勵充分的情況下T值能取大些;另在判別準則的取值上還是取決于經驗，如第(1)步中設置閥值時CMI_O的取75%;第(4)步中MAC的取70%，都是根據經驗取值的，而且在一定程度上取決于模態試驗的效果。在這些問題上需要工作者能進行一些經驗判斷。

通過頻率識別和索力識別結果，可知主梁和索在0～3 Hz范圍內有共振現象;二次擬合法比頻率間距法能更準確地識別出較長拉索的低階頻率;崖門大橋的中跨斜拉索索力分布比較對稱均勻;但在如何精確地修正安裝減震器后的索長問題上還得進一步研究。

［1］王衛鋒，韓大建.斜拉橋的索力測試及其參數識別［J］.華南理工大學學報(自然科學版)，2001，29(1)：18－21.

［2］續秀忠，華宏星，陳兆能.基于環境激勵的模態參數辨識方法綜述［J］.振動與沖擊，2002，21(3)：2－5.

［3］王 彤，張令彌.運行模態分析的頻域空間域分解法及其應用［J］. 航空學報，2006，27(1)：62－66.

［4］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.

［5］常 軍，孫利民，張啟偉.基于兩階段穩定圖的隨機子空間識別結構模態參數［J］.地震工程與工程振動，2008，28(3)：47－51.

［6］ Zui H，Shinke T，Namita Y.Practical formulas for estimation of cable tension by vibration method［J］.Journal of Structural Engineering，ASCE，1996，122(6)：651 －656.

［7］王衛鋒.PC斜拉橋的施工監測［D］.廣州：華南理工大學建筑學院，2000.

［8］James G H，Carne T G，Lauffer J P.The natural excitation technique(NExT)for modal parameter extraction from operating wind turbines［R］.Albuquerque：Sandia National Laboratory report，SAND92 －1666.UC －261，1993.

［9］王 濟，胡 曉.MATLAB在振動信號處理中的應用［M］.北京：中國水利水電出版社，2006.

［10］ Juang J N，Pappa R S.An eigensystem realization algorithm for modal parameter identification and model reduction［J］.Journal of Guidance，Control and Dynamics，1985，8(5)：620－627.

［11］ Pappa R S，Elliott K B.Consistent－mode indicator for the eigensystem realization algorithm［J］.Journal of Guidance Control and Dynamics ，1993，16(5)：852 －858.

［12］ Yun G J.Modal identification and damage detection for structural health monitoring under ambient vibration environment［C］//Structures Congress 2009.Austin，USA：American Society of Civil Engineers，2009.

［13］ Pappa R S，James G H，Zimmerman D C.Autonomous modal identification of the space shuttle tail rudder［J］.Journal of Spacecraft and Rockets ，1998，35(2)：163 －169.

［14］ Illinois SHM Project：http：//shm.cs.uiuc.edu.

［15］侯立群，歐進萍，趙雪峰，等.哈爾濱四方臺斜拉橋模態參數和索力識別［J］. 振動與沖擊，2009，28(5)：106－110.

［16］胡利平.崖門大橋主橋斜拉橋－工程檢測試驗報告［R］.廣東省交通建設工程質量檢測中心，2002.

Modal identification and cable tension estimation of
long span cable-stayed bridge based on ambient excitation

YE Xi-jun1，YAN Quan-sheng1，LI Jian2，WANG Wei-feng1，ZHU Tian-feng1，LIU Ming-hui1
(1.School of Engineering and Transportation，South China University of Technology，Guangzhou 510641，China;2.College of Environmental Engineering，University of Illinois at Urbana－Champaign，Champaign 61801，USA)

For health monitoring of long span cable－stayed bridge，modal parameters and cable tension are the key parameters to assess the condition of the bridge.A finite element model of Yamen bridge was built.Modal test and cable vibration test of the bridge were performed under ambient excitation.An improved multiple reference DOFs stabilization diagram algorithm based on ERA(Eignsystem Realization Algorithm)was presented.By setting different reference DOFs in each group of data，NExT(Natural Excitation Technique)－ERA was used to identify modal parameters.Damping ratio，consistent mode indicator from observability(CMI－O)and modal assurance criterion(MAC)were used as thresholds to identify the most accurate modal parameters.Lower order frequencies were estimated by quadratic fit method，and the cable tension was estimated by two different methods.Based on the analysis of deck and cable vibration，it is evident that the vertical vibration of the bridge deck is tightly coupled with the cable vibrations within the frequency range of 0～3Hz.

ambient excitation;cable－stayed bridge;modal parameter identification;multiple reference DOFs stabilization diagram;cable tension estimation

TU311;O327

A

廣東省交通運輸廳科技項目(2010－02－015)

2011－03－30 修改稿收到日期：2011－09－09

葉錫鈞 男，博士生，1984年10月生