庫水下降聯(lián)合降雨作用下樹坪滑坡流固耦合分析

，世梅,

(三峽大學(xué) 三峽庫區(qū)地質(zhì)災(zāi)害教育部重點(diǎn)實驗室，湖北 宜昌 443002)

1 研究背景

三峽水庫蓄水以后，沿岸地區(qū)的自然條件發(fā)生顯著性變化,一方面庫水上升引起滑坡體地下水位抬升，從而改變滑坡體內(nèi)巖土體的應(yīng)力狀態(tài)，減小滑坡體的有效應(yīng)力；另一方面，庫水的周期性波動又通過孔隙改變坡體的物理力學(xué)特性，兩方面的相互作用構(gòu)成了滑坡滲流場及應(yīng)力場的耦合作用問題[1-4]。庫區(qū)滑坡經(jīng)過多年的周期性庫水漲落，巖土體物理力學(xué)特性發(fā)生很大的變化，一些滑坡已近失穩(wěn)邊緣[5]，在庫水急劇下降時，由于滑坡體內(nèi)地下水來不及滲出，滑坡體內(nèi)外形成高水頭差，引起向外的滲流力，極大程度上影響滑坡體的穩(wěn)定性，特別是庫水下降期的后期為三峽庫區(qū)降雨集中期[6]，庫水下降聯(lián)合強(qiáng)降雨將可能引起庫區(qū)滑坡失穩(wěn)。

本文應(yīng)用大型通用有限元軟件ABAQUS對三峽庫區(qū)秭歸縣樹坪滑坡在庫水下降聯(lián)合降雨作用下的滲流場與應(yīng)力場進(jìn)行三維耦合分析，并結(jié)合監(jiān)測資料，研究庫水下降及降雨作用下滑坡變形與破環(huán)影響規(guī)律。

2 樹坪滑坡概況

圖1 樹坪滑坡全貌圖

樹坪滑坡位于秭歸縣沙鎮(zhèn)溪樹坪村三峽水庫右岸庫岸斜坡(如圖1所示)。該滑坡屬古崩滑堆積體,分布高程為65～500 m，滑坡體前緣突入長江，剪出口高程約65～68 m，滑體南北縱長約800 m，東西寬約700 m，面積約55萬m2，厚約30～70 m，總體積約為2 750萬m3。樹坪滑坡屬多期性巨型滑坡，物質(zhì)組成較復(fù)雜，根據(jù)地質(zhì)勘察和鉆探資料，從物質(zhì)組成大致可分為以下幾類：①耕植土層,主要為黃褐色、灰黃色以及紫紅色的粉質(zhì)黏土夾碎石；②坡積物(Qdl),主要為粉質(zhì)黏土層夾碎塊石，呈灰黃、淺褐黃色，粉質(zhì)黏土呈可塑狀態(tài)，碎塊石含量4%～10%，碎塊石多呈次棱角狀，以粉砂巖、泥質(zhì)粉砂巖為主，結(jié)構(gòu)松散～稍密，透水性較差；③滑坡堆積物(Qdel),碎塊石黏土層,呈棕黃色、紫紅色，土石比大小隨位置差異較大，碎塊石成分主要為泥質(zhì)粉砂巖、泥灰?guī)r和灰?guī)r。碎石呈次棱角狀，個別粒徑達(dá)15 cm。土的成分為壤土、粉質(zhì)黏土，呈硬塑-可塑狀，填充于碎塊石中，結(jié)構(gòu)不均勻，稍密，滑體堆積物結(jié)構(gòu)松散，透水性較好；④滑帶土,紫紅色角礫石土層，較濕，結(jié)構(gòu)緊密，土可塑。碎石呈次棱角狀～次圓狀，碎石上可見擦痕，土層中可見明顯揉皺、光滑鏡面?；瑤裆钶^大，層厚一般在10～20 cm；⑤基巖,滑坡地段為三迭系中統(tǒng)巴東組地層[7]。

為研究庫水及降雨作用對滑坡變形的影響，樹坪滑坡上布置有多種監(jiān)測儀器，現(xiàn)僅介紹與本文研究相關(guān)的GPS監(jiān)測、地下水位監(jiān)測。整個滑坡布設(shè)有8個GPS監(jiān)測點(diǎn)、2個基準(zhǔn)點(diǎn)，對滑坡體進(jìn)行大地變形監(jiān)測。并在前緣2個鉆孔中各安裝有地下水位監(jiān)測儀。樹坪滑坡平面圖如圖2所示。

圖2 樹坪滑坡平面圖

3 ABAQUS流固耦合求解方法

3.1 基本原理

庫水聯(lián)合降雨作用下滑坡體中孔隙流體壓力的變化會引起多孔介質(zhì)骨架有效應(yīng)力發(fā)生變化，這些變化又會反過來影響孔隙流體的流動和壓力的分布。滲流場與應(yīng)力場這種相互影響稱為流固耦合，在水庫水位波動和降雨條件下滑坡體中具有普遍性、基礎(chǔ)性、動態(tài)性和交叉性的特點(diǎn)。應(yīng)力平衡方程[8]和滲流連續(xù)方程為

IN-PN=0 ;

(1)

(2)

式中：IN為內(nèi)力矩陣；PN為外力矩陣；δuW為孔隙水壓力的變分；J表示土體體積的變化；ρW為水的密度；nW為土體孔隙率；V為體積；k為多孔介質(zhì)的滲透系數(shù)，一般由試驗測定；?為測壓管水頭，是位置水頭與壓力水頭之和： ??/?x為滲透梯度。

3.2 邊界條件

流固耦合方程的求解，還需給出相應(yīng)的定解條件[9]，本例滲流邊界條件為混合邊界條件，即一部分邊界水頭給定，另一部分流量給定，即降雨在坡面接觸部分的表面邊界為流量邊界。降雨引起坡面入滲及產(chǎn)流條件的確定是由降雨強(qiáng)度和坡面巖土體的入滲率決定的。當(dāng)降雨強(qiáng)度小于坡面巖土體的入滲率時，降雨產(chǎn)生的水流全部入滲至坡體；當(dāng)降雨強(qiáng)度大于坡面巖土體的入滲率時，坡面上將產(chǎn)生徑流。雨水在坡面的滲流流速v(由坡面巖土體的滲透系數(shù)k決定)和降雨強(qiáng)度q存在式(3)和式(4)的關(guān)系，根據(jù)該關(guān)系就可確定坡面入滲模擬的邊界條件：

v>q,取入滲速度為q， 單位為mm/s ；

(3)

v≤q,取入滲速度為v， 單位為mm/s。

(4)

固體應(yīng)力場邊界采用混合邊界條件，即固體骨架表面部分邊界已知應(yīng)力，另一部分邊界已知位移。

4 數(shù)值分析模型及初始條件

4.1 計算參數(shù)

根據(jù)樹坪滑坡地質(zhì)勘查報告中提供的有關(guān)巖土體的物理力學(xué)性質(zhì)試驗數(shù)據(jù)及參數(shù)建議值為基本依據(jù)，后選取樹坪滑坡具有代表性的2-2′剖面(如圖3)進(jìn)行智能位移反演，獲得有限元分析計算參數(shù)取值范圍如表1所示。

圖3 樹坪滑坡2-2′工程地質(zhì)剖面圖

表1樹坪滑坡有限元計算物理力學(xué)參數(shù)取值范圍

Table1Physico-mechanicalparametersofShupinglandslideforfiniteelementcalculation

部位重度γ/(kN·m-3)變形模量E/MPa泊松比μ黏聚力c/MPa內(nèi)摩擦角Φ/(°)滲透系數(shù)k/(m·d-1)滑體22.0524.750.1800.250311.388 880滑帶22.1021.600.2750.361250.025 920滑床26.6048 0000.2503.100410.002 592

4.2 計算模型

根據(jù)樹坪滑坡的地質(zhì)條件和地形地貌特征，選取樹坪滑坡三維數(shù)值計算模型的范圍：沿水流方向為725 m，垂直水流方向為840 m，模型底面高程為25 m(計算模型在滑坡平面圖上的位置如圖2所示)。計算域包含滑體、滑帶和基巖，整個計算域剖分了14 478個六面體單元，共計17 542個節(jié)點(diǎn)，三維計算模型與網(wǎng)格見圖4。

圖4 計算模型網(wǎng)格

4.3 計算初始條件

樹坪滑坡坡內(nèi)地下水位線往往隨外界因素(庫水、降雨)的變化而變化。根據(jù)地質(zhì)模型確定的縱剖面地下水分布及在滑坡體前緣182 m處鉆孔地下水位監(jiān)測數(shù)據(jù)來初步計算確定初始地下水位(圖5所示)，初始狀態(tài)孔壓分布見圖6，前緣地下水位基本與實測水位接近。計算域的底部應(yīng)力邊界為法向約束，前后兩側(cè)及左右兩側(cè)均采用法向約束，坡體表面為自由邊界。

圖5 樹坪滑坡前緣地下水位(庫水穩(wěn)定時)

圖6 樹坪滑坡三維計算模型初始狀態(tài)孔隙水壓力分布

4.4 計算工況

水庫水位變化和降雨作用是誘發(fā)庫區(qū)滑坡的最主要因素[10-11]，根據(jù)三峽庫區(qū)庫水調(diào)度及本項研究目標(biāo)制定計算工況如下。

工況1：庫水在175 m水位穩(wěn)定150 d；

工況2：庫水在175 m水位時約以0.186 m/d的速度下降至145 m水位，歷時161 d；

工況3：庫水在175 m水位約以0.186 m/d的速度下降至145 m水位的開始2 d時疊加100年一遇降雨(由文獻(xiàn)[6]及不同周期降雨強(qiáng)度計算得100年一遇降雨強(qiáng)度為180.24 mm/d)，歷時161 d；

工況4：庫水在175 m水位約以0.186 m/d的速度下降至145 m水位的最后2 d時疊加100年一遇降雨，歷時161 d。

5 計算結(jié)果及分析

為了分析樹坪滑坡在庫水下降及降雨作用下的穩(wěn)定性，通過能進(jìn)行巖土體飽和-非飽和流固耦合分析的有限元軟件(ABAQUS)對樹坪滑坡在4種工況下的滑坡變形破壞過程進(jìn)行模擬，獲得了該滑坡孔隙水壓力、位移場、應(yīng)力場和塑性區(qū)的變化過程、變化特征和變化規(guī)律。

圖7 滑坡體在各工況下穩(wěn)定后孔隙水壓力分布

5.1 孔隙水壓力結(jié)果

圖7為滑坡體在各工況下穩(wěn)定后孔隙水壓力分布圖。如圖7所示，工況2、工況3前緣后緣孔壓基本一致，說明滑坡體孔壓場受降雨影響較小，主要原因是對于整個滑坡體來說，總降雨量較小，不足以改變穩(wěn)定后的滲流場；對比圖7(a)與圖7(b)、圖7(c)，可見，庫水下降后的后緣孔壓由175 m穩(wěn)態(tài)時的1.957 MPa增大至2.521 MPa，說明庫水下降對滑坡地下水滲流場有較為明顯的影響，降雨的影響主要集中在降雨結(jié)束時滑坡體前緣的淺表層。

圖8 各工況下水平方向位移

5.2 位移場變化規(guī)律

圖8為滑坡體在各工況條件下水平方向位移圖，對比圖8中的(a)，(b)，(c)可見，庫水下降初期疊加降雨較庫水下降對滑坡位移的影響大，工況2和3的最大位移分別為25.44 cm和18.67 cm，由圖8中(c)，(d)可知，庫水下降后期時疊加降雨對位移的影響大于初期時疊加降雨，工況4的水平位移達(dá)到了44.05 cm，分析其原因為庫水下降后坡體內(nèi)應(yīng)力場及巖土體力學(xué)特性發(fā)生很大變化，此時遇到強(qiáng)降雨時可能導(dǎo)致滑坡局部失穩(wěn)。由上述4種工況計算結(jié)果可得到如下結(jié)論：強(qiáng)降雨對樹坪滑坡位移的影響較庫水作用大。

圖9 各工況下水平位移和庫水位變化趨勢圖

圖10 滑坡GPS變形圖

為進(jìn)一步分析庫水和降雨對樹坪滑坡變形的趨勢，由各工況計算結(jié)果繪制如圖9所示水平位移和庫水位變化關(guān)系圖。圖9(a)趨勢圖表明庫水下降初期滑坡水平位移突然增大，后逐漸減緩，但仍保持繼續(xù)增大的趨勢；圖9(b)趨勢圖表明庫水下降初期時疊加降雨對滑坡水位位移影響很大，后水平位移有減小的趨勢，但總體位移變形仍舊很大；圖9(c)趨勢圖表明工況4前期同工況2趨勢相同，僅在庫水下降末期時疊加降雨后有突然增大的位移趨勢，工況4為最危險工況。為比較數(shù)值模擬結(jié)果與實際位移趨勢，選取樹坪滑坡滑體東部GPS監(jiān)測點(diǎn)ZG86，ZG87數(shù)據(jù)與滑體西區(qū)中部點(diǎn)ZG88數(shù)據(jù)(如圖10)同計算結(jié)果進(jìn)行比較，發(fā)現(xiàn)滑坡體模擬位移變形趨勢同GPS基本吻合，如圖10所示，2012年6月初滑坡區(qū)受到強(qiáng)降雨突然發(fā)生大變形。工況3和工況4位移模擬結(jié)果也表明坡體發(fā)生大變形的部位與實際情況基本類似。

圖11 各工況下滑體第一主應(yīng)力圖

5.3 應(yīng)力場及塑性區(qū)

為分析庫水及降雨對樹坪滑坡應(yīng)力場的影響，選取各工況下滑體第一主應(yīng)力圖進(jìn)行分析(如圖11)，各工況下滑體第一主應(yīng)力均有明顯變化。圖11(a)表明滑體僅在庫水下降作用下只有滑體后部應(yīng)力發(fā)生較大變化，最大主應(yīng)力為0.436 9 MPa，滑體前緣應(yīng)力變化較小。圖11(b)、圖11(c)表明在庫水和降雨共同作用下滑體應(yīng)力值有較大變化，工況3最大值為0.389 6 MPa，主要集中在中后部，工況4最大主應(yīng)力為2.565 MPa，主要集中在滑體中部。由工況計算結(jié)果可知，影響滑坡應(yīng)力場的主要因素為降雨作用，次要因素為坡體結(jié)構(gòu)及庫水作用。

圖12所示為各工況下塑性區(qū)計算結(jié)果，分別比較各圖發(fā)現(xiàn)，各工況下滑體上均沒有出現(xiàn)明顯塑性區(qū)，僅在滑體后緣有局部塑性區(qū)，不影響滑坡整體穩(wěn)定，沒有出現(xiàn)塑性貫通現(xiàn)象。

圖12 各工況下塑性區(qū)分布

以上分析了樹坪滑坡在庫水下降聯(lián)合降雨作用下滲流場、位移場及應(yīng)力場變化規(guī)律，計算結(jié)果雖與實際情況接近，但三維流固耦合分析仍需取得新的發(fā)展。數(shù)值模擬結(jié)果的真實性取決于滑坡的初始狀態(tài)(初始應(yīng)力分布、初始孔壓分布等)，本文利用勘察資料及地下水位監(jiān)測資料初步獲得滑坡體初始孔壓分布，接近滑坡體真實狀態(tài)，但初始應(yīng)力場的準(zhǔn)確獲得仍較困難，尤其是庫水的周期性漲落作用，必將引起滑坡土體的蠕變和松弛，初始應(yīng)力狀態(tài)應(yīng)該結(jié)合滑坡土體的蠕變和松弛共同獲得。

6 結(jié) 論

樹坪滑坡屬于多期型巨型滑坡，一旦失穩(wěn)滑坡將會對長江通航及三峽大壩產(chǎn)生重大影響，本文通過研究樹坪滑坡監(jiān)測數(shù)據(jù)及應(yīng)用ABAQUS軟件對樹坪滑坡進(jìn)行數(shù)值模擬，并得到以下結(jié)論：

(1) 庫水和降雨作用對樹坪滑坡的變形影響較大，尤其是強(qiáng)降雨作用可能導(dǎo)致樹坪滑坡發(fā)生滑動。

(2) 庫水下降對樹坪滑坡位移變形的趨勢有先陡增后緩慢增加的特點(diǎn)。

(3) 庫水下降導(dǎo)致樹坪滑坡滑坡體應(yīng)力場及土體物理力學(xué)特性發(fā)生較大的變化，特別是庫水下降后期，此時若疊加強(qiáng)降雨將導(dǎo)致滑坡體的大變形，這也與CPS監(jiān)測所顯示樹坪滑坡常在庫區(qū)枯水期時經(jīng)歷強(qiáng)降雨發(fā)生大位移的結(jié)果相同。

參考文獻(xiàn)：

[1] 王 緩,徐志英,速寶玉.復(fù)雜裂隙巖體滲流與應(yīng)力彈塑性全耦合分析[J].巖石力學(xué)與工程學(xué)報,2000,19(2):177-181. (WANG Yuan, XU Zhi-ying, SU Bao-yu. Complete Coupled Analysis of Fluid Flow and Elastoplastic Stress in Complicated Fractured Rock Masses[J]. Chinese Journal of Rock Mechanics and Engineering, 2000,19(2):177-181. (in Chinese))

[2] 高海鷹,夏頌佑.三維裂隙巖體滲流場與應(yīng)力場耦合模型研究[J].巖土工程學(xué)報,1997,19(2):102-105.(GAO Hai-ying, XIA Song-you. Research on 3-D Model of Coupled Seepage and Stress Fields in Fracture Rock Mass[J]. Chinese Journal of Geotechnical Engineering, 1997, 19(2): 102-105. ( in Chinese))

[3] 仵彥卿,柴軍瑞.裂隙網(wǎng)絡(luò)巖體三維滲流場與應(yīng)力場禍合分析[J].西安理工大學(xué)學(xué)報,2000,16(1):1-5. (WU Yan-qing, CHAI Jun-rui. Analysis of 3-D Model of Coupled Seepage and Stress Fields in Fracture Network of Rock Mass[J]. Journal of Xi’an University of Technology, 2000, 16(1): 1-5. (in Chinese))

[4] 楊林德,楊志錫.各向異性飽和土體的滲流耦合分析和數(shù)值模擬[J].巖石力學(xué)與工程學(xué)報,2002,21(10):1447-1451.(YANG Lin-de，YANG Zhi-xi. Coupling Analysis and Numeric Simulations on Seepage Flow in Anisotropic Saturated Soils[J]. Chinese Journal of Rock Mechanics and Engineering，2002,21(10):1447-1451. (in Chinese))

[5] 黃潤秋，許 強(qiáng)，戚國慶.降雨及水庫誘發(fā)滑坡的評價與預(yù)測[M].北京：科學(xué)出版社，2007.(HUANG Run-qiu, XU Qiang， QI Guo-qing. Evaluation and Prediction of Landslide Induced by Rainfall and Reservoir[M]. Beijing: Science Press, 2007.(in Chinese))

[6] 王世梅，鄭 宏，童富果.水庫蓄水后三峽庫區(qū)重大水庫復(fù)活型滑坡空間預(yù)測評價研究[R].宜昌：三峽大學(xué)三峽庫區(qū)地質(zhì)災(zāi)害教育部重點(diǎn)實驗室,2013. (WANG Shi-mei, ZHENG Hong, TONG Fu-guo. Spatial Prediction Evaluation of Major Revived Landslides in Reservoir Area after the Impoundment of Three Gorges Reservoir[R].Yichang: Key Laboratory of Geological Hazards on Three Gorges Reservoir Area under Ministry of Education, China Three Gorges University,2013.(in Chinese))

[7] 易 武，尚 敏，黃海峰,等. 秭歸譚家河滑坡地質(zhì)模型研究報告[R]. 宜昌：三峽大學(xué)三峽庫區(qū)地質(zhì)災(zāi)害教育部重點(diǎn)實驗室，2010. (YI Wu, SHANG Min, HUANG Hai-feng,etal. Research Report of Geological Model for Tanjiahe Landslide in Zigui County[R].Yichang: Key Laboratory of Geological Hazards on Three Gorges Reservoir Area under Ministry of Education,China Three Gorges University, 2010. (in Chinese))

[8] 張 欣. 基于ABAQUS流固耦合理論的庫岸滑坡穩(wěn)定性分析[D].濟(jì)南：山東大學(xué),2005. (ZHANG Xin. Stability Analysis of Reservoir Landslide Based on the ABAQUS Seepage-Stress Coupling Theory[D]. Jinan: Shandong University, 2005. (in Chinese))

[9] 費(fèi) 康，張建偉.ABAQUS在巖土工程中的應(yīng)用[M].北京：中國水利水電出版社，2010.(FEI Kang, ZHANG Jian-wei. Application of ABAQUS in Geotechnical Engineering[M]. Beijing: China Water Power Press, 2010.(in Chinese))

[10] 張學(xué)年,盛祝平.長江三峽工程庫區(qū)順層岸坡研究[M].北京：地震出版社，1993. (ZHANG Xue-nian, SHENG Zhu-ping. Study on Layered Bank Slope along the Reservoir Area of the Yangtze River Three Gorges Project[M]. Beijing: Seismological Press,1993. (in Chinese))

[11] 杜榮桓，劉新民.長江三峽工程庫區(qū)滑坡與泥石流研究[M].成都：四川科學(xué)技術(shù)出版社，1991. (DU Rong-huan, LIU Xin-min. Study on Landslides and Debris-flows in the Reservoir Area of the Three Gorges Project[M]. Chengdu: Sichuan Science and Technology Press, 1991. (in Chinese))