考慮分區(qū)各向異性和滲流作用的邊坡穩(wěn)定性研究

張幫鑫 ,賈劍青 ,賴遠明 ,王宏圖 ,辛成平

（1.蘭州交通大學(xué) 交通運輸學(xué)院，蘭州 730070；2.中國科學(xué)院 西北生態(tài)環(huán)境資源研究院，蘭州 730000；3.重慶大學(xué) 資源與安全學(xué)院，重慶 400044）

土體強度的各向異性和滲流作用是影響邊坡穩(wěn)定性的重要因素[1-2]；強度折減法是邊坡穩(wěn)定性分析的重要方法之一[3]，近年來，學(xué)者們對此開展了一系列研究。蘇永華等[1]建立了考慮土體含水率分布與坡面滲流作用的降雨入滲分析模型，并確定了邊坡穩(wěn)定性系數(shù)計算公式；Sun 等[3]提出了一種基于特征點位移突變的殘差位移增量準(zhǔn)則，將最大平均殘余位移增量對應(yīng)的強度折減系數(shù)作為邊坡的安全系數(shù)；Rao 等[4]研究了非均質(zhì)和各向異性黏性土邊坡的抗滑樁加固技術(shù)和穩(wěn)定性；Xu 等[5]探討了非均質(zhì)系數(shù)和各向異性系數(shù)與邊坡體加固效果的關(guān)系；羅凌暉等[6]研究了土體強度各向異性對成層土邊坡穩(wěn)定性的影響，并結(jié)合強度折減法分析了其穩(wěn)定性；曾鈴等[7]試驗研究了不同降雨條件下砂土和粉質(zhì)黏土體積含水率與基質(zhì)吸力的變化規(guī)律；Zeng 等[8]研究了降雨入滲條件下，邊坡體內(nèi)瞬態(tài)飽和區(qū)的形成、發(fā)展和消散過程；何曉瑩[9]考慮滲流正交異性的耦合作用，研究了膠凝砂礫石壩的滲流及應(yīng)力變化特性；王正成等[10]考慮彈性模量、泊松比等參數(shù)的變化，建立了流固耦合模型，分析了某壩基滲流場與應(yīng)力場特征；史卜濤等[11]提出了物質(zhì)點強度折減法，并分析了某邊坡穩(wěn)定性；Nie 等[12]研究了強度折減法的破壞準(zhǔn)則，并與極限平衡法進行了比較；劉亞棟[13]分區(qū)考慮土體強度各向異性，計算分析了邊坡穩(wěn)定性，但該分區(qū)方法僅視覺性地將應(yīng)力偏轉(zhuǎn)方向大體一致的區(qū)域進行集中分區(qū)。目前，鮮有文獻考慮土體分區(qū)各向異性和滲流作用對邊坡穩(wěn)定性的影響。筆者以某邊坡工程為例，考慮土體分區(qū)各向異性和雨水的入滲作用，建立流固耦合模型，并采用強度折減法計算分析邊坡穩(wěn)定性。

1 流固耦合理論與強度折減法

1.1 流固耦合理論

流固耦合的基本思想是采用有限元方法將邊坡體離散為有限個土體單元，并通過求解由達西理論和有限元理論所建立的有限元方程，即可得到單元體的節(jié)點位移和孔隙水壓力，然后結(jié)合滲流場和應(yīng)力場相互作用機理，推出兩場耦合的有限元模型[14-15]。流固耦合中，滲流場影響應(yīng)力場模型可表達為[14-15]

式中：k為滲透系數(shù)；σ和p為k的影響參數(shù)；H為形函數(shù)；q為 流量；n2為邊界Γ2的法線 方向；n3為邊界Γ3的法線方向。

應(yīng)力場影響滲流場模型為[14-15]

式中：L為微分算子矩陣；n為面力邊界Sσ法線方向的方向余弦矩陣，即為滲流場的水頭分布函數(shù)；X為節(jié)點外荷載矩陣；{u}為位移場；B為應(yīng)力插值矩陣；D為彈性矩陣；Sσ為已知面力邊界；tˉ為已知面力邊界Sσ上的面力分布，也為水頭分布H（x,y）的函數(shù)；SM為位移邊界；{uˉ}為位移邊界SM上的位移矢量。

則流固耦合有限元模型為[14-16]

式中：K為滲透系數(shù)相關(guān)矩陣；{f}為滲流場水頭分布函數(shù)；{H}為形函數(shù)；M為整體剛度矩陣；X為節(jié)點外荷載矩陣；F為滲透力矩陣。

1.2 強度折減法

強度折減法是在外荷載不變的情況下，將土體抗剪強度參數(shù)的黏聚力和內(nèi)摩擦角進行折減，得到新的土體強度參數(shù)，并通過不斷增大折減系數(shù)Fr直至邊坡達到極限平衡狀態(tài)，此時折減系數(shù)即為邊坡體的安全系數(shù)，可表示為[15,17]

2 工程背景

研究的邊坡位于陜西省咸陽市旬邑縣境內(nèi)，為中國國家高速公路菏澤至寶雞聯(lián)絡(luò)線某合同段項目；邊坡體地質(zhì)構(gòu)造為發(fā)育中、新生代陸相地層，地下水主要由降雨及周邊河谷暗流組成；項目區(qū)地處內(nèi)陸，年平均降雨量為589.4 mm。邊坡橫斷面如圖1 所示（為便于分析，下文圖表中①為水泥土；②為素土；③為淺黃色粉土；④為黃褐色粉質(zhì)黏土）。各土層物理力學(xué)參數(shù)如表1 所示。

表1 物理力學(xué)參數(shù)Table 1 The physical and mechanical parameters

圖1 邊坡橫斷面圖（單位：m）Fig.1 Cross-sectional view of the slope (Unit: m)

3 邊坡穩(wěn)定性分析

3.1 模型建立和參數(shù)確定

3.1.1 模型建立 建立的邊坡有限元模型如圖2所示。根據(jù)邊坡成層土體地質(zhì)分布特征并結(jié)合其幾何形狀，將邊坡體劃分為I 區(qū)、Ⅱ區(qū)、Ⅲ區(qū)和Ⅳ區(qū)4個區(qū)域。

圖2 有限元模型Fig.2 Finite element model

流固耦合分析中，在模型gh邊施加水平和垂直位移約束，fg和ah邊施加水平位移約束，其他邊均為自由邊界；考慮降雨入滲的影響，在abcdef邊施加降雨荷載，根據(jù)該地區(qū)年平均降雨量的統(tǒng)計結(jié)果，取雨水入滲強度為6.365×10-5m/h；在fghi區(qū)域內(nèi)考慮孔隙水壓力變化；模型網(wǎng)格單元類型選用孔隙流體單元（CPE4P）。

3.1.2 參數(shù)確定 土層抗剪強度參數(shù)會隨土體最大主應(yīng)力方向角α的變化而變化[14,18-19]。當(dāng)角α分別為0°、30°、60°和90°時，取土體各向異性系數(shù)K=0.7，采用Casagrande 法[13]計算所得各土層的強度參數(shù)，如表2 所示。利用最小二乘法并結(jié)合表2，推導(dǎo)土體①、②、③、④各向異性強度參數(shù)方程，分別為式（5）、式（6）、式（7）、式（8）。

表2 抗剪強度參數(shù)變化Table 2 Variation of shear strength parameters

式中：Ac=-0.068 05×(1-3 cos2α) ；Aφ=0.061 32× (1-3 cos2α)，Ac和Aφ分別為黏聚力和內(nèi)摩擦角的各向異性狀態(tài)參數(shù)；b為中間主應(yīng)力參數(shù)[20-21]；α為最大主應(yīng)力方向角，對于平面應(yīng)變問題，角α可通過式（9）確定[13]，即

式 中：σx、σy和τxy分別為該單元在xy平 面內(nèi)水平和垂直方向的正應(yīng)力和剪應(yīng)力；r表示應(yīng)力Mohr 圓半徑。

滲透系數(shù)隨基質(zhì)吸力的變化為[14-15,22]

式中：Kws為土體飽和時的滲透系數(shù)，取Kws=5×10-6m/s，即Kws=0.018 m/h；uw為土體中的水壓；aw、bw和cw為材料系數(shù)，取aw=1 000，bw=0.01，cw=1.7。

飽和度隨基質(zhì)吸力的變化為[14-15]

式中：Sr為飽和度；Si為殘余飽和度，取Si=0.08；Sn為最大飽和度，取Sn=1；as、bs和cs為材料系數(shù)，分別取as=1、bs=1.5×10-5、cs=3.5。

3.2 分級加載及分區(qū)方法

邊坡初始應(yīng)力狀態(tài)對模型單元體最大主應(yīng)力、最小主應(yīng)力及剪應(yīng)力等均產(chǎn)生一定影響；隨著自重應(yīng)力的逐漸增加，邊坡體一定范圍內(nèi)的土體逐步進入塑性狀態(tài)并產(chǎn)生一定的塑性變形。為考慮塑性應(yīng)變在坡體局部累積過程中對最大主應(yīng)力方向角α的影響，通過分級加載方式將各層土體的重力施加于模型。通過分級加載，實現(xiàn)邊坡體不同區(qū)域二級分區(qū)的動態(tài)調(diào)整；與此同時，計算各二級分區(qū)各向異性強度參數(shù)并對其進行賦值。該加載共有3 級，各級荷載施加值如表3 所示。

表3 分級荷載施加值Table 3 Load application values for graded loads

為使土體強度的各向異性更切合實際，在分級荷載施加過程中采用了分區(qū)技術(shù)。分區(qū)采用初分、細分和精分3 級控制。1）初分，每級荷載施加完成后，計算所有土體單元角α的大小，并以0°～30°、30°～60°和60°～90°為劃分標(biāo)準(zhǔn)，對I、Ⅱ、Ⅲ和Ⅳ 4 個大區(qū)進行初次區(qū)域劃分。2）細分，初次區(qū)域劃分后，由于同一區(qū)域內(nèi)可能包含不同性質(zhì)的土層，為保證后續(xù)劃分結(jié)果的準(zhǔn)確性，對包含不同性質(zhì)土層的初次劃分區(qū)域進行二次細分，即以土層界面為分界線，將該分區(qū)劃分為與土層數(shù)相等的若干二級分區(qū)。3）精分，細分完成后，結(jié)合數(shù)值模擬結(jié)果并利用式（5）～式（9）各向異性強度參數(shù)計算公式，計算3 次分級加載下各區(qū)域內(nèi)土體的c、φ值，與此同時，對細分后發(fā)生變化的區(qū)域即時進行調(diào)整。第1級加載時，取α=0°時對應(yīng)的各土層強度參數(shù)。各級加載的分區(qū)演化如圖3 所示，各參數(shù)的演變?nèi)鐖D4所示。

圖3 各級加載的分區(qū)演化Fig.3 Partition evolution for each level of loading

圖4 強度參數(shù)調(diào)整Fig.4 Strength parameter adjustment

圖3 表明：在3 次荷載施加過程中，I、Ⅱ、Ⅲ和Ⅳ 4 個初始分區(qū)的二級分區(qū)數(shù)量沒有增加或減少，但二級分區(qū)的分布域發(fā)生了一定變化，且分區(qū)演化范圍逐漸減??；分區(qū)演化主要發(fā)生在α為30°～60°的范圍內(nèi)。第2 次加載后，I2、I3、Ⅱ2、Ⅱ3和Ⅳ2五個二級分區(qū)范圍有所擴大，如圖3（b）；第3 次加載后，Ⅱ2和Ⅱ3范圍進一步擴大，如圖3（c）。從圖4 可以看出，各分區(qū)α平均值和強度參數(shù)的變化范圍隨著荷載的逐級施加而逐漸減小。

3.3 基于流固耦合的邊坡穩(wěn)定性分析

建立考慮土體分區(qū)各向異性（工況1）、未分區(qū)各向異性（工況2）和未分區(qū)各向同性（工況3）3 種工況的流固耦合模型，模擬分析了邊坡應(yīng)力場、位移場和滲流場的變化規(guī)律。模擬所得各工況平均應(yīng)力云圖如圖5 所示，塑性應(yīng)變云圖如圖6 所示，流速 云圖如圖7 所示。

圖5 平均應(yīng)力云圖Fig.5 Average stress cloud diagram

圖6 塑性應(yīng)變云圖Fig.6 Plastic strain cloud diagram

圖7 流速云圖Fig.7 Flow velocity cloud diagram

由圖5～圖7 可知：工況1 的平均應(yīng)力最大值較工況2 下降了6.4%，較工況3 上升了15.8%。工況1、2、3 的邊坡滑動帶寬依次增大，平均滑動帶寬分別為2、3、5 m；3 種工況的坡體滑動深度、滲流域和流速大小[22-23]逐漸增大，其中，滑動深度分別為18.0、19.0、19.7 m。模擬過程中還發(fā)現(xiàn)：3 種工況的邊坡塑性區(qū)貫通時間逐漸延長；工況1 在第10 個折減增量步時塑性區(qū)完全貫通，而工況2、3 分別在第18、第29 個折減增量步時塑性區(qū)才完全貫通，說明工況1的邊坡滑動時間先于工況2、工況3。

為確定邊坡的穩(wěn)定性狀態(tài)，以二級臺階坡肩頂點d點為追蹤研究點（如圖2 所示），分別以有限元計算結(jié)果是否收斂（標(biāo)準(zhǔn)A）、特征點位移出現(xiàn)拐點（標(biāo)準(zhǔn)B）和是否形成連續(xù)的貫通區(qū)（標(biāo)準(zhǔn)C）作為評價標(biāo)準(zhǔn)[4,24]，計算所得3 種工況的邊坡穩(wěn)定性系數(shù)變化曲線如圖8 所示（圖中A、B和C為判斷標(biāo)準(zhǔn)，1、2和3 表示工況）。

圖8 穩(wěn)定性系數(shù)變化曲線Fig.8 Stability coefficient variation curve

研究表明，標(biāo)準(zhǔn)A和標(biāo)準(zhǔn)B均可作為邊坡失穩(wěn)的判斷依據(jù)，而標(biāo)準(zhǔn)C為邊坡失穩(wěn)的必要非充分條件[3]；與其他兩種標(biāo)準(zhǔn)相比，標(biāo)準(zhǔn)B與采用極限平衡法計算所得邊坡穩(wěn)定性系數(shù)的誤差最小[14,24]。由圖8 可以看出：以標(biāo)準(zhǔn)B為判斷標(biāo)準(zhǔn)計算所得邊坡的穩(wěn)定性系數(shù)相對最小，3 種工況的邊坡穩(wěn)定性系數(shù)分別為1.109、1.141 和1.409；工況1 的穩(wěn)定性系數(shù)相比于工況2、工況3 分別下降了2.8%和21.3%。

4 結(jié)論

1）方向角α的大小沿坡體外自由邊界向坡體內(nèi)部逐漸增大，且隨著荷載的逐級施加，該變化過程逐漸趨于穩(wěn)定。

2）考慮分區(qū)土體強度各向異性時，與其他兩種工況相比，邊坡體塑性區(qū)貫通時間最短，滑動帶深度最淺，滲流域和流速最小。

3）采用3 種判斷標(biāo)準(zhǔn)確定的邊坡穩(wěn)定性系數(shù)中，考慮分區(qū)土體強度各向異性時的邊坡穩(wěn)定性系數(shù)值最小，其次為未考慮分區(qū)土體強度各異性時，而未分區(qū)各向同性時的邊坡穩(wěn)定性系數(shù)值最大。