凍融循環(huán)作用對(duì)土質(zhì)邊坡穩(wěn)定性的影響

宋彥琦， 鄭俊杰， 李向上， 李小龍， 邵志鑫

(中國(guó)礦業(yè)大學(xué)(北京)力學(xué)與建筑工程學(xué)院，北京 100083)

隨著季節(jié)的交替，凍土區(qū)會(huì)反復(fù)發(fā)生凍融循環(huán)作用，在此過(guò)程中，土壤理化性質(zhì)的改變使邊坡穩(wěn)定狀態(tài)逐漸遭受破壞，最終發(fā)生失穩(wěn)。據(jù)調(diào)查，中國(guó)高海拔地區(qū)多年凍土面積居世界首位[1]，因此研究?jī)鋈谘h(huán)作用對(duì)邊坡穩(wěn)定性的影響具有重要且現(xiàn)實(shí)的意義。

中外對(duì)凍融循環(huán)作用下邊坡穩(wěn)定性的研究取得了大量成果。實(shí)驗(yàn)方面，常丹等[2]對(duì)青藏粉砂土進(jìn)行了常規(guī)三軸實(shí)驗(yàn)，結(jié)果表明隨著凍融循環(huán)次數(shù)的增加，土樣的彈性模量和靜強(qiáng)度先減小后增大，黏聚力逐漸減小，內(nèi)摩擦角先減小后增加。楊更社等[3]對(duì)白堊系飽和砂巖進(jìn)行了實(shí)驗(yàn)研究，發(fā)現(xiàn)凍結(jié)作用對(duì)飽和砂巖的力學(xué)性質(zhì)的影響較大。崔廣芹等[4]研究了內(nèi)蒙古包頭凍土地區(qū)的粉砂土，發(fā)現(xiàn)對(duì)于不同干密度和不同含水率的土樣，隨著凍融次數(shù)的增加，黏聚力會(huì)減小，內(nèi)摩擦角變化不大，黏聚力對(duì)抗剪強(qiáng)度影響較大。葛琪[5]對(duì)凍融界面的抗剪強(qiáng)度進(jìn)行了研究，發(fā)現(xiàn)隨著凍融循環(huán)次數(shù)和含水量的增加，凍融界面土體的抗剪強(qiáng)度逐漸降低，且含水量的影響更大，但其影響隨凍融次數(shù)增加而遞減并最終穩(wěn)定。梁波等[6]對(duì)青藏鐵路清水河、北麓河的地基土進(jìn)行取樣進(jìn)行室內(nèi)模擬試驗(yàn)，發(fā)現(xiàn)隨著凍融循環(huán)次數(shù)的增加，土體密度逐漸趨于某一定值，融沉總變形在增加的同時(shí)呈現(xiàn)指數(shù)遞減的趨勢(shì)。通過(guò)這些實(shí)驗(yàn)可以看出，凍融循環(huán)作用是影響巖土材料強(qiáng)度的一個(gè)重要因素。數(shù)值模擬方面，鄭穎人等[7-8]和趙尚毅等[9]通過(guò)有限元軟件利用強(qiáng)度折減法對(duì)邊坡安全系數(shù)進(jìn)行求解，發(fā)現(xiàn)有限元強(qiáng)度折減法完全可以用來(lái)解決工程實(shí)際問(wèn)題，并對(duì)有限元強(qiáng)度折減法的計(jì)算精度和影響因素進(jìn)行了探討和完善。張永明等[10]和陳印東等[11]通過(guò)強(qiáng)度折減法利用有限元軟件對(duì)邊坡進(jìn)行模擬，其計(jì)算結(jié)果與極限平衡法基本一致，表明強(qiáng)度折減法在邊坡穩(wěn)定性的判定上具有合理性與可行性。史俊濤等[12]利用容重增加法和強(qiáng)度折減法對(duì)某土質(zhì)邊坡進(jìn)行穩(wěn)定性分析，發(fā)現(xiàn)強(qiáng)度折減法相比容重增加法使用更加合理。

綜合上述研究可以看出，雖然中外在凍融循環(huán)作用對(duì)巖土力學(xué)性能的實(shí)驗(yàn)研究以及邊坡穩(wěn)定性的數(shù)值模擬研究已經(jīng)取得了一定的成果，但針對(duì)凍融循環(huán)這個(gè)影響因素對(duì)邊坡整體穩(wěn)定性的數(shù)值模擬研究很少，且強(qiáng)度折減法對(duì)于邊坡穩(wěn)定性的評(píng)價(jià)具有合理性和實(shí)用性。本文利用強(qiáng)度折減法，通過(guò)大型有限元軟件ABAQUS計(jì)算獲得邊坡的穩(wěn)定安全系數(shù)，與簡(jiǎn)化Bishop法獲得的理論安全系數(shù)進(jìn)行比較，并以特征部位位移發(fā)生突變作為破壞判據(jù)，得出了邊坡安全系數(shù)隨坡角的減小而增大，隨凍融循環(huán)次數(shù)增加而逐漸降低最終趨于平穩(wěn)的變化趨勢(shì)，該文可為邊坡穩(wěn)定性影響因素的探究提供參考。

1 邊坡穩(wěn)定性的數(shù)值模擬

1.1 強(qiáng)度折減法的基本原理

強(qiáng)度折減法是在理想彈塑性的有限元計(jì)算中，逐漸降低邊坡巖土體的抗剪強(qiáng)度參數(shù)，直至達(dá)到極限破壞狀態(tài)，同時(shí)得到邊坡的強(qiáng)度儲(chǔ)備安全系數(shù)Fs[13]。

對(duì)于摩爾-庫(kù)侖材料，關(guān)于剪應(yīng)力τ有：

(1)

所以有：

(2)

(3)

式中，σ為法向應(yīng)力，c和φ為有效黏聚力和有效內(nèi)摩擦角；c′和φ′為折減后的有效黏聚力和有效內(nèi)摩擦角。

目前有限元中邊坡整體失穩(wěn)的判據(jù)主要有三種[8]：①以有限元計(jì)算是否收斂作為破壞判據(jù)；②以特征部位位移發(fā)生突變作為破壞判據(jù)；③以是否形成連續(xù)的貫通區(qū)作為破壞判據(jù)。

1.2 計(jì)算模型和材料參數(shù)

某邊坡模型尺寸如圖1所示，土樣的力學(xué)參數(shù)采用文獻(xiàn)[14]中尾礦土的數(shù)據(jù)，其中土體密度ρ=1.492 g·cm-3，泊松比取v=0.35，彈性模量E采取圍壓σ3=0.2 MPa時(shí)測(cè)得的數(shù)值，不同凍融循環(huán)次數(shù)下土體物理力學(xué)參數(shù)如表1所示。

圖1 模型尺寸Fig.1 Model size

表1 不同凍融循環(huán)次數(shù)下土體物理力學(xué)參數(shù)

按照二維平面應(yīng)變建立有限元模型，計(jì)算采用Mohr-Coulomb準(zhǔn)則，模型邊界條件為：左右兩側(cè)限制水平方向的位移，底部限制水平和垂直兩個(gè)方向的位移，上部為自由邊界。模型只受重力荷載作用，整個(gè)模型包括349個(gè)單元，單位類型為CPE4，網(wǎng)格劃分如圖2所示。

圖2 有限元網(wǎng)格劃分Fig.2 Finite element mesh

1.3 計(jì)算結(jié)果

1.3.1 塑性應(yīng)變分析

塑性應(yīng)變等值云圖反映了邊坡塑性應(yīng)變的分布，圖3所示為1次凍融循環(huán)作用下不同折減系數(shù)時(shí)三種不同情況下邊坡的塑性應(yīng)變等值云圖，邊坡變形顯示比例設(shè)置為0。

圖3 不同折減系數(shù)下的塑性區(qū)Fig.3 Continuous contour of plastic strain under different reduction factors

由圖3塑性應(yīng)變的分布可直觀判斷出邊坡失穩(wěn)的發(fā)展過(guò)程，即塑性區(qū)最先從坡腳處產(chǎn)生，隨著邊坡發(fā)生屈服的區(qū)域逐漸擴(kuò)大并向坡頂延伸，直至到達(dá)坡頂，即坡腳與坡頂間的塑性區(qū)發(fā)生貫通，邊坡發(fā)生失穩(wěn)。

1.3.2 滑動(dòng)面確定

滑動(dòng)面的位置與形狀可以通過(guò)位移等值線分布圖來(lái)確定。

圖4為1次凍融循環(huán)作用下有限元計(jì)算不收斂時(shí)的位移等值云圖。圖4清楚地定位了滑動(dòng)面的位置，形狀呈大致的圓弧狀。

圖4 用位移等值云圖表示的滑動(dòng)面位置和形狀Fig.4 Position and shape of the sliding surface represented by the displacement equivalent cloud

1.3.3 安全系數(shù)分析

在計(jì)算過(guò)程中，將場(chǎng)變量取為強(qiáng)度折減系數(shù)，在0.5～7變化，線性增加場(chǎng)變量的大小，增量為0.25，特征點(diǎn)取為邊坡坡面左上角頂點(diǎn)，見圖1紅點(diǎn)。圖5為1次凍融循環(huán)作用下特征點(diǎn)處折減系數(shù)FV隨水平位移U1的變化關(guān)系圖。

圖5 折減系數(shù)隨位移的變化關(guān)系Fig.5 Variation of the reduction factor with displacement

如圖5所示，A(-0.001 3,4.954 2)處特征點(diǎn)位移發(fā)生突變，B(-0.629 2,5.044 3)處計(jì)算不收斂。折減系數(shù)從0.5～4.954 2這個(gè)階段，特征點(diǎn)位移幾乎一直為0，說(shuō)明這個(gè)階段坡頂面幾乎沒(méi)有產(chǎn)生任何變形，由圖3(a)～圖3(c)可看出當(dāng)折減系數(shù)未達(dá)到4.954 2時(shí)塑性區(qū)并未蔓延至特征點(diǎn)，因此在這個(gè)階段，邊坡尚處于穩(wěn)定階段。直至折減系數(shù)大于4.954 2之后，水平位移絕對(duì)值迅速增加，說(shuō)明這個(gè)階段特征點(diǎn)周圍土體受到了很大的影響，導(dǎo)致坡頂面迅速產(chǎn)生大量變形，位移為負(fù)值說(shuō)明該部位處于向下滑動(dòng)的變形狀態(tài)，由圖3(d)可看出當(dāng)折減系數(shù)達(dá)到4.954 2時(shí)塑性區(qū)延伸至特征點(diǎn)，邊坡產(chǎn)生了一個(gè)完整的塑性貫通區(qū)，當(dāng)位移達(dá)到一定值時(shí)，邊坡發(fā)生失穩(wěn)。

1.3.4 安全系數(shù)確定

建立相同的邊坡計(jì)算模型，如圖1和圖2所示，分別采用凍融循環(huán)次數(shù)為0、1、5、10、15次時(shí)土體的力學(xué)參數(shù)，如表1所示，計(jì)算邊坡的安全系數(shù)，選取不同判據(jù)得出的安全系數(shù)如表2所示。其中，安全系數(shù)F1、F2、F3分別為以特征部位位移發(fā)生突變作為破壞判據(jù)、以是否形成連續(xù)的貫通區(qū)作為破壞判據(jù)、以有限元計(jì)算是否收斂作為破壞判據(jù)所得出的結(jié)果。

表2 不同凍融循環(huán)次數(shù)下不同破壞判據(jù)求出的邊坡安全系數(shù)

由表2可看出采用有限元中邊坡整體失穩(wěn)的三種判據(jù)得到的安全系數(shù)相差不大，其中以特征部位位移發(fā)生突變作為破壞判據(jù)與以是否形成連續(xù)的貫通區(qū)為破壞判據(jù)得出的安全系數(shù)更為接近。

2 邊坡穩(wěn)定性的理論分析

2.1 邊坡穩(wěn)定分析的簡(jiǎn)化Bishop法

采用文獻(xiàn)[15]中簡(jiǎn)化Bishop的積分表達(dá)求解邊坡的安全系數(shù)。簡(jiǎn)化Bishop的積分表達(dá)具體如下：邊坡高度為H，坡角為β，邊坡危險(xiǎn)滑移線為圓弧，以坡角為坐標(biāo)原點(diǎn)，水平方向?yàn)閤軸，垂直方向?yàn)閥軸，滑移線的圓心為O(a,b)(a，b為圓心在x，y軸的坐標(biāo))，半徑為R，將邊坡滑體部分分為若干條塊，條塊寬度為積分變量dx，條塊高度h的表達(dá)式為

(4)

式(4)中：

邊坡安全系數(shù)F為

(5)

式(5)中：γ為滑體的容重，α為內(nèi)摩擦角，c為有效黏聚力。

2.2 滑動(dòng)圓弧圓心范圍的確定

簡(jiǎn)單土質(zhì)邊坡最危險(xiǎn)圓弧圓心的范圍，可按潘家錚[16]的方法確定，該方法如圖6所示。

圖6 滑動(dòng)圓弧圓心范圍的確定Fig.6 Determination of the center of the sliding arc center

邊坡坡度為1：m，高為H，長(zhǎng)為L(zhǎng)=mH，作出邊坡的中垂線和中法線，以邊坡中點(diǎn)O為圓心，以R1=L/2和R2=3L/4為半徑作圓弧，在邊坡上側(cè)分別與中垂線和中法線圍成一個(gè)四邊形，最危險(xiǎn)圓心常常位于這個(gè)四邊形范圍內(nèi)。

在搜索最危險(xiǎn)圓心時(shí)，為方便編程，將圓心范圍擴(kuò)大為四邊形的外接矩形CDEF，如圖6所示。

2.3 計(jì)算結(jié)果

通過(guò)MATLAB編程實(shí)現(xiàn)邊坡的安全系數(shù)的求解，其中計(jì)算步長(zhǎng)設(shè)置為0.1 m，計(jì)算結(jié)果如表3所示。

表3 不同凍融循環(huán)次數(shù)下邊坡安全系數(shù)的理論解

3 結(jié)果分析

3.1 不同破壞判據(jù)與理論解的對(duì)比分析

表4為凍融循環(huán)次數(shù)分別為0、1、5、10、15次時(shí)，邊坡采取強(qiáng)度折減法選取不同破壞判據(jù)和簡(jiǎn)化Bishop法得出的理論解的誤差對(duì)比，不同破壞判據(jù)結(jié)果的誤差取5次凍融循環(huán)情況的誤差平均值。其中破壞判據(jù)一為特征部位位移發(fā)生突變，破壞判據(jù)

表4 不同破壞判據(jù)求出的安全系數(shù)與理論解的誤差對(duì)比

二為形成連續(xù)的貫通區(qū)，破壞判據(jù)三為有限元計(jì)算不收斂。

由表4可看出，分別采用三種破壞判據(jù)得出的不同凍融循環(huán)條件下邊坡的安全系數(shù)與簡(jiǎn)化Bishop法得到的理論值的誤差都很小，即通過(guò)有限元模擬與簡(jiǎn)化Bishop法得到的安全系數(shù)基本一致，說(shuō)明這兩種方法在本例中都是可行的。

通過(guò)判據(jù)一得出的安全系數(shù)與理論解的誤差和通過(guò)判據(jù)二得出的安全系數(shù)與理論解的誤差相差甚小，且都小于通過(guò)判據(jù)三得出的安全系數(shù)與理論解的誤差，即以特征部位位移發(fā)生突變和是否形成連續(xù)的貫通區(qū)作為破壞判據(jù)相比以有限元是否收斂作為破壞判據(jù)更為準(zhǔn)確，分析其原因是當(dāng)塑形變形區(qū)從坡腳貫通到坡頂時(shí)，若無(wú)其他約束條件比如抗滑樁、錨桿等邊坡支護(hù)設(shè)施，塑性區(qū)上方部位會(huì)迅速產(chǎn)生大量變形，位移自然會(huì)大幅度增加，但此時(shí)由于載荷步步長(zhǎng)設(shè)置或其他原因?qū)е掠邢拊?jì)算不一定不收斂，即有限元計(jì)算收斂時(shí)不一定表明邊坡處于穩(wěn)定狀態(tài)。因此在實(shí)際工程應(yīng)用中，可根據(jù)實(shí)際情況，將特征部位位移是否發(fā)生突變和塑性區(qū)是否貫通兩個(gè)判據(jù)綜合考慮判斷邊坡的穩(wěn)定性狀態(tài)。

在接下來(lái)對(duì)邊坡安全系數(shù)的求解過(guò)程中，通過(guò)有限元強(qiáng)度折減法進(jìn)行計(jì)算且全部以特征部位位移發(fā)生突變?yōu)槠茐呐袚?jù)確定安全系數(shù)數(shù)值。

3.2 凍融循環(huán)次數(shù)對(duì)邊坡穩(wěn)定性的影響

為更直觀地分析得到凍融循環(huán)次數(shù)對(duì)邊坡安全系數(shù)的影響，將表2中以特征部位位移發(fā)生突變?yōu)槠茐呐袚?jù)得到的安全系數(shù)繪成折線圖，如圖7所示。

圖7 安全系數(shù)隨凍融循環(huán)次數(shù)的變化關(guān)系Fig.7 Variation of safety factor with the number of freeze-thaw cycles

由圖7可看出，在本例中，當(dāng)坡比固定時(shí)，凍融循環(huán)作用后，邊坡穩(wěn)定性逐漸下降，從10次凍融循環(huán)到15次凍融循環(huán)，有限元計(jì)算的邊坡的安全系數(shù)下降比僅為7.494 3%，而經(jīng)1次凍融循環(huán)作用后，安全系數(shù)下降比達(dá)24.959 1%，下降幅度最大。

3.3 邊坡坡比對(duì)邊坡穩(wěn)定性的影響

為研究邊坡坡比對(duì)邊坡穩(wěn)定性的影響，分別選取坡比為1∶0.5、1∶0.75、1∶1、1∶1.25建立不同坡比的邊坡模型，如圖8所示，計(jì)算相應(yīng)邊坡在凍融循環(huán)次數(shù)分別為0、1、5、10、15次時(shí)的穩(wěn)定系數(shù)，結(jié)果如圖9所示。

圖8 不同坡比邊坡尺寸示意Fig.8 Slope size at different slope ratios

圖9 安全系數(shù)隨坡比的變化關(guān)系Fig.9 Variation of safety factor with slope ratio

由圖9可清楚地看出，當(dāng)坡比一定時(shí)，隨著凍融循環(huán)次數(shù)的增加，四種不同坡比情況下的邊坡的安全系數(shù)下降的幅度都會(huì)逐漸減小，并最終趨于穩(wěn)定；當(dāng)凍融次數(shù)一定時(shí)，四種不同坡比情況下的邊坡的安全系數(shù)都會(huì)隨著坡角的減小逐漸增大，即邊坡越平緩，安全系數(shù)越大，邊坡穩(wěn)定性越好。

4 結(jié)論

(1)通過(guò)強(qiáng)度折減法，采用有限元中邊坡整體失穩(wěn)的三種判據(jù)求出一定坡比下不同凍融循環(huán)次數(shù)下的安全系數(shù)，三種判據(jù)得到的結(jié)果相差不大。將其與簡(jiǎn)化Bishop法得到的理論安全系數(shù)進(jìn)行對(duì)比，兩者計(jì)算結(jié)果基本一致。

(2)以特征部位位移是否發(fā)生突變和是否形成連續(xù)的貫通區(qū)作為破壞判據(jù)相比以有限元是否收斂作為破壞判據(jù)更為準(zhǔn)確，在實(shí)際工程應(yīng)用中，可根據(jù)實(shí)際情況，將特征部位位移是否發(fā)生突變和塑性區(qū)是否貫通兩個(gè)判據(jù)綜合考慮判斷邊坡的穩(wěn)定性狀態(tài)。

(3)通過(guò)強(qiáng)度折減有限元方法建立ABAQUS力學(xué)模型，研究邊坡的穩(wěn)定性，相比簡(jiǎn)化Bishop法而言，更加客觀簡(jiǎn)便，不需要假定危險(xiǎn)滑動(dòng)面，就能夠直接動(dòng)態(tài)地展示不同邊坡的失穩(wěn)破壞過(guò)程，該方法可為邊坡穩(wěn)定性影響因素的探究提供參考。

(4)討論了凍融循環(huán)以及邊坡波比對(duì)邊坡穩(wěn)定性的影響，但是邊坡研究還應(yīng)考慮現(xiàn)場(chǎng)的環(huán)境溫度，濕度等因素的影響，另外力學(xué)模型的簡(jiǎn)化程度等因素對(duì)于模擬結(jié)果都是有影響的，這些都是今后研究中需要進(jìn)一步考慮的問(wèn)題。