實測模態下三跨整體橋精細化基準有限元模型

張志勇，何樂洽，薛俊青，陳夢圓，黃福云，Giuseppe Carlo Marano，Bruno Briseghella

(福州大學土木工程學院，福建 福州 350108)

0 引言

基于物理模型的損傷識別，是土木工程結構健康監測的一項重要研究內容. 通過比較結構無損狀態和損傷狀態下動力特性的差異，可發現定位和定量結構的損傷程度，從而為工程結構的智能化管養提供有力的決策支持[1]. 有限元技術是現代橋梁與建筑結構設計和分析的重要支撐. 通過建立合理的設計模型，采用必要的設計假定和簡化，能夠解決結構設計階段所面臨的主要問題[2]. 大量的結構長期健康監測結果表明，工作狀態下的梁橋，其動力特性(如自振頻率、阻尼等)受到各種因素的顯著影響[3]. 前人的研究工作揭示，根據結構無損/初始狀態下實測模態信息建立的有限元基準模型可為結構的長期健康監測提供精確、可靠的參考信息[4]. 作為能夠全面、正確反映結構真實行為的基準有限元模型，除了應保證與實測結果較好的吻合度外，還應保證模型的各項參數真實準確的物理意義.

整體式橋臺橋梁(簡稱整體橋)屬于無伸縮縫橋梁，常見于中小跨徑橋梁[5-6]. 整體橋將傳統橋梁的上下部結構連成一個整體，為高階超靜定結構. 本文以一座三跨整體橋為工程背景，考慮結構-土相互作用，基于環境振動試驗實測模態信息，開展了該橋基于梁格法設計簡化模型與梁-殼單元精細化模型的模型修正. 將與試驗頻率和振型的吻合度為目標，通過調整模型的物理參數，建立該橋初始狀態的基準動力有限元模型[7]. 最后為整體橋動力分析提出了建議.

1 工程背景

馬巒山橋位于深圳市馬巒山郊野公園北側，于2019年年底建成通車(見圖1). 該橋分南、北兩幅，雙向八車道. 由于兩幅橋的結構完全一致，故僅以南幅橋作為本文的研究對象. 該橋是一座三跨預應力混凝土箱梁整體橋，采用先成橋后填筑臺后土的方式建造. 橋臺后填土采用中等密實砂性土. 在樁基深度方向，樁周土分別為素填土、粉質粘土、粘土、殘積土與不同風化程度的粉砂巖. 橋梁跨度為87.92 m(28.96 m+30 m+28.96 m)，橋寬為17 m，主梁采用單箱三室等高箱梁. 全橋取消了橋臺和橋墩處的伸縮縫、伸縮裝置和支座，采用整體式橋臺將主梁梁端與橋臺固結在一起，墩梁也采用固結. 墩下和橋臺下均采用單排混凝土樁基礎. 由于溫度變形產生的梁體伸縮量引起橋臺和橋墩的縱橋向變形，并同時通過柔性樁和橋頭搭板引起臺后填土和樁周土的相應變形，從而產生結構和土的相互作用.

圖1 深圳馬巒山整體橋Fig.1 Maluan mountain integral abutment bridge in Shenzhen

2 環境激勵下的模態試驗

圖2 試驗測點布置圖(單位： cm)Fig.2 Layout of test points (unit: cm)

為了獲取該橋實際的動力特性，在風和少量車輛為主的環境激勵下進行不同測點布置的兩次模態試驗. 試驗測試儀器為4臺瑞士GEOSIG公司生產的三向無線加速度傳感器GMSPLUS，可同時記錄縱橋向、橫橋向和豎向的加速度響應. 試驗一與試驗二中，分別在橋面上等間距的布置了13個與39個測點. 測點位于道路兩側和主梁橫斷面中間，測點的布置如圖2. 采用隨機子空間算法，基于實測加速度信號識別到了該橋前五階試驗模態，包括一階豎彎、二階豎彎、三階豎彎、一階扭轉和二階扭轉. 試驗結果詳見表1，其中f為自振頻率，ξ為阻尼比.

表1 實測模態與修正后計算模態對比

3 有限元模型的建立

3.1 基于梁格法的有限元模型

3.1.1上部結構模擬

圖3 箱梁截面圖與梁格法劃分(單位： cm) Fig.3 Cross section of box girder and division of grillage method (unit: cm)

圖3為馬巒山整體橋的箱型主梁橫截面. 該橋主梁寬17 m，高1.6 m，縱向預應力鋼束布置在腹板內部. 主梁基于梁格法設計，采用MIDAS-civil軟件建立全橋三維有限元模型. 圖3中虛線所示為橫斷面梁格劃分. 虛擬橫梁沿縱橋向的間距為1 m. 虛擬橫梁截面僅計入主梁的頂板和底板對于剛度的貢獻，不計入其質量. 橋面混凝土調平層與瀝青鋪裝層以及防撞護欄僅作為附加質量考慮，忽略其剛度影響. 模型中采用的是基于Timoshenko梁理論的三維梁單元，采用分布質量且不考慮鋼筋的影響.

3.1.2下部結構模擬

橋臺采用梁格法模擬，主梁與橋臺采用共節點的方式連接，橋墩、承臺和樁基均用梁單元模擬，相互之間采用主從節點的剛性連接.

模型中考慮了橋臺-臺后土與樁-樁周土的相互作用，分別采用非線性與線性拉壓彈簧模擬橋臺-土與樁-土的相互作用. 其中橋臺處僅在縱橋向設置了土彈簧，而樁周在縱、橫橋兩個方向上均設置了相同的土彈簧，如圖4所示. 下面簡要介紹橋臺土彈簧與樁周土彈簧剛度的計算方法[5, 8].

該橋臺后土為中等密實砂性土，橋臺與臺后土僅考慮縱橋向的相互作用，忽略橋臺與土之間的豎向和水平向摩擦力. 模型中，將臺后土彈簧作用在橋臺梁格節點上，按照該節點與相鄰節點的間距和土層的性質來確定其剛度. 例如對于橋臺上某節點，當橋臺位移為零時受到靜止土壓力Po=Koγzbt的作用.模型中，將其作為附加力單獨施加到橋臺上； 當橋臺擠壓填土時，土壓力變化至被動土壓力Pp=Kpγzbt，此時橋臺土彈簧受壓； 當橋臺遠離填土時，土壓力變化至主動土壓力Pa=Kaγzbt，此時橋臺土彈簧受拉.其中Ko、Kp、Ka分別為靜止、被動、主動土壓力系數，γ為土的容重，z為該節點處臺后填土深度，b為該節點處左右相鄰單元長度和的1/2，t為該節點處上下相鄰單元長度和的1/2.由此，該節點處土彈簧受壓時的剛度為k1=(Pp-Po)/Δ1； 受拉時的剛度為k2=(Pa-Po)/Δ2.其中，Δ1與Δ2為橋臺節點位移. 上述各參數的取值均可參考NCHRP規范曲線的相應數值[9]. 相應地，橋臺梁格模型上各節點的土壓力與橋臺位移關系可用圖5來表示.

圖4 馬巒山整體橋結構模型南立面圖(半幅，單位： cm)Fig.4 Structural model of Maluan mountain bridge: south elevation view(half span, unit: cm)

圖5 橋臺各節點土壓力與橋臺位移關系(H為橋臺高度)Fig.5 Relationship of earth pressure at each node of abutment and its displacement(H: abutment height)

整體橋中樁-土的相互作用不可忽略，通常可按照溫克爾假定，把樁視作彈性地基上的梁，將樁側土離散為線性彈簧. 該橋樁周土包括素填土、粉質粘土、粘土、殘積土與不同風化程度的粉砂巖(從上到下). 模型中，樁基梁單元按土層厚度進行取整劃分. 在樁基節點處分別設置橫橋向與縱橋向的土彈簧. 對樁基某節點處的土彈簧剛度，取為Kz=Cz·az·bz.其中，地基系數Cz=m·z，m為地基土的比例系數，z為樁基節點所在的地基土深度，az為該節點處土層的厚度，bz為該節點處樁基計算寬度. 上述各參數的取值均可參考《公路橋涵地基與基礎設計規范(JTG 3363—2019)》[10]的附錄P. 樁周土彈簧施加在樁基單元的節點處. 考慮到該樁為摩擦樁，且環境激勵條件下豎向位移微小，與結構位移比較可以忽略不計. 故將樁基單元節點豎向固結.

3.2 基于梁-殼單元的精細化有限元模型

3.2.1上部結構模擬

利用通用有限元軟件ANSYS建立了全橋精細化的三維梁-殼單元模型. 橋梁上部結構采用了殼單元SHELL181，將箱型主梁橫截面分上頂板、中腹板、下底板三個部分進行模擬. 單元厚度分別取為相應截面厚度的平均值. 主梁在與橋臺、橋墩連接處設置有橫隔，橫隔單元厚度按實際厚度取值.

將橋面混凝土調平層與瀝青鋪裝層作為附加質量加到主梁頂板上，忽略其剛度貢獻. 采用梁單元BEAM188模擬防撞護欄，考慮到護欄結構的非連續性和不確定性，根據參考文獻[11]將其初始彈性模量做一定量的折減. 與梁格模型中一樣，所有單元均采用了分布質量且不考慮鋼筋的影響.

3.2.2下部結構模擬

橋梁下部結構中，橋臺采用殼單元SHELL181模擬，橋墩、承臺和樁基均采用梁單元BEAM188模擬. 主梁與橋臺采用共節點的方式連接. 主梁與橋墩基于約束方程建立連接. 橋臺與承臺、橋墩與承臺、承臺與樁基均采用主從節點的剛性連接.

與梁格模型一樣，梁-殼單元模型基于土彈簧考慮了橋臺-臺后土與樁-樁周土的相互作用. 用離散的非線性彈簧單元COMBIN39，模擬臺后填土對橋臺的作用. 同前，彈簧的剛度按NCHRP所給出的土壓力系數曲線推算確定，具體計算方法和梁格模型相同. 同時，用離散的線性彈簧單元COMBIN14模擬樁周土對樁的作用，其剛度計算方法同前. 與梁格模型相同，將樁基單元節點采用豎向固結.

最終，簡化的梁格模型與精細化的梁-殼單元模型分別見圖6與圖7. 其中，梁格模型總計有1 356個節點，1 918個單元，而梁-殼單元模型總計有6 662個節點，6 975個單元. 兩類模型方位與圖1實橋一致. 兩側橋臺不等高，左側為6.97 m，右側為6.46 m，且臺下與墩下的樁基長度也不等. 梁格模型與梁-殼單元模型各部分的材料特性取值參考了《公路鋼筋混凝土及預應力混凝土橋涵設計規范(JTG 3362—2018)》[12].

圖6 馬巒山整體橋有限元梁格模型Fig.6 Finite element grillage model of Maluanshan bridge

圖7 精細化的ANSYS梁-殼單元模型Fig.7 Refined ANSYS beam-shell element model

4有限元模型的參數化修正

4.1 目標函數的確定

采用基于靈敏度有限元模型修正，將有限元模型的計算模態與結構實測模態的頻率和振型的殘差平方和作為目標函數. 通過最小化該目標函數，可以獲得較為精確的動力有限元模型. 該過程可視作一個數值優化求解問題，即如下所示非線性最小二乘問題：

(1)

式中：J(θ)為目標函數，其中θ為模型參數；Wε為權重因子對角矩陣用來平衡目標函數中各項殘差的貢獻.出于簡化考慮，此處Wε取為單位矩陣；εz為模態數據z的計算值與實測值之間的殘差.對于無阻尼特征值z=λ，此時殘差ελ, i(θ)可表達為

(2)

式中：λi=ωi2=(2πfi)2，i∈{1, 2, …,nλ}.其中上波浪號表示實測值，f或ω表示結構固有頻率或圓頻率.

對于特征向量(振型)，采用模態置信準則(modal assurance criterion, MAC)計算，公式為：

其中:φa表示計算振型向量;φe表示實測振型向量.此時z=MAC，殘差εMAC, i(θ)可表達為

εMAC, i(θ)=1-MAC

(3)

式中：i∈{1, 2, …,nMAC}.

對于梁格模型采用僅考慮頻率的目標函數：

(4)

對于精細化的梁-殼單元模型，采用聯合頻率與振型的目標函數，最終的目標函數如下：

(5)

4.2 參數靈敏度分析與選取

為了選取對計算模態影響較大的模型參數，從而確定優化問題的設計變量，需要進行模型的參數靈敏度分析. 針對計算頻率ω參數靈敏度矩陣定義為式(6)，則針對計算振型與實測振型MAC的參數靈敏度矩陣同理可得.

(6)

由無阻尼廣義的特征值方程Kφi=λiMφi可知，同時修正計算模型的質量與剛度將無法得到唯一解. 故靈敏度分析時，不考慮模型參數的質量密度，僅分析模型參數的彈性模量. 針對表1中前五階計算頻率與振型MAC值，并考慮到一階縱向剛體模態對整體橋的特殊性(結構-土相互作用)，分別對梁格模型與精細化梁-殼單元模型進行參數靈敏度分析.

分析結果發現： 在梁格模型中，自振頻率對主梁與橋臺梁格的縱梁與橫梁的剛度均較為敏感，而對樁周與臺后土彈簧剛度不敏感. 據此，考慮將主梁縱梁彈性模量、主梁與橋臺虛擬橫梁彈性模量、橋臺縱梁彈性模量3個參數作為梁格模型待修正參數.

在梁-殼單元模型中，針對主梁彈性模量、橋墩彈性模量、橋臺彈性模量、橫隔梁彈性模量、護欄彈性模量、樁周土彈簧剛度、臺后土彈簧剛度7個參數進行靈敏度分析，結果(取絕對值)見表2. 由表2可知，實測的前五階頻率對主梁剛度、橋墩剛度、橋臺剛度、橫隔梁剛度較為敏感，對樁周與臺后土彈簧剛度靈敏度較小，且均小于2%. 樁周與臺后土彈簧剛度僅對一階縱向剛體頻率敏感，約10%. 但由于在環境振動條件下未實測到該橋縱向剛體模態，不宜對其進行修正. 同時，由表1中各階實測振型圖結果獲知，橋臺底部振型各向位移和橋墩底部振型各向位移基本上均在0～0.03之間. 據此認為，在試驗條件環境激勵情況下，橋臺和橋墩承臺下的樁基基礎剛度對結構整體剛度影響較小. 由于本文關注重點在于整體橋橋面的動力特性，故不再將臺后土彈簧剛度與樁周土彈簧剛度作為修正參數. 且文中已將橋臺彈性模量和橋墩彈性模量列為模型修正參數，認為可以部分補償將樁基模型節點豎向固結造成的誤差. 綜合以上分析，同時考慮到非結構構件的影響，增加護欄剛度作為待修正參數. 最終確定以下5個參數作為梁-殼單元模型待修正參數，分別是： ① 主梁彈性模量； ② 橋墩彈性模量； ③ 橋臺彈性模量； ④ 主梁橫隔彈性模量； ⑤ 護欄彈性模量.

表2 梁-殼單元模型參數靈敏度矩陣

4.3 修正過程及結果4.3.1 梁格模型修正

參考式(4)定義的目標函數，采用手動調整4.2節中確定的3個模型參數的方式，進行梁格模型的修正. 修正前后梁格模型的計算模態與實測模態的對比見表3. 修正前后各參數的調整情況見表4. 由表3可知修正后的模型各階頻率的計算值與實測值絕對誤差均在5%以內. 這表明選取的參數與修正的過程是基本合理的.

表3 修正前后梁格模型模態與試驗模態的對比

表4 梁格模型修正前后參數變化

4.3.2梁-殼單元模型的修正

圖8 梁-殼單元模型目標函數的迭代收斂曲線Fig.8 Iterative convergence curve of objective function in beam shell element model

根據式(5)定義的目標函數，利用MATLAB和ANSYS的交互式訪問，實現了梁-殼單元模型的參數化修正. 其中，式(1)的求解采用了MATLAB優化工具箱提供的非線性最小二乘函數“lsqnonlin”. 該函數采用了基于Trust-region-reflective(TRR)的優化算法.

表5 修正前后梁-殼單元模型模態與試驗模態對比

表6 梁-殼單元模型修正前后參數變化

上述梁格模型與梁-殼單元模型的修正均以試驗一的模態結果為參考. 實測模態與修正后的計算模態對比見表1.

4.4 修正結果分析

由表3和表5可知，經過參數化修正后，無論是梁格模型或精細化梁-殼單元模型與實測頻率的吻合度均較為良好. 但由表4和表6可知，兩種模型修正前后參數的變化量差異顯著. 具體分析如下.

針對梁格模型，其修正參數特別是主梁及橋臺的虛擬橫梁彈性模量增大了69倍. 這表明梁格法中橫梁采用基于真實截面與材料特性換算得到的抗彎剛度EI，難以準確地模擬結構的扭轉剛度[13]. 針對梁-殼單元模型： ①修正后，主梁彈性模量增大了18%，可能補償了初始模型中將箱型主梁板厚用均值代替引起的誤差； ②修正后，橋墩和橋臺彈性模量分別增大了34.2%和56.6%，可能補償了采用土彈簧簡化方法模擬結構-土相互作用引起的誤差； ③修正后，橫隔梁彈性模量降低了78%，可能因為真實結構橫隔與圖紙存在差異，而修正后的橫隔梁彈性模量值補償了該誤差； ④修正后護欄彈性模量增大了42.9%，考慮到初始值為材料真實彈性模量的70%，修正后數值接近結構的真實彈性模量，說明該橋護欄連續性較好、剛度貢獻較為突出[7, 14-15].

5 結語

1) 基于梁格法建立的動力有限元模型，對扭轉模態存在顯著的固有模型結構誤差. 經過模型參數修正后，雖然能獲得與實測頻率較為接近的結果，但虛擬橫梁的物理參數顯著不合理，難以作為考慮環境因素影響的長期健康監測有限元基準模型.

2) 精細化梁-殼單元模型能較為合理的擬合實測模態結果，包括自振頻率和振型. 其上部結構物理參數基本合理，可考慮引入溫度效應對材料屬性的影響，以此得到消除溫度影響的高精度健康監測基準有限元模型. 相關問題有待深入研究.

3) 精細化基準模型的建立，需要綜合考慮模型的精度和效率. 對于該整體橋，建立的梁-殼單元模型，通過選取合理的修正參數，采用基于參數靈敏度分析的模型修正方法，其計算結果精確、可靠.