候 超 ，靳曉光 ,2,3，何 杰 ，張 馳

（1.重慶大學土木工程學院，重慶 400044；2.重慶大學山地城鎮(zhèn)建設(shè)與新技術(shù)教育部重點實驗室，重慶 400044；3.重慶大學煤礦災害動力學與控制國家重點實驗室，重慶 400044）

據(jù)相關(guān)數(shù)據(jù)統(tǒng)計，永久性凍土區(qū)和季節(jié)性凍土區(qū)約占我國國土面積的3/4，當環(huán)境溫度低于0 ℃時寒區(qū)巖體內(nèi)水分凍結(jié)，溫度高于0 ℃時則消融.隨著季節(jié)變換及晝夜更替環(huán)境溫度在0 ℃上下震蕩，巖體內(nèi)水分發(fā)生周期性的凍結(jié)和消融，即凍融循環(huán)[1-2].隨著“一帶一路”倡議的實施，大量的工程建設(shè)在西部寒區(qū)興起，越來越多的巖土工程遭遇到凍融問題，其中巖石在凍融循環(huán)作用下的損傷模型一直是研究者關(guān)注的問題，對寒區(qū)巖石的損傷劣化機理及巖土工程穩(wěn)定控制具有重要的理論意義[3-6].

目前，對凍融循環(huán)引起的巖石損傷劣化機理已經(jīng)有了較為一致的認識：巖石礦物顆粒在低溫作用下體積收縮，而孔隙和裂隙內(nèi)的水凍結(jié)成冰產(chǎn)生約9.08%的體積膨脹，由于不同礦物顆?？s脹率有所差異，導致顆粒間的縮脹變形不協(xié)調(diào)，礦物顆粒與微孔隙之間產(chǎn)生巨大的凍脹力，對巖石造成一定的損傷[7]；溫度升高時，巖石內(nèi)孔隙、裂隙中的冰體融化，隨著凍脹力的釋放及水分的遷移，進一步促進巖石的損傷；在長期的凍融循環(huán)作用下巖石內(nèi)部微裂紋不斷萌生和發(fā)展，造成巖石承載力下降[8-9].

迄今為止，國內(nèi)外學者從宏觀、微觀及細觀等多方面對凍融巖石的損傷劣化機理進行了大量研究：在凍融巖石損傷模型方面，張慧梅等[10]基于巖石內(nèi)部缺陷分布的隨機性，運用損傷力學理論建立了溫度-荷載耦合作用下巖石的損傷模型，并通過試驗驗證模型的合理性；Huang 等[11]假設(shè)巖石微元體強度服從Weibull 分布，以彈性模型作為損傷變量結(jié)合最大拉應變準則建立了凍融巖石的損傷模型，對模型進行了驗證，并應用于寒區(qū)隧道的穩(wěn)定性分析；Jia 等[12]將凍融作用視為疲勞損傷，建立了凍融巖石的疲勞損傷模型，揭示了反復凍融下巖石力學強度衰減機理；為實現(xiàn)對凍融巖石變形全過程的預測，袁超等[13]選用Drucker-Prager 準則，建立了能反映凍融巖石各變形階段特征的損傷模型；張峰瑞等[14]對經(jīng)歷不同凍融循環(huán)次數(shù)的花崗巖開展細觀特征分析和剪切蠕變試驗，基于試驗結(jié)果提出了凍融巖石損傷黏性元件，構(gòu)建了花崗巖凍融剪切蠕變本構(gòu)模型.

總體看來，對巖石的凍融損傷模型研究較少，模型多基于莫爾-庫倫準則及Drucker-Prager 準則，不能準確反映巖石的張拉破壞特性和峰后強度特征.眾所周知，巖石的抗拉強度遠低于抗壓強度，巖石微孔隙內(nèi)水結(jié)成冰形成的凍脹力可視為拉應力，巖石在凍融和受荷作用下易發(fā)生張拉破壞.因此，本文基于最大拉應變破壞準則和統(tǒng)計損傷理論建立了考慮凍融和荷載耦合作用的巖石損傷模型，推導了模型參數(shù)的理論解，利用已有研究結(jié)果對本文建立的損傷模型進行了驗證，探討了不同凍融次數(shù)下巖石的總損傷曲線演化規(guī)律，對模型參數(shù)進行了分析，并采用數(shù)值模擬方法計算了凍融作用對隧道工程穩(wěn)定性的影響規(guī)律.

1 損傷本構(gòu)模型建立

1.1 損傷變量及本構(gòu)方程

基于Lemaitre 應變等價準則[15]，巖石在荷載作用下的有效應力為

式中： σi為名義應力；D為受荷損傷變量.

天然巖石中含有大量微孔隙和裂隙，為了表征巖石微元體的強度引入Weibull 分布函數(shù)，則巖石微元體的強度概率密度函數(shù)可表示為

受荷損傷變量D為失效單元數(shù)Nf與總單元數(shù)N的比值，如式（3）.

則可得到受荷損傷變量為

參照張慧梅等[16]提出的凍融損傷模型，將凍融作用下承受荷載的巖石劃分為凍融損傷部分、受荷損傷部分和未損傷部分，則

式中： σr為損傷部分承受的殘余應力；A為受荷面積；A1為未損傷部分面積；An為凍融損傷部分面積；A2為受荷損傷部分面積.

巖石的凍融損傷變量、受荷損傷變量和總損傷變量如式（8）～（10）.

聯(lián)立式（8）～（10）可得

聯(lián)立式（1）和式（6）、（7）、（11）可得

假定巖石只在主應力方向發(fā)生損傷，且未損傷部分服從廣義胡克定律，結(jié)合巖石變形協(xié)調(diào)關(guān)系得到凍融作用下承受三軸壓縮荷載巖石在主應力方向的損傷模型為

式中：En為巖石凍融循環(huán)n次后的彈性模量； μn為巖石凍融循環(huán)n次后的泊松比； σ3為圍壓； σ1為主應力方向應力；ε1為主應力方向應變.

由宏觀損傷力學可知，巖石的凍融損傷變量為

式中：E0為未受凍融作用巖石的彈性模量.

將式（5）、（14）代入式（11）、（13）得

假設(shè)巖石微元體強度服從最大拉應變準則，如式（16）.

將式（16）代入式（15）得

1.2 模型參數(shù)確定

參數(shù)f0、m可用峰值應力法確定[11]，峰值點處應力、應變?yōu)?/p>

式中： σ1為軸向應力.

峰值點處應力對應變的偏微分等于0，則

將式（18）、（19）代入式（17）得

將式（21）代入式（20）得

由式（22）可得

則

將式（24）代入式（21）得

由式（25）及Weibull 分布函數(shù)參數(shù)定義可得

由式（24）得

最后，引入修正系數(shù)k對式（17）進行修正得

k的取值為

式（28）及式（29）即為本文所建立的基于最大拉應變破壞準則并考慮凍融和荷載耦合作用的寒區(qū)凍融巖石損傷模型.

2 模型驗證

2.1 實例一：低孔隙率硬質(zhì)花崗巖

Tan 等[17]對采自西藏嘎隆拉山高速公路隧道中的花崗巖開展了飽和狀態(tài)下，溫度變化范圍為-40～+ 40 ℃，循環(huán)次數(shù)高達150 次的凍融循環(huán)試驗.文獻[17]所采用的花崗巖孔隙率為0.67%，平均單軸抗壓強度為135.73 MPa，屬于致密的硬質(zhì)巖石.為驗證本文所建立損傷本構(gòu)模型的適用范圍及合理性，利用文獻中圍壓為10 MPa，凍融循環(huán)次數(shù)為0、50、100、150 次的試驗曲線與本文損傷模型計算的理論曲線進行對比分析.同時，基于Huang 等[11]建立的損傷本構(gòu)模型計算了低孔隙率硬質(zhì)花崗巖的應力-應變曲線，將其作為對比曲線之一.

表1 給出了不同凍融次數(shù)下低孔隙率硬質(zhì)花崗巖的力學參數(shù).表1 中，N為凍融次數(shù).由于凍融作用對硬質(zhì)巖石的泊松比影響程度較低，因此，本文在計算過程中，假設(shè)不同凍融次數(shù)下低孔隙率硬質(zhì)花崗巖的泊松比數(shù)值不變，均為0.15.試驗與理論模型的對比結(jié)果如圖1 所示，圖中：理論模型1 代表本文所建立的損傷本構(gòu)模型；理論模型2 代表文獻[11]中的損傷本構(gòu)模型.

表1 低孔隙率硬質(zhì)花崗巖力學參數(shù)Tab.1 Mechanical parameters of hard granite with low porosity

從圖1 中可以看出：當花崗巖凍經(jīng)歷0、50、100 次凍融循環(huán)作用后，試驗數(shù)據(jù)曲線與理論模型1、理論模型2 在峰值強度前均具有較好的對應關(guān)系；不同的是，在峰后階段，理論模型1 可較好反映出巖石的殘余強度，而理論模型2 的應力在峰值強度后跌落速率較快，無法體現(xiàn)巖石的殘余強度特征；當凍融循環(huán)次數(shù)達到150 次時，理論模型2 在峰值強度前與試驗數(shù)據(jù)吻合度較理論模型1 高，此時，雖然理論模型1 的峰值強度與試驗數(shù)據(jù)曲線之間存在一定的差距，但曲線整體趨勢較為一致.可見，本文所建立的凍融巖石損傷本構(gòu)模型可以適用于低孔隙率硬質(zhì)巖石，且能較好地體現(xiàn)巖石的殘余強度特征.

2.2 實例二：高孔隙率軟質(zhì)紅砂巖

張慧梅等[16]對未受凍融循環(huán)作用及受凍融循環(huán)后的紅砂巖開展了三軸力學試驗，分析了不同凍融循環(huán)次數(shù)和圍壓下紅砂巖的力學特性.文獻[16]所采用的紅砂巖孔隙度為14.26%，單軸抗壓強度約為4 MPa，屬于高孔隙率軟巖[18-19].為驗證本文所建立的損傷本構(gòu)模型在高孔隙率軟巖中的適用性，利用文獻[16]中的試驗數(shù)據(jù)曲線與本文計算的理論曲線進行對比分析.類似的，此處也加入了與文獻[11]中的損傷本構(gòu)模型的對比.表2 給出了不同凍融次數(shù)和圍壓下高孔隙率軟質(zhì)紅砂巖的力學參數(shù).試驗與理論模型的對比結(jié)果如圖2 所示.

圖2 高孔隙率軟質(zhì)紅砂巖試驗與理論曲線對比Fig.2 Comparison of experimental and theoretical curves of soft red sandstone with high porosity

表2 高孔隙率軟質(zhì)紅砂巖力學參數(shù)Tab.2 Mechanical parameters of soft red sandstone with high porosity

由圖2 可知：對于經(jīng)歷不同凍融作用的高孔隙率軟質(zhì)紅砂巖，在峰值應力前，試驗數(shù)據(jù)與理論模型1、理論模型2 的吻合性均較好；在峰值強度后階段，理論模型2 應力快速跌落，與試驗數(shù)據(jù)中巖石具有一定的殘余強度形變規(guī)律不符；理論模型1 則可較好地反映高孔隙率軟質(zhì)紅砂巖的峰后形變特征.以上現(xiàn)象進一步證明了理論模型1 的合理性，也證明了本文所建立的損傷本構(gòu)可適用于孔隙率較高的軟巖.

由圖2、3 可知：本文所建立的凍融巖石損傷本構(gòu)模型具有較廣泛的適用性，不僅適用于低孔隙率硬巖，也適用于高孔隙率軟巖；同時，本文所建立的模型還可較好地反映出不同類別巖石的峰后形變特征.

3 損傷變量演化規(guī)律

由式（5）、（11）、（14）及文獻[16-17]中的試驗數(shù)據(jù)計算得到兩類巖石在凍融循環(huán)下的總損傷演變曲線，如圖3 和圖4 所示.

圖3 低孔隙率硬質(zhì)花崗巖總損傷演變曲線Fig.3 Total damage evolution curves of hard granite with low-porosity

圖4 高孔隙率軟質(zhì)紅砂巖總損傷演變曲線Fig.4 Total damage evolution curves of soft red sandstone with high porosity

從圖3 中可以看出：不同凍融循環(huán)次數(shù)下花崗巖的總損傷演變曲線均呈“S”形；未受凍融影響的花崗巖初始階段總損傷近似為0，受凍融循環(huán)作用的花崗巖在初始階段總損傷為凍融損傷量，且初始階段的總損傷量隨凍融次數(shù)的增加而增加；隨著外力的加載，巖石內(nèi)部微裂紋和缺陷不斷發(fā)展擴張，受荷損傷開始快速增加，花崗巖的總損傷曲線也呈現(xiàn)出加速損傷的趨勢，巖石逐漸進入屈服階段；當巖石強度逐漸喪失后，總損傷曲線斜率開始變緩，并最終趨于穩(wěn)定值1，此時巖石已完全破壞，總損傷達到最大值.由圖4 可知：紅砂巖的總損傷曲線也呈現(xiàn)出近似“S”形的趨勢，其演變規(guī)律與花崗巖類似，不再展開描述.

基于以上分析，將低孔隙率硬質(zhì)花崗巖和高孔隙率軟質(zhì)紅砂巖總損傷演變曲線分為初始損傷階段、加速階段和完全損傷階段（圖5）：在初始損傷階段，巖石內(nèi)的總損傷量較低，且隨著應變的增加，總損傷量近乎不變；在加速損傷階段，總損傷量隨應變的增加而急劇增長，巖石內(nèi)產(chǎn)生大量的裂紋；在完全損傷階段，總損傷的增長趨勢逐漸變緩并最終達到穩(wěn)定，巖石在此階段發(fā)生宏觀破壞.

圖5 凍融受荷巖石總損傷演變規(guī)律Fig.5 The evolution of total damage of rock under loading and freeze-thaw action

4 模型參數(shù)分析

4.1 參數(shù)m 和 f0 對應力-應變曲線的影響

為確定本文所建立的損傷模型中參數(shù)m和f0的物理意義及對應力-應變曲線的影響，以凍融循環(huán)50 次，圍壓為10 MPa 的低孔隙率硬質(zhì)花崗巖以及凍融循環(huán)5 次，圍壓為4 MPa 的高孔隙率軟質(zhì)紅砂巖為例進行參數(shù)分析，如圖6、7 所示.

圖6 參數(shù) m 對應力-應變曲線的影響Fig.6 Influence of parameter m on the stress-strain curve

從圖6 中可以看出：隨著m的增加，低孔隙率硬質(zhì)花崗巖和高孔隙率軟質(zhì)紅砂巖的峰值強度均有所增加；在峰值應力后，應力跌落速率隨m的增大而增大，說明巖石的脆性特征逐漸增強.從圖7 中可以看出：隨著f0的增加，兩類巖石的應力-應變曲線峰值點升高且向右移動，達到峰值應力所需的應變量也隨之增加；巖石的塑性變形特征隨著f0的增加而增加.由此可以推斷出，參數(shù)m和f0的物理意義分別為巖石的脆性變形特征和塑性變形特征.

圖7 參數(shù)f0對應力-應變曲線的影響Fig.7 Influenceofparameterf0on the stress-strain curve

4.2 參數(shù) m 和 f0 對總損傷曲線的影響

由式（5）、（11）和式（14）計算得到總損傷量，并分析參數(shù)m和f0對總損傷量的影響，如圖8、9 所示.由圖8 可知：隨著參數(shù)m的增加，2 類巖石的初始損傷階段皆有所延長，即，巖石進入加速損傷階段的臨界應變值變大；加速損傷階段的斜率隨m的增加而增加，且?guī)r石進入完全損傷階段的臨界應變值隨m的增加而減小，說明了巖石提前進入完全損傷階段.由圖9 可知：隨著參數(shù)f0的增加，2 類巖石的初始損傷階段皆有所延長，這與參數(shù)m對總損傷曲線的影響規(guī)律一致；不同的是，隨著f0的增加，2 類巖石加速損傷階段的斜率均逐漸減小，進入完全損傷階段的臨界應變值增大，即巖石將逐漸滯后進入完全損傷階段.以上現(xiàn)象進一步表明參數(shù)m可反映巖石的脆性形變特征，而參數(shù)f0則代表巖石的塑性形變特征.

圖8 參數(shù) m 對總損傷量的影響Fig.8 Influence of parameter m on the total damage

圖9 參數(shù) f0 對總損傷量的影響Fig.9 Influence of parameter f0 on the total damage

5 工程算例分析

為分析凍融循環(huán)對大型工程穩(wěn)定性的影響，本節(jié)將以隧道工程為算例進行數(shù)值模擬分析.圖10 給出了工程算例的計算模型.如圖所示，該隧道斷面為三心拱形，拱頂埋深24.5 m，隧道圍巖為軟質(zhì)紅砂巖.該隧道采用全斷面開挖法，開挖循環(huán)進尺為2 m，計算開挖總長度為20 m.隧道支護方式為錨噴支護，噴射混凝土等級為C20，錨噴襯砌厚度0.15 m，錨桿截面面積為490.63 mm2，長度為4 m，錨桿間弧長1.8 m.

圖10 隧道工程算例計算模型（單位：m）Fig.10 Calculation model of tunnel engineering example （unit: m）

基于文獻[16]中所給出不同凍融次數(shù)下的紅砂巖力學參數(shù)，在大型巖土工程有限元軟件MIDAS GTS 中分別計算了凍融循環(huán)次數(shù)為0、5、10 次時，隧道開挖后的位移場和應力場分布.不同工況下隧道圍巖的力學參數(shù)見表3，在計算過程中假設(shè)圍巖的容重和內(nèi)摩擦不隨凍融次數(shù)而改變.計算得到的圍巖豎向位移場、圍巖最大主應力矢量場以及襯砌結(jié)構(gòu)最大主應力場分別如圖11～13 所示.

圖11 不同凍融次數(shù)下隧道圍巖豎向位移Fig.11 Vertical displacement of surrounding rock for tunnel under different freeze-thaw cycles

表3 不同凍融次數(shù)下隧道圍巖力學參數(shù)Tab.3 Mechanical parameters of surrounding rock for tunnel under different freeze-thaw cycles

從圖11 中可以看出：隧道開挖后拱頂和底部的圍巖垂直位移量較大；當圍巖經(jīng)歷0、5、10 次凍融循環(huán)后，隧道拱頂?shù)淖畲笙鲁亮糠謩e為3.75 、3.99 、4.42 mm，經(jīng)歷10 次凍融循環(huán)后，隧道拱頂?shù)膰鷰r下沉量增加了17.87%；類似地，當圍巖經(jīng)歷0、5、10 次凍融循環(huán)后，隧道底部的最大隆起量分別為4.73 、5.05 、5.64 mm，經(jīng)歷10 次凍融循環(huán)后，隧道底部的圍巖隆起量增加了19.24%.

從圖12 中可以看出：隧道開挖后拱頂與底部分布拉應力積聚區(qū)，隧道兩側(cè)圍巖分布壓應力積聚區(qū)；當圍巖經(jīng)歷0、5、10 次凍融循環(huán)后，隧道圍巖的最大拉應力極值分別為0.297 、0.300、0.305 MPa，經(jīng)歷10 次凍融循環(huán)后，隧道圍巖最大拉應力極值增加了2.7%；當圍巖經(jīng)歷0、5、10 次凍融循環(huán)后，隧道圍巖的最大壓應力極值分別為0.448、0.450、0.457 MPa，經(jīng)歷10 次凍融循環(huán)后，隧道圍巖的最大壓應力極值增加了2.01%.隧道圍巖的最大主應力隨凍融循環(huán)變化較小，與文獻[11]中得出的結(jié)論一致，這是由于圍巖內(nèi)主要是自重應力，自重應力在凍融循環(huán)下基本保持不變.因此，圍巖的最大主應力變化量很小，而由于巖石力學性質(zhì)的改變，圍巖位移隨凍融循環(huán)的變化量則較大.

此外，隧道襯砌結(jié)構(gòu)的應力分布也受到凍融作用的影響，從圖13 中可以看出：隧道襯砌拱頂和底部主要為拉應力積聚區(qū)，襯砌兩側(cè)分布壓應力積聚區(qū)，隨著凍融循環(huán)次數(shù)的增加，隧道襯砌拉、壓應力極值也逐漸增大；當圍巖經(jīng)歷0、5、10 次凍融循環(huán)后，隧道襯砌最大拉應力極值分別為4.88、5.25、5.93 MPa，經(jīng)歷10 次凍融循環(huán)后，隧襯砌最大拉應力極值增加了21.52%；當圍巖經(jīng)歷0、5、10 次凍融循環(huán)后，隧道襯砌最大壓應力極值分別為3.75、3.99、4.42 MPa，經(jīng)歷10 次凍融循環(huán)后，隧襯砌最大壓應力極值增加了17.87%.凍融作用后隧道圍巖的力學性質(zhì)劣化，導致圍巖壓力傳遞至襯砌結(jié)構(gòu)，引起襯砌應力的增加.由此可知，雖然凍融作用對隧道圍巖應力場影響不大，但隨著圍巖的軟化也可引起隧道結(jié)構(gòu)的破壞，影響隧道的整體穩(wěn)定性.因此，在寒區(qū)進行巖土工程建設(shè)時，應特別注意對凍融循環(huán)對巖土性質(zhì)的影響并采取適當?shù)姆乐未胧?

圖13 不同凍融次數(shù)下隧道襯砌最大主應力Fig.13 Maximum principal stress in tunnel lining under different freeze-thaw cycles

6 結(jié) 論

本文基于統(tǒng)計損傷力學理論，假設(shè)巖石微元體強度服從Weibull 分布，結(jié)合最大拉應變破壞準則建立了考慮凍融和荷載耦合作用的寒區(qū)凍融巖石損傷本構(gòu)模型；基于試驗數(shù)據(jù)及已有的損傷本構(gòu)模型驗證了所建立的損傷本構(gòu)模型的正確性，并對總損傷演變過程及模型參數(shù)進行了分析；最后，在數(shù)值模擬軟件中分析了凍融循環(huán)對隧道工程穩(wěn)定性的影響，得到的主要結(jié)論如下：

1） 以低孔隙率硬質(zhì)花崗巖和高孔隙率軟質(zhì)紅砂巖為驗證對象，結(jié)合試驗數(shù)據(jù)及已有的損傷本構(gòu)模型確定所建立模型的合理性和適用范圍，結(jié)果表明，本文所建立的模型可較好地重現(xiàn)2 類巖石的應力-應變?nèi)^程曲線，且能較好地反映巖石的峰后強度.

2） 不同凍融循環(huán)次數(shù)下低孔隙率硬質(zhì)花崗巖和高孔隙率軟質(zhì)紅砂巖的總損傷演變曲線均呈“S”形，且總損傷演變曲線可分為初始損傷階段、加速階段和完全損傷階段.

3） Weibull 分布參數(shù)m和f0分別代表巖石的脆性特征和塑性特征；隨著m增加，2 類巖石的應力峰值逐漸增加，峰后應力跌落速率逐漸增大，總損傷曲線加速損傷階段的斜率也逐漸增加，巖石提前進入完全損傷階段；隨著f0的增加，巖石達到峰值應力所對應的應變逐漸增加，總損傷曲線加速損傷階段斜率逐漸減小，巖石滯后進入完全損傷階段.

4） 采用數(shù)值模擬方法，計算了凍融循環(huán)對隧道圍巖位移、應力及隧道襯砌應力的影響，隨著凍融次數(shù)的增加，隧道圍巖的垂直位移、最大主應力及襯砌最大主應力均逐漸增加；凍融循環(huán)對隧道圍巖的應力分布影響較小，但由于圍巖物理力學性質(zhì)的劣化，圍巖壓力傳遞至襯砌可導致隧道結(jié)構(gòu)產(chǎn)生損傷.