混凝土溫度自應力數值計算方法

彭曉宇

(石家莊鐵道大學土木工程系 石家莊 050043)

為了保證混凝土橋梁結構的耐久性，控制由溫度應力產生的裂縫十分重要。對于溫度應力的計算可以分為2個部分：構件內部溫度不同的各部分之間相互約束而引起的應力，稱為溫度自應力；超靜定結構中，構件由于溫度變化產生的變形受到其他構件或外部支承的約束而引起的應力，稱為溫度次應力[1]。

對于溫度應力的計算方法，目前有理論分析法、經驗公式法和有限元法[2-3]，其中經驗公式法簡單易算，但是計算精度較差；有限元法雖然計算精度高，但是計算工作量大[4]。而早在1986年劉興發[5]就提出采用有限元差分法來計算溫度場，但是該方法受限于其求解性。

為了能在工程設計中高效準確地計算出溫度自應力，特別是在太陽輻射[6]和氣溫變化等環境因素的影響下，內部溫度呈非線性分布的橋梁結構的溫度自應力，本文對溫度自應力的計算公式進行推導，根據推出的公式編寫出相關計算程序，利用軟件計算以提高設計效率，并對計算程序進行校核，對比有限元軟件Abaqus計算結果驗證程序計算的正確性。

1 溫度自應力推導

1.1 截面圖示與假定

混凝土橋梁大都采用矩形、工字形截面或箱形截面，考慮到箱形截面可以等化為工字形截面，而矩形截面是工字形截面的一個特例，因此公式推導和編程只考慮工字形截面的情況，其截面示意見圖1。

圖1 計算截面示意圖

推導假定如下。

1) 截面為“工”字形梁或箱梁。

2) 在溫度力作用下橫截面保持為平面，即符合平截面假定。

3) 橫截面上各點處的縱向線段均處于單軸應力狀態。

1.2 公式推導

截面在不均勻溫度荷載作用下，其內部軸向力與彎矩的總和均為零。據此可得式(1)、式(2)：

(1)

(2)

式中：σ為截面正應力；dA為積分面域；Ec為混凝土彈性模量；Δε為截面應變的變化量。

采用《鐵路橋涵設計基本規范》附錄B中溫度荷載，溫度在梁各高度處產生總應變表示為

εt=αty=α·t0·e-ay

(3)

式中：α為混凝土的線膨脹系數；a為溫差曲線的指數；t0為溫差。

梁各高度處的實際應變值為

(4)

令

Δε=(εt-ε)

(5)

可得應力：σ=Ec·Δε

對式(1)

(6)

對式(2)

(7)

將式(3)～式(5)代入式(6)、式(7),簡化得

解得

(8)

式中：

A=b1h1-b3h1+b3(h0-h2)-

b2(h0-h2)+b2h0

b3e-a(h0-h2)+b2e-a(h0-h2)-b2e-ah0]

b3e-ah1(ah1+1)-b3e-a(h0-h2)[a(h-h2)+1]-

b2e-a(h0-h2)[a(h0-h2)+1]}

σ=EcΔε=Ec(εt-ε)=

(9)

將式(8)代入式(9)求出

式中：ε為截面應變；A，B，C，D，E，F為計算參數。

2 編寫程序及校核

將上述推導公式用計算機語言編寫成程序，用以高效計算溫度自應力。為了驗證推導公式及編寫程序的準確性，從截面在不均勻溫度荷載作用下，其截面上的軸向力與彎矩的總和均為0出發，求解程序所計算截面的軸向力及截面彎矩，判斷其是否為0。因此，校核計算可利用原程序，將整個橫截面劃分為若干段(設截面被分為NS段)，用溫度自應力計算程序分別計算每段的軸向力及彎矩，并進行求和，若其結果趨近于0，則證明溫度自應力計算程序正確，計算過程如下。

假設有一工字形截面，截面對應的參數如下：b1=1.5 m，h1=0.9 m，b2=2.3 m，h2=0.7 m，b3=0.88 m，h0=3.5 m。并取：α=10-5，a=7，t0=16 ℃，Ec=3.4×104MPa。

將上述驗證思路編寫成程序并進行校核計算。整理程序計算的結果，并探討不同等分塊數對應截面軸向力及彎矩的關系，不同等分塊數與對應的軸向力之和關系見表1及圖2，不同等分塊數與對應的截面彎矩之和關系見表2及圖3。

表1 不同等分塊數對應的軸向力之和

圖2 不同等分塊數對應的軸向力之和關系曲線

NS(等分塊數)∑M=∫Ωσy dA/(kN·m) 20155.8320018.1812 0001.7485 0000.79610 0000.366100 0000.035

圖3 不同等分塊數對應的彎矩之和關系曲線

由表1、2、圖2、3可見，當截面分段數越多，其內力越趨近于0。計算結果與前面分析的截面在不均勻溫度荷載作用下，截面上的軸向力與彎矩的總和均為0相符。這表明工字截面溫度自應力計算程序是正確的。

3 結合有限元軟件計算驗證計算程序

為了進一步證明工字截面溫度自應力推導公式及計算程序的正確性。用Abaqus大型通用有限元軟件建立對應實體模型。并分層賦予溫度以模擬溫度場的分布，然后對比計算結果。

假設一截面b1=b2=b3=1.5 m，h0=2 m，則工字形截面變為寬1.5 m、高2 m的矩形截面，取α=10-5，a=7，t0=16 ℃，Ec=3.4×104MPa，泊松比μ=0.2。

先分別取y=0.5，1 m，代入工字截面溫度自應力計算程序，計算對應位置的溫度自應力，應力計算結果分別為714，374 kPa。

然后用大型有限元軟件Abaqus建立三維有限元實體模型。不同于殼單元和梁單元可以在其厚度或梁高方向設置非線性溫度場，Abaqus實體單元只能賦予一個溫度值，因此，為了能夠模擬非線性溫度場，可以將建立的實體有限元模型沿梁高方向分為20層，然后將每層賦予對應的溫度值。計算的溫度自應力云圖見圖4。

圖4 Abaqus計算應力云圖

提取y=0.5，1 m處的溫度自應力，得應力計算結果為710.2，376.9 kPa。

將程序計算的溫度自應力結果與Abaqus建立實體模型計算的結果進行對比，對比結果見表3。

表3 計算結果對比表

根據以上的校核計算結果，在y=0.5 m和y=1 m位置，本程序計算結果與Abaqus計算結果誤差在1%以下，因而可以認為溫度自應力計算程序正確且具有良好的精度[7]。

4 結論

1) 本文中基于3個假設下推導的溫度自應力計算公式是正確的，其求解結果符合截面在不均勻溫度荷載作用下的受力特性。

2) 本程序計算的溫度應力與Abaqus有限元軟件建立的三維有限元模型計算結果誤差在允許范圍之內。從而證明該計算程序具有良好的精度，可應用于工程實際計算中。

[1] 彭友松,強士中.混凝土橋梁結構溫度自應力計算方法探討[J].西南交通大學學報,2006(4):452-455.

[2] 朱伯芳.大體積混凝土溫度應力與溫度控制[M].北京：中國水利水電出版社，2012．

[3] 朱伯芳,宋敬廷.混凝土溫度場及溫度徐變應力的有限元分析[G]//水利水電工程應用電子計算機資料選編.北京：水利電力出版社，1977.

[4] 李國興.基Matlab的混凝土平面溫度自應力差分解法[J].水電能源科學,2016(4):78-81，27.

[5] 劉興法.預應力混凝土箱梁溫度應力計算方法[J].土木工程學報,1986(1):44-54.

[6] 徐長武,任志剛,霍凱成.太陽輻射作用下鋼管膨脹混凝土界面性能試驗與分析[J].工程力學,2015(8):201-210.

[7] 徐豐，王波，張海龍.混凝土連續箱梁橋溫度應力對比分析[J].交通科技,2008(4):9-12.