蔣望濤 姜海強 馬勤國? 李永東 李治國

(1.華南理工大學 土木與交通學院/亞熱帶建筑科學國家重點實驗室/華南巖土工程研究院，廣東 廣州 510640；2.中鐵隧道勘察設計研究院有限公司，廣東 廣州 511400)

寒區(qū)隧道開挖后含水圍巖在低溫環(huán)境中會凍結膨脹，凍結圍巖邊界受到襯砌及未凍結圍巖的約束而產(chǎn)生凍脹力，隧道開挖與圍巖凍脹會造成圍巖凍脹損傷并導致圍巖應力重分布，進而演變形成不同的應力區(qū)，增大寒區(qū)隧道結構與圍巖的變形，嚴重影響寒區(qū)隧道的穩(wěn)定性[1- 5]。因此，開展寒區(qū)隧道彈塑性解答的研究對于寒區(qū)隧道的穩(wěn)定性具有重要意義，而圍巖凍脹損傷及凍脹特性是計算寒區(qū)隧道圍巖彈塑性解的重要基礎。

目前考慮損傷的隧洞彈塑性解在礦井、地下巷道等得到了廣泛應用[6- 9]，但在寒區(qū)隧道中研究較少。代表性研究成果有：劉紅巖等[10- 11]考慮巖體凍融循環(huán)損傷和彈性本構得到寒區(qū)隧道圍巖凍脹力的理論計算方法；譚賢君等[12]基于熱動力學、連續(xù)介質(zhì)力學以及分凝勢理論建立了圍巖溫度-滲流-應力-損傷(THMD)多場耦合模型，并開發(fā)相關數(shù)值計算程序求得多場耦合條件下的寒區(qū)隧道圍巖凍脹力。寒區(qū)隧道圍巖的凍脹特性是彈塑性求解的另一重要因素，而已有研究基于均勻凍脹特性，結合M-C準則、D-P準則以及粘彈性理論得到寒區(qū)隧道應力、應變以及位移的解析式[3- 4，13- 14]。然而寒區(qū)隧道由于冷空氣流入沿著圍巖徑向方向表現(xiàn)為單向凍結狀態(tài)，夏才初等[15]利用對比試驗驗證了飽和巖石在單向凍結條件下具有不均勻凍脹性，且不均勻凍脹性表現(xiàn)為橫觀各向同性凍脹；Lv等[16]推導了考慮橫觀各向同性凍脹的寒區(qū)隧道圍巖M-C彈塑性解；Feng等[17]基于橫觀各向同性特性推導出5類D-P準則下的寒區(qū)隧道圍巖彈塑性解。因此，寒區(qū)隧道圍巖不僅有凍脹損傷特征，也有不均勻凍脹的特點，但鮮見文獻將三者結合起來進行相關理論研究。

此外，已有寒區(qū)隧道圍巖的彈塑性解大部分基于Mohr-Coulomb(M-C)和Druck-Prager(D-P)準則而建立，無法詳細描述中間主應力對凍結圍巖強度的發(fā)揮效果，所得結果偏向保守或偏向危險。統(tǒng)一強度理論在考慮中間主應力效應的基礎上，通過調(diào)節(jié)中間主應力效應反映巖土類各種強度準則[18]，具有很好的工程應用性。代表性的研究有：張常光等[19]考慮橫觀各向同性凍脹建立了寒區(qū)隧道位移及應力的塑性統(tǒng)一解；曾開華等[20]考慮中間主應力效應和圍巖剪脹特性建立了深埋圓形隧道彈塑性統(tǒng)一解；張長安[21]基于統(tǒng)一強度理論對比分析了M-C準則、廣義M-N準則及外接圓D-P準則的強度效應。

因此，本研究基于損傷理論推導出考慮損傷的寒區(qū)隧道統(tǒng)一強度理論，建立了彈塑性損傷力學計算模型，考慮中間主應力效應及不均勻凍脹特性求解得到寒區(qū)隧道彈塑性應力場和位移場的統(tǒng)一解；其次，探討圍巖損傷特征、不均勻凍脹特征以及中間主應力效應條件下的圍巖應力場及位移場的分布規(guī)律；最后，分析各影響因素條件下寒區(qū)隧道襯砌上徑向應力及塑性半徑的變化趨勢。本研究成果可為寒區(qū)隧道工程設計提供理論參考。

1 基本假設

基于損傷理論和統(tǒng)一強度理論，考慮寒區(qū)隧道圍巖不均勻凍脹特性影響，求得寒區(qū)隧道應力與位移彈塑性統(tǒng)一解，并做出如下基本假設：

(1)襯砌視為彈性介質(zhì)，圍巖視為滿足統(tǒng)一強度理論的彈塑性介質(zhì)。

(2)隧道縱向長度遠大于隧道橫斷面尺寸，故將隧道視為平面應變問題。

(3)隧道處于靜水壓力作用下，忽略襯砌和圍巖重力的影響。

(4)將圍巖視為不均勻凍脹介質(zhì)，圍巖不均勻凍脹系數(shù)保持常數(shù)不變，隧道縱向的線凍脹率和垂直方向的線凍脹率不相等，即為橫觀各向同性凍脹，但是在凍結過程中各自保持不變。

2 寒區(qū)隧道彈塑性損傷力學模型

寒區(qū)隧道彈塑性損傷力學計算模型如圖1所示，計算模型包括4個區(qū)域：襯砌結構(Ⅰ區(qū)，rL>r≥ra，r為隧道徑向半徑)，凍結塑性損傷區(qū)(Ⅱ區(qū)，rp>r≥rL)，凍結彈性區(qū)(Ⅲ區(qū)，rf>r≥rp)，未凍結彈性區(qū)(Ⅳ區(qū)，r≥rf)。凍結塑性損傷區(qū)(Ⅱ區(qū))和凍結彈性區(qū)(Ⅲ區(qū))圍巖在低溫環(huán)境下發(fā)生凍結膨脹，分別受到襯砌(Ⅰ區(qū))和未凍結彈性區(qū)(Ⅳ區(qū))的約束作用，產(chǎn)生凍脹反力作用在襯砌和未凍結彈性區(qū)上。

從圖1得知：ra為隧道襯砌內(nèi)半徑；rL為隧道襯砌外半徑；rp為凍結塑性損傷區(qū)半徑；rf為凍結彈性區(qū)半徑；σf為作用在襯砌與凍結塑性損傷區(qū)交界處的徑向應力；σp為作用在凍結塑性損傷區(qū)與凍結彈性區(qū)交界處的徑向應力；σh為作用在凍結彈性區(qū)與未凍結彈性圍巖交界處的徑向應力；p0為初始地應力。

圖1 寒區(qū)隧道彈塑性損傷力學計算模型Fig.1 Elasto-plastic damage mechanics calculation model of cold regions tunnel圖1中計算模型的各區(qū)域滿足變形幾何條件：

(1)

式中，εr和εθ分別為徑向應變和切向應變，u為隧道徑向位移，r為隧道徑向半徑。

基于不均勻凍脹特性和彈性理論，凍結圍巖本構方程可表示為[14，16- 17]：

(2)

隧道徑向的線凍脹率ξr和環(huán)向的線凍脹率ξθ滿足下式關系[16，22]：

(3)

式中，k為不均勻凍脹系數(shù)，k>1代表著圍巖的不均勻凍脹性，ξV為圍巖體積凍脹率。

3 寒區(qū)隧道彈塑性統(tǒng)一解

如圖1所示為彈塑性損傷力學計算模型，模型中有3個交界面，交界面處滿足位移及應力連續(xù)邊界方程，得到6個連續(xù)邊界方程如下：

(4)

3.1 考慮損傷的統(tǒng)一強度理論

平面應變狀態(tài)下最大和最小的主應力分別為σθ和σr，中間主應力為σz，滿足σ1=σθ，σ2=σz，σ3=σr，統(tǒng)一強度理論[18- 20]如下：

σθ=Mσr+N

(5)

(6)

式中：M和N為統(tǒng)一強度方程參數(shù)；c為圍巖黏聚力；φ為圍巖內(nèi)摩擦角；b為統(tǒng)一強度理論參數(shù)，取值范圍為0～1，當b=0時，式(5)表示為M-C準則；當b=1時，式(5)退化為雙剪應力準則；當1>b>0時，式(5)代表系列新強度準則。

基本假設條件(4)中假定圍巖為橫觀各向同性介質(zhì)，圍巖損傷采用有效應力表示：

(7)

將式(7)代入式(5)，統(tǒng)一強度理論可表示為

σθ=Mσr+(1-D)N

(8)

圍巖單向壓縮條件下的巖石應力應變曲線可表示為雙線性[6- 9]，如圖2所示。圖中：AB段為彈性階段，彈性模量為E；B點為峰值強度點；BC段為塑性損傷軟化階段，降模量為λ。為了描述超過峰值后的損傷程度，定義損傷程度比值λ/Ef。

圖2 巖石單軸壓縮應力-應變曲線Fig.2 Stress-strain curve of rock under uniaxial compression

因此，圍巖三維損傷演化方程為[6- 9]：

(9)

式中，εc為峰值應變，εi為等效應變。

等效應變εi可由ε1、ε2、ε3三向主應變表示：

(10)

考慮圍巖剪脹特性，由已有研究成果[7]得知塑性損傷區(qū)內(nèi)有效應變εi表達式：

(11)

將式(11)代入式(9)求得圍巖三維損傷演化方程：

(12)

將式(12)代入式(8)求得考慮損傷的統(tǒng)一強度方程為

(13)

3.2 襯砌支護區(qū)

寒區(qū)隧道襯砌支護區(qū)視為厚壁圓筒彈性介質(zhì)，襯砌僅受σf外壓作用，利用彈性理論[23]，求得應力如下：

(14)

襯砌支護結構彈性本構方程[23]為：

(15)

結合幾何方程(1)得到襯砌支護區(qū)的徑向位移表達式：

(16)

3.3 凍結塑性損傷區(qū)

平面應力應變狀態(tài)下，凍結圍巖塑性損傷區(qū)微分平衡方程為[23]：

(17)

將式(13)代入式(17)，利用邊界條件：當r=rL時，σr=σf，得到凍結圍巖塑性損傷區(qū)應力表達式：

(18)

凍結塑性損傷區(qū)圍巖應變包括塑性應變和彈性應變，由式(2)得到如下表達式：

(19)

在凍結塑性損傷區(qū)內(nèi)，圍巖應變分量滿足協(xié)調(diào)方程：

(20)

將式(19)代入式(20)，得到式(21)：

(21)

由線性非相關流動性法則[24]，得圍巖塑性損傷區(qū)內(nèi)塑性切向應變與徑向應變之間滿足以下關系：

(22)

根據(jù)式(22)線性非相關流動法則，結合式(21)，得到

(23)

將式(18)代入式(2)得到凍結圍巖塑性損傷區(qū)彈性應變分量表達式：

(24)

將式(24)代入式(23)得到凍結圍巖塑性損傷區(qū)切向應變?nèi)缦拢?/p>

(25)

(26)

將式(26)代入式(25)得到凍結圍巖塑性損傷區(qū)塑性切向應變表達式為

(27)

將式(24)和(27)代入式(1)，以此得到凍結圍巖塑性損傷區(qū)徑向位移表達式為：

(28)

3.4 凍結彈性區(qū)

采用彈性力學理論[23]求得凍結圍巖彈性區(qū)應力表達式：

(29)

將式(29)代入式(2)得到凍結圍巖彈性區(qū)應變?yōu)椋?/p>

(30)

利用式(30)和幾何方程(1)，得到凍結圍巖彈性區(qū)徑向位移：

(31)

3.5 未凍結彈性區(qū)

將未凍結圍巖區(qū)視為厚壁圓筒彈性介質(zhì)，利用彈性理論[23]，求得應力為：

(32)

未凍結圍巖彈性本構方程[23]表示為：

(33)

同理，相應的徑向位移為

(34)

3.6 求解步驟

為了求解σf、σp、σh以及rp4個未知量，根據(jù)力學計算模型中6個邊界條件，采用以下求解步驟：

(2)在r=rL、r=rp及r=rf交界處，式(4)中剩余4個邊界條件，可化簡表示包含4個方程的方程組。

(3)利用MATLAB軟件對步驟(2)中的方程組進行求解，求出σf、σp、σh以及rp4個未知量的值。

4 算例分析

4.1 統(tǒng)一解驗證

基于損傷理論和統(tǒng)一強度理論，考慮圍巖不均勻凍脹特性建立寒區(qū)隧道彈塑性統(tǒng)一解。所得的彈塑性統(tǒng)一解包含統(tǒng)一強度理論參數(shù)b、不均勻凍脹系數(shù)k和圍巖體積凍脹率ξV、損傷降模量λ以及位移釋放系數(shù)η等因素的共同影響。隧道相關計算參數(shù)見表1。

表1 寒區(qū)隧道幾何和材料特性表[14，17]1)Table 1 Geometrical and material properties of cold regions tunnel[14，17]

當b=0、k=1、λ=0以及η=0時，彈塑性統(tǒng)一解退化為文獻[14]中考慮均勻凍脹的M-C準則統(tǒng)一解，將本研究推導的應力彈塑性統(tǒng)一解與文獻[14]的應力解進行對比，如圖3所示。

圖3 本文解析解與文獻[14]解析解的對比Fig.3 Comparison between results from this study and those of the reference[14]

從圖3可以看出，當參數(shù)b=0、k=1、λ=0以及η=0時，本研究推導的應力彈塑性統(tǒng)一解與文獻[14]的應力解完全吻合，驗證了本研究推導的考慮損傷的寒區(qū)隧道彈塑性統(tǒng)一解的正確性。

4.2 圍巖應力場及位移場分析

為分析寒區(qū)隧道應力場及位移場的分布規(guī)律，考慮損傷降模量λ、圍巖體積凍脹率ξV、不均勻凍脹系數(shù)k以及統(tǒng)一強度理論參數(shù)b等4個因素的影響，得到不同影響因素條件下的寒區(qū)隧道應力及位移分布曲線。

4.2.1 圍巖損傷特征

圖4示出了不同損傷程度下寒區(qū)隧道應力分布及彈塑性區(qū)位移分布曲線。隨著λ/Ef的增大，支護區(qū)環(huán)向及徑向應力以及凍結彈性區(qū)的切向應力逐漸增大，而凍結塑性損傷區(qū)的切向應力及凍結區(qū)的徑向應力卻不斷減少，其余區(qū)域內(nèi)的應力變化較小。隨著λ/Ef比值增大，相同距離處的位移不斷增大，λ/Ef比值從0.0增至1.5時塑性區(qū)域外壁面的位移保持不變，但塑性半徑向外擴大了0.221 m，塑性區(qū)域內(nèi)壁面位移增大了0.948 mm，且凍結塑性區(qū)內(nèi)零位移半徑增大0.090 m。這表明損傷程度變大不僅能增大寒區(qū)隧道襯砌及圍巖變形，也能擴大凍結圍巖塑性區(qū)半徑，在寒區(qū)隧道工程設計中需要考慮圍巖凍結損傷程度，提高襯砌支護及時性和采取保溫防凍措施能有效減小圍巖損傷程度。

(a)應力分布曲線

(b)彈塑性區(qū)位移分布曲線圖4 損傷程度對應力和位移的影響(k=2，b=0.5，ξV=0.016 5)Fig.4 Influence of damage degree on stress and displacement(k=2，b=0.5，ξV=0.016 5)

4.2.2 圍巖不均勻凍脹特性

圖5示出了不同體積凍脹率下寒區(qū)隧道應力分布及彈塑性區(qū)位移分布曲線。圍巖體積凍脹率的增加會顯著提高寒區(qū)隧道圍巖各區(qū)域的應力，對隧道襯砌及圍巖的安全十分不利。當ξV從0.010 5增大至0.022 5時，凍結塑性區(qū)內(nèi)、外壁面位移分別相對增大了128.36%和172.07%，塑性區(qū)域半徑增大了0.138 m，但塑性區(qū)內(nèi)位移零點保持不變。這表明在寒區(qū)隧道工程中需要詳細考慮凍脹作用，通過采取保溫防凍措施提高寒區(qū)隧道的抗凍性能。

(a)應力分布曲線

(b)彈塑性區(qū)位移分布曲線圖5 體積凍脹率對應力和位移的影響(λ/Ef=0.5，k=2，b=0.5)Fig.5 Influence of frost heave ratio of volume on stress and displacement(λ/Ef=0.5，k=2，b=0.5)

圖6示出了不同不均勻凍脹系數(shù)時的寒區(qū)隧道應力分布及彈塑性區(qū)位移分布曲線?？梢钥闯?，隨著參數(shù)k從1增大至3，作用于襯砌上的凍脹力增大了0.16 MPa，凍結彈性區(qū)的切向應力得到有效減緩，其余區(qū)域應力變化不大。同時，不均勻凍脹系數(shù)k的增大將塑性區(qū)域半徑縮小了0.25 m，其余位置處位移曲線較接近。這說明考慮不均勻凍脹系數(shù)相比均勻凍脹系數(shù)對于寒區(qū)隧道設計更有利，且圍巖不均勻凍脹系數(shù)需要通過試驗測試得到。

(a)應力分布曲線

(b)彈塑性區(qū)位移分布曲線圖6 不均勻凍脹系數(shù)對應力和位移的影響(λ/Ef=0.5，b=0.5，ξV=0.016 5)Fig.6 Influence of coefficients of non-uniform frost heave on stress and displacement(λ/Ef=0.5，b=0.5，ξV=0.016 5)

4.2.3 中間主應力效應

圖7示出了不同統(tǒng)一強度理論參數(shù)b下的寒區(qū)隧道應力分布及彈塑性區(qū)位移分布曲線。隨著參數(shù)b從0.0增加到1.0，襯砌區(qū)及凍結塑性損傷區(qū)的環(huán)向應力不斷減少，其余區(qū)域應力差距不大。這表明考慮中間主應力效應能有效發(fā)揮凍結圍巖的承載潛能。同時，隨著參數(shù)b的增大，凍結塑性損傷區(qū)內(nèi)同一位置處的位移逐漸減少，塑性區(qū)域內(nèi)、外壁位移分別相對減少了42.88%和17.78%，塑性區(qū)域半徑縮小了0.278 m，且塑性區(qū)域內(nèi)零位移半徑減小了0.090 m。這說明考慮中間主應力系數(shù)能有效縮小寒區(qū)隧道的位移及塑性半徑，在寒區(qū)隧道工程設計時不容忽視。

(a)應力分布曲線

(b)彈塑性區(qū)位移分布曲線圖7 統(tǒng)一強度理論參數(shù)對應力和位移的影響(λ/Ef=0.5，k=2，ξV=0.016 5)Fig.7 Influence of parameter of unified strength theory on stress and displacement(λ/Ef=0.5，k=2，ξV=0.016 5)

4.3 參數(shù)分析

從上文分析得知，影響寒區(qū)隧道的彈塑性解的因素包括損傷程度、不均勻凍脹特性以及中間主應力效應等。為更進一步分析本文推導出來的彈塑性統(tǒng)一解參數(shù)構成，詳細分析下述參數(shù)對襯砌支護上的徑向應力σf和塑性半徑rP的影響，參數(shù)包括：損傷程度比值λ/Ef、體積凍脹率ξV、不均勻凍脹系數(shù)k、統(tǒng)一強度理論參數(shù)b、初始地應力p0以及位移釋放系數(shù)η，如圖8-13所示。

圖8 不同λ/Ef的應力及塑性半徑曲線Fig.8 Stress and plastic radius curves with different λ/Ef

圖9 不同ξV的應力及塑性半徑曲線Fig.9 Stress and plastic radius curves with different ξV

圖10 不同k的應力及塑性半徑曲線Fig.10 Stress and plastic radius curves with different k

圖11 不同b的應力及塑性半徑曲線Fig.11 Stress and plastic radius curves with different b

圖12 不同p0的應力及塑性半徑曲線Fig.12 Stress and plastic radius curves of different p0

圖13 不同η的應力及塑性半徑曲線Fig.13 Stress and plastic radius curves of different η

圖8示出了不同損傷程度比值λ/Ef的應力σf和塑性半徑rp曲線，隨著比值λ/Ef從0.00增大至2.00時，應力σf增大了0.904 MPa，相對增加28.10%，塑性半徑從6.041 m變化至6.499 m，相對增加7.58%。這表明凍結塑性損傷區(qū)λ/Ef比值越大，圍巖凍結塑性區(qū)半徑越大，襯砌支護結構越不安全。

體積凍脹率ξV代表著寒區(qū)隧道圍巖的凍結膨脹程度，如圖9所示。當參數(shù)ξV從0.004 5增加到0.028 5時，應力σf從0.802 MPa增大至5.897 MPa，相對增加了6.353倍，塑性半徑相對增加了0.489 m。這說明圍巖體積凍脹率的增加對襯砌支護結構及圍巖凍脹區(qū)域有著明顯影響，當隧道內(nèi)圍巖具有較高的體積凍脹率時，隧道結構由于凍脹而引發(fā)的變形也會越大。

圖10示出了不同不均勻凍脹系數(shù)k對應的應力σf和塑性半徑rp曲線，當參數(shù)k從1.0變化至3.0時，應力σf相對增大了0.163 MPa，塑性半徑卻相對減小了0.250 m。這表明考慮均勻凍脹對于寒區(qū)隧道工程設計偏于危險，而考慮不均勻凍脹是十分有必要的。

圖11示出了不同統(tǒng)一強度理論參數(shù)b對應的應力σf和塑性半徑rp曲線，當參數(shù)b從0.0增加至1.0時，應力σf和塑性半徑rp逐漸減小，其中應力σf相對減小了30.06%，塑性半徑相對減小了0.273 m。進一步說明考慮中間主應力系數(shù)能有效發(fā)揮寒區(qū)隧道圍巖的承載潛能。

圖12示出了不同初始地應力p0對應下的應力σf和塑性半徑rp曲線，隨著p0從2.5 MPa增加至12.5 MPa，應力σf從3.303 MPa增大至6.093 MPa，相對增加了0.844倍，而塑性半徑卻相對減小了0.241 m。這說明初始地應力p0越大，寒區(qū)隧道需要更高強度的襯砌支護。

位移釋放系數(shù)η代表著隧道開挖后支護的及時性。如圖13所示，當參數(shù)η從0.0增加至0.8時，應力σf相對減小了25.12%，而塑性半徑卻相對增大了0.496 m。這表明η越大隧道開挖圍巖初始變形增大，后期襯砌支護凍脹變形量變小，同時圍巖凍結區(qū)更容易出現(xiàn)塑性損傷現(xiàn)象。

5 結論

(1)基于損傷理論模型，考慮中間主應力效應和寒區(qū)隧道圍巖不均勻凍脹特性，建立了寒區(qū)隧道圍巖彈塑性應力場和位移場的統(tǒng)一解。通過與相關文獻寒區(qū)隧道M-C準則應力解對比分析，驗證了本研究推導的統(tǒng)一解的正確性，結果可為寒區(qū)隧道工程設計提供理論參考，具有很好的應用性。

(2)塑性損傷特征、圍巖不均勻凍脹特性以及中間主應力效應對寒區(qū)隧道的應力場和位移場有著顯著的影響。不考慮圍巖塑性損傷及凍脹特征對于寒區(qū)隧道是偏危險的，考慮不均勻凍脹特性更能代表寒區(qū)隧道圍巖實際的凍脹特性，考慮中間主應力效應能有效發(fā)揮凍結圍巖的承載潛能。由此可見，考慮寒區(qū)隧道圍巖的損傷、凍脹特性及中間主應力效應是非常有必要的。

(3)寒區(qū)隧道襯砌支護上徑向應力σf和塑性半徑rp受很多因素影響，增大λ/Ef和ξV會加快塑性半徑的變化，圍巖更容易達到塑性狀態(tài)；增大k和p0對襯砌支護不利，但凍結塑性區(qū)半徑反而變??；增大b能發(fā)揮凍結圍巖的承載潛能，有效縮減凍結塑性區(qū)半徑；增大η意味支護拖延，圍巖更容易轉(zhuǎn)為塑性狀態(tài)，圍巖變形增大。