基于扭轉-滑移耦合約束的梁橋抗傾覆滑移整體穩定分析

童上航, 肖汝誠, 莊冬利, 孫 斌

(同濟大學 橋梁工程系, 上海 200092)

0 引 言

梁橋傾覆倒塌事故屢有發生,造成了嚴重的生命財產損失．車輛超載是引起梁橋傾覆的直接原因,據筆者不完全統計,近年來我國因車輛超載而倒塌的橋梁見表1,這些事故中的橋梁均為獨柱墩連續梁橋,并且事故車輛都是靠邊緣行駛導致傾覆扭矩過大,梁體發生滑移傾覆而落地．

為防止此類事故發生,眾多學者致力于提高梁橋抗傾覆的整體穩定性能,并取得了一定成果．然而梁橋傾覆滑移機理尚未明確,計算方法仍不完善．目前梁橋傾覆計算方法主要分為剛體法、變形體法和精細化有限元法．剛體法將上部結構視為剛體,支座視為理想點支撐,選取荷載側梁橋最外側兩支座的連線作為傾覆軸,驗算抗力矩和傾覆力矩大小[1]．剛體法忽略了梁體變形和支座變形,會高估梁橋的抗傾覆能力[2]．變形體法認為梁橋會發生小變形,利用支反力驗算抗傾覆能力．比如美國橋梁規范(AASHTO, 2007)通過控制最小支反力來保證不發生傾覆[3];國內交通運輸部2018版新規范(JTG 3362—2018)將傾覆劃分為兩個特征狀態,第一個狀態控制支反力,第二個狀態采用安全系數法[4]．彭衛兵等[5]基于變形體法提出了考慮支座尺寸影響的實用計算方法．然而,變形體法不能精確考慮梁橋大位移引起的幾何非線性問題和支座與上部結構復雜相互作用引起的約束非線性問題．一些學者采取精細化有限元法模擬了傾覆過程[6-9],但精細化有限元建模復雜,計算成本高,不適合工程應用．可見上述三種計算方法各有不足,因此,建立一種兼顧分析精度和計算效率的改進方法顯得尤為重要．

一些學者已嘗試建立改進方法來計算部分支座失效后梁橋的平衡狀態,例如Shi等[10]將等效慣性矩應用于Euler-Bernoulli梁并假定支座不發生滑移,建立了一種簡化模型去分析梁體與支座脫空后的平衡狀態．事實上,支座滑移是梁橋整體失穩破壞的控制條件之一,并且支座摩擦因數直接關系到主梁的抗滑移能力[11-13]．由此可知,建立改進方法應明確支座與上部結構的相互作用機理,著重考慮支座脫空滑移對梁橋整體失穩的影響．鑒于此,本文提出了一種基于曲梁單元考慮扭轉-滑移耦合約束的梁橋整體穩定性分析方法,該方法可計算各種支座失效情況下的梁橋平衡狀態,并準確判斷梁橋整體失穩模式．具體來說,本文工作如下:推導了七自由度曲梁單元剛度矩陣,建立了各種支座失效情況下的支座約束方程,并將支座約束方程引入有限元總方程參與求解．將梁橋整體失穩分為傾覆和滑移兩種情況,提出了以各支座受力狀態判斷梁橋整體是否失穩的計算流程,并編制了相應的計算程序．以簡支超靜定曲梁為例,對比解析解、ANSYS中含翹曲自由度的Beam 188單元計算結果和本文所提出七自由度曲梁單元的計算結果,檢驗了本文所提出曲梁單元的精度．利用所提出方法評估了某匝道橋事故案例,分析了在不同荷載放大率下的支反力大小變化及支座失效情況,還原了失穩發展歷程．通過對比本文所提出方法和傳統桿系有限元模型的支反力結果,驗證了本文方法能更精確地模擬各種支座失效情況下的梁橋平衡狀態．

1 曲梁單元剛度矩陣

文獻[14]推導了包含豎向撓度和扭轉角的曲梁單元剛度矩陣,本文所提出的支座約束方程中包含梁的橫向位移,而曲梁縱向位移與橫向位移的剛度是相互耦合的．因此本文在文獻[14]的基礎上引入橫向位移和縱向位移,推導了七自由度曲梁單元剛度矩陣．

在曲梁的剪切中心軸線上建立直角流動坐標系,其中z軸沿軸線切線方向,x軸指向軸線曲率中心,y軸垂直xOz平面向下．u,v,w分別為x,y,z方向的位移,φ為截面扭轉角．截面內力分量和流動坐標系如圖1所示．設梁的截面面積為A,彈性模量為E,剪切模量為G,繞x軸的慣性矩為Ix,繞y軸的慣性矩為Iy,扭轉慣性矩為Id,扇形慣性矩為Iω．

圖1 曲梁截面內力分量與流動坐標系Fig. 1 The internal force components and the co-rotational coordinate system of the curved beam

曲梁單元一個結點有七個自由度,包括縱向位移w、橫向位移u及其導數u′(對縱坐標z的導數)、豎向位移v及其導數v′、扭轉角φ及其導數φ′．設單元長度為2λ,曲率半徑為R,如圖2所示．

圖2 七自由度曲梁單元Fig. 2 The curved beam element with 7 degrees of freedom

令自然坐標ξ=z/λ,縱向位移w采用線性位移模式,橫向位移u、豎向位移v和扭轉角φ都采用三次位移模式:

(1)

以上位移模式結合結點位移δi=[wiuiu′iviv′iφiφ′i]T和δj=[wjuju′jvjv′jφjφ′j]T,可解得a0,a1,b0,b1,b2,b3,c0,c1,c2,c3,d0,d1,d2,d3等系數．將這些系數代回方程(1),可推導出單元形函數矩陣N:

(2)

其中形函數

曲梁的幾何方程、截面內力與變形的關系如下[14-16]:

(3)

(4)

其中N為截面軸力,My和Mx分別是沿y方向和沿x方向的彎矩,Td為自由扭轉矩,Tω為約束扭轉矩,B為雙力矩．將式(2)、(3)代入式(4),解得廣義應力與結點位移δ的關系:

[NMyMxTdBTω]T=DBδ,

(5)

其中彈性矩陣D為

廣義應變矩陣B為

B=

其中N′1表示N1對z求導．

由虛功原理和變分原理可解得曲梁單元剛度矩陣為

(6)

當曲梁單元的半徑趨于無窮大時,其單元剛度矩陣便退化為直梁的單元剛度矩陣．

記外力集中荷載為F=[PxPyPzMxMyMz]T,分布荷載為f=[qxqyqzmxmymz]T,所有外力荷載分量的正方向與流動坐標系正方向一致,可得單元結點荷載為

(7)

其中

N′3表示N3對z求導．

2 支座約束方程與梁橋失穩判斷

當忽略支座尺寸對梁橋傾覆的影響時,支座視為一個點支撐上部結構:對于單支座墩,其失效模式可分為脫空和滑移;對于雙支座墩,其失效模式可分為單側脫空不滑移、單側脫空且滑移和雙支座都脫空．本文提出的迭代求解思路是先計算所有支座未失效情況下的有限元解,然后檢查每個支座處的支反力和扭轉角．如果支反力受拉,則支座脫空,去除該支座;如果扭轉角大于摩擦角,則支座發生滑移,應引入扭轉-滑移耦合約束方程．然后求解改變支座約束條件后的有限元總方程．

2.1 支座約束

2.1.1 支座未脫空且未滑移

對于單支座墩,梁橋橫向位移和豎向位移都被支座限制,因此約束條件為

(8)

而對于雙支座墩,不僅約束橫向位移和豎向位移,還會約束梁的扭轉以及曲梁平面內彎曲,即

(9)

2.1.2 單墩支座脫空和雙支座墩完全脫空

支座未脫空未滑移時,為實現約束方程(8)、(9),可采用對角元素改1法[17],并將結點載荷列陣相應元素改為0．支座僅承壓,如果支座支反力計算結果為受拉,說明支座發生脫空,應去除該位置的支座約束條件．直接在總剛度矩陣和結點載荷列陣中還原相應的元素即可．

2.1.3 需要引入約束方程的支座失效情況

需引入約束方程的失效情況包括單墩支座滑移、雙支座墩單側脫空不滑移和雙支座墩單側脫空且滑移．以單墩支座發生滑移為例,當梁橋在單墩支座結點處的扭轉角大于摩擦角時,梁橋在該支座處會發生滑移．設支座給梁橋底面法向力為N,摩擦力為f,忽略摩擦力在縱橋向的分量,摩擦因數為μ,梁橋截面剪心到底面的距離為hb,剪心OSC在截面內的位移r=(u,v)．為方便在有限元求解方程中引入扭轉-滑移耦合約束,規定Px為支座反力沿x方向的分量,Py為沿y方向的分量,mz為支座反力給發生位移后剪心位置的扭矩,Px,Py和mz將作為待求量參與有限元總方程求解,如圖3所示．

① 扭轉角大于0(即與流動坐標系正方向一致)

梁橋在支座處結點的橫向位移u,豎向位移v,扭轉角φ滿足幾何約束條件為

g1=utanφ-v+hb(1-secφ)=0,

(10)

支反力給剪心的扭矩為

mz=Px(v-hb)-Pyu,

(11)

因此第二個約束方程為

g2=-mz+Px(v-hb)-Pyu=0．

(12)

摩擦力f和法向力N滿足

(13)

由此可得第三個約束方程為

g3=Py(sinφ-μcosφ)+Px(cosφ+μsinφ)=0．

(14)

(a) φ>0 (b) φ<0

② 扭轉角小于0

同理可得三個約束方程為

(15)

對于其他需要引入約束方程的支座失效情況,其推導過程與單支座墩滑移類似．所有的約束方程匯總見表2．

表2 支座約束方程

2.2 梁橋失穩判斷

在求解過程中,先計算所有支座都有效的位移解,根據支反力正負和扭轉角大小判斷各支座是否失效(脫空或滑移),然后改變各支座約束條件,在有限元求解總方程中引入支座約束方程重新求解,再檢查各支座的失效模式,如果失效模式發生改變則修改約束方程,則再次重新求解,直至各支座失效模式不改變為止．最后根據各支座的失效情況判斷梁橋是否失穩(整個上部結構滑移或傾覆)．

2.2.1 利用Newton-Raphson迭代法求解含約束條件的有限元方程

對于支座未滑移未脫空、單支座脫空和雙支座完全脫空這三種支座失效模式,只需修改總剛度矩陣和結點載荷列陣即可．而對于單支座滑移、雙支座單側脫空不滑移和雙支座單側脫空且滑移,則需要引入約束方程,支反力Px,Py和支反力對剪心扭矩mz將作為待求量引入結點位移列陣共同參與求解．當有k個支座需要引入約束方程時,有限元求解總方程為

(16)

其中

利用Newton-Raphson迭代法求解方程(16),即迭代求解

a(k)=a(k-1)-J-1(a(k-1))F(a(k-1)),

(17)

直至增量小于精度ε:

‖J-1(a)F(a)‖<ε,

(18)

其中J(a)為F(a)的Jacobi矩陣,

2.2.2 梁橋失穩判斷方法

上部結構可能發生的失穩模式有滑移和傾覆．當所有支座都發生了脫空或滑移時,上部結構將發生整體滑移．如果沒發生整體滑移,但承壓支座僅剩一個或兩個,那么上部結構將發生傾覆．實際工程中支座的容許扭轉角較小,摩擦角較小,滑移往往會先于傾覆發生．

綜上所述,梁橋失穩判斷流程為:

我國體育教師教育培訓主要由職前教育系統和在職培訓系統兩大系統組成。體育專業院校和師范院校的體育院系主要承擔體育教師的職前教育，各級各類的教育學院和教師進修學校主要承擔體育教師的在職培訓。目前，我國高校中體育教師的培養還存在著諸多問題，其中“訓練”比“教育”重要的思想使得體育教師職前教育的課程體系缺乏“師范性”。而體育教師正要通過體育教育學類的課程來學習、掌握和運用來實現其專業的壟斷性和不可替代性。因此，建立和完善體育教育學科的課程體系是提高體育教師教育質量、保障體育教師專業地位的重要條件。

1) 假定所有支座有效,求解第一次有限元方程

Kδ=P．

2) 根據各支座的扭轉角和支反力判斷各支座的失效模式．對于單支座墩,如果支反力受拉,則支座脫空;如果支反力受壓但扭轉角大于摩擦角,則支座發生滑移．對于雙支座墩,如果兩個支座的支反力都受拉,則兩個支座都脫空;如果只有一個支座受壓,則是單側支座脫空,然后根據扭轉角是否大于摩擦角判斷是否發生滑移．

3) 如果所有支座都發生脫空或滑移,則上部結構發生滑移,退出計算．

4) 如果承壓的支座只剩一個或兩個,則上部結構發生傾覆,退出計算．

5) 如果上部結構未發生滑移或傾覆,則判斷與上一次求解相比支座失效模式是否發生改變:如果未發生改變則退出計算,梁橋未發生整體失穩;反之則根據支座失效模式修改支座約束方程G和總剛矩陣Keq,利用Newton-Raphson迭代法求解有限元方程

然后返回步驟2)．

3 算 例 分 析

3.1 簡支超靜定曲梁分析

為檢驗曲梁單元的精度,以簡支超靜定曲梁為例,對比理論解析解、ANSYS中含翹曲自由度的Beam 188單元計算結果和七自由度曲梁單元計算結果．在均布豎向荷載和均布扭矩作用下,簡支超靜定曲梁的扭轉角φ和豎向撓度v解析解為[14]

(19)

表3 各參數取值

圖4 撓度的理論結果與有限元結果比較

3.2 某匝道橋事故分析

某匝道橋發生倒塌,事故發生時有3輛重約110 t的貨車行駛在曲梁靠外側(近傾覆側),1輛較輕的車在靠內側(計算中忽略)．倒塌曲梁段長75 m,曲線半徑220 m,平面布置和箱梁橫斷面如圖6、7所示．單墩支座摩擦因數取0.3[18],聯端雙支座摩擦因數取0.01,主梁材料為C50混凝土．

圖6 平面布置圖(單位: m)

定義荷載放大率為施加汽車荷載的比例系數．車輛荷載的荷載放大率從0.1開始逐漸增大,計算所有支座的豎向支反力,如圖8所示．當荷載放大率為0.45時,聯端內側支座P1-2和P4-2發生脫空;放大率為0.55時,聯端外側支座P1-1和P4-1發生滑移;放大率達到0.93時,箱梁最大扭轉角16.8°,已超過單墩支座摩擦角,單墩支座P2和P3發生滑移,此時所有未脫空支座都發生滑移,箱梁在支座處受到的水平力指向外側,箱梁整體往外側滑移,還沒有達到事故車輛荷載大小就已經倒塌．因此,在事故車輛荷載作用下,聯端內側支座發生了脫空,并且其他所有支座發生了滑移,箱梁整體發生滑移然后沿P2和P3連線傾覆,扭轉角發散．由式(15)第三個約束方程可知,墩柱頂受到的水平力會隨著扭轉角發散急劇增大,墩柱底彎矩過大發生強度破壞,事故現場P2和P3墩柱也確實發生墩底破壞而倒塌[8]．

圖8 不同荷載放大率下各支座的豎向支反力Fig. 8 The vertical support reaction of each bearing under different load amplification ratios

規范方法基于傳統的桿系模型驗算抗傾覆能力[4],但傳統桿系模型假定所有支座有效,不能模擬上部結構與支座的復雜相互作用,而本文所提出的方法可以模擬支座脫空和支座滑移帶來的影響．根據圖8,選取荷載放大率為0.1,0.5和0.8三種工況,分別對應所有支座有效、聯端內側支座脫空和聯端外側支座滑移三種情況,對比本文方法和傳統模型的豎向支反力計算結果,如圖9所示．在0.1倍車載作用下,兩者的豎向支反力最大相差2.80%;在0.5倍車載作用下,本文方法判斷出聯端內側支座P1-2和P4-2發生脫空,支反力為0,而傳統模型算得支座P1-2和P4-2承受拉力;在0.8倍車載作用下,本文方法結果中的各支反力與傳統模型相比都有顯著差異,可見支座脫空滑移對支反力影響較大．傳統模型不能準確計算上部結構與支座的相互作用,無法模擬實際事故中支座脫空和支座摩擦滑移的現象,而本文提出的模型克服了其不足之處,能夠更精確地模擬各種支座失效情況下的梁橋平衡狀態．從計算效率來看,以0.1倍車載工況為例,本文方法計算耗時0.11 s(采用MATLAB計算),傳統桿系模型耗時0.52 s(采用MIDAS/Civil計算),可以預計本文方法對于更復雜的梁橋結構具有更高的計算效率．

圖9 所提出計算模型與傳統模型的豎向支反力對比Fig. 9 Comparison of the vertical support reactions given by the proposed model with the traditional model

4 結 論

1) 本文對比了簡支超靜定曲梁的解析解、Beam 188單元有限元解和本文所提出七自由度曲梁單元的有限元解．從位移的計算結果來看,本文所提出的曲梁單元比Beam 188單元更接近解析解,驗證了曲梁單元的高精度性．

2) 本文分析了某匝道橋倒塌事故．在不同荷載放大率條件下通過計算各支座的支反力大小和判斷各支座失效情況,還原了該匝道橋失穩發展歷程:在事故車輛荷載作用下,聯端內側支座脫空,其他所有支座發生滑移,因此該梁橋先發生滑移然后傾覆落地．通過對比分析本文方法與傳統桿系有限元模型在所有支座有效、聯端內側支座脫空和聯端外側支座滑移三種情況下的計算結果,驗證了本文方法能更精確地模擬各種支座失效情況下的梁橋平衡狀態．