库水下降联合降雨作用下树坪滑坡流固耦合分析

，世梅,

(三峡大学 三峡库区地质灾害教育部重点实验室，湖北 宜昌 443002)

1 研究背景

三峡水库蓄水以后，沿岸地区的自然条件发生显著性变化,一方面库水上升引起滑坡体地下水位抬升，从而改变滑坡体内岩土体的应力状态，减小滑坡体的有效应力；另一方面，库水的周期性波动又通过孔隙改变坡体的物理力学特性，两方面的相互作用构成了滑坡渗流场及应力场的耦合作用问题[1-4]。库区滑坡经过多年的周期性库水涨落，岩土体物理力学特性发生很大的变化，一些滑坡已近失稳边缘[5]，在库水急剧下降时，由于滑坡体内地下水来不及渗出，滑坡体内外形成高水头差，引起向外的渗流力，极大程度上影响滑坡体的稳定性，特别是库水下降期的后期为三峡库区降雨集中期[6]，库水下降联合强降雨将可能引起库区滑坡失稳。

本文应用大型通用有限元软件ABAQUS对三峡库区秭归县树坪滑坡在库水下降联合降雨作用下的渗流场与应力场进行三维耦合分析，并结合监测资料，研究库水下降及降雨作用下滑坡变形与破环影响规律。

2 树坪滑坡概况

图1 树坪滑坡全貌图

树坪滑坡位于秭归县沙镇溪树坪村三峡水库右岸库岸斜坡(如图1所示)。该滑坡属古崩滑堆积体,分布高程为65～500 m，滑坡体前缘突入长江，剪出口高程约65～68 m，滑体南北纵长约800 m，东西宽约700 m，面积约55万m2，厚约30～70 m，总体积约为2 750万m3。树坪滑坡属多期性巨型滑坡，物质组成较复杂，根据地质勘察和钻探资料，从物质组成大致可分为以下几类：①耕植土层,主要为黄褐色、灰黄色以及紫红色的粉质黏土夹碎石；②坡积物(Qdl),主要为粉质黏土层夹碎块石，呈灰黄、浅褐黄色，粉质黏土呈可塑状态，碎块石含量4%～10%，碎块石多呈次棱角状，以粉砂岩、泥质粉砂岩为主，结构松散～稍密，透水性较差；③滑坡堆积物(Qdel),碎块石黏土层,呈棕黄色、紫红色，土石比大小随位置差异较大，碎块石成分主要为泥质粉砂岩、泥灰岩和灰岩。碎石呈次棱角状，个别粒径达15 cm。土的成分为壤土、粉质黏土，呈硬塑-可塑状，填充于碎块石中，结构不均匀，稍密，滑体堆积物结构松散，透水性较好；④滑带土,紫红色角砾石土层，较湿，结构紧密，土可塑。碎石呈次棱角状～次圆状，碎石上可见擦痕，土层中可见明显揉皱、光滑镜面。滑带埋深较大，层厚一般在10～20 cm；⑤基岩,滑坡地段为三迭系中统巴东组地层[7]。

为研究库水及降雨作用对滑坡变形的影响，树坪滑坡上布置有多种监测仪器，现仅介绍与本文研究相关的GPS监测、地下水位监测。整个滑坡布设有8个GPS监测点、2个基准点，对滑坡体进行大地变形监测。并在前缘2个钻孔中各安装有地下水位监测仪。树坪滑坡平面图如图2所示。

图2 树坪滑坡平面图

3 ABAQUS流固耦合求解方法

3.1 基本原理

库水联合降雨作用下滑坡体中孔隙流体压力的变化会引起多孔介质骨架有效应力发生变化，这些变化又会反过来影响孔隙流体的流动和压力的分布。渗流场与应力场这种相互影响称为流固耦合，在水库水位波动和降雨条件下滑坡体中具有普遍性、基础性、动态性和交叉性的特点。应力平衡方程[8]和渗流连续方程为

IN-PN=0 ;

(1)

(2)

式中：IN为内力矩阵；PN为外力矩阵；δuW为孔隙水压力的变分；J表示土体体积的变化；ρW为水的密度；nW为土体孔隙率；V为体积；k为多孔介质的渗透系数，一般由试验测定；ø为测压管水头，是位置水头与压力水头之和： ∂ø/∂x为渗透梯度。

3.2 边界条件

流固耦合方程的求解，还需给出相应的定解条件[9]，本例渗流边界条件为混合边界条件，即一部分边界水头给定，另一部分流量给定，即降雨在坡面接触部分的表面边界为流量边界。降雨引起坡面入渗及产流条件的确定是由降雨强度和坡面岩土体的入渗率决定的。当降雨强度小于坡面岩土体的入渗率时，降雨产生的水流全部入渗至坡体；当降雨强度大于坡面岩土体的入渗率时，坡面上将产生径流。雨水在坡面的渗流流速v(由坡面岩土体的渗透系数k决定)和降雨强度q存在式(3)和式(4)的关系，根据该关系就可确定坡面入渗模拟的边界条件：

v>q,取入渗速度为q， 单位为mm/s ；

(3)

v≤q,取入渗速度为v， 单位为mm/s。

(4)

固体应力场边界采用混合边界条件，即固体骨架表面部分边界已知应力，另一部分边界已知位移。

4 数值分析模型及初始条件

4.1 计算参数

根据树坪滑坡地质勘查报告中提供的有关岩土体的物理力学性质试验数据及参数建议值为基本依据，后选取树坪滑坡具有代表性的2-2′剖面(如图3)进行智能位移反演，获得有限元分析计算参数取值范围如表1所示。

图3 树坪滑坡2-2′工程地质剖面图

表1树坪滑坡有限元计算物理力学参数取值范围

Table1Physico-mechanicalparametersofShupinglandslideforfiniteelementcalculation

部位重度γ/(kN·m-3)变形模量E/MPa泊松比μ黏聚力c/MPa内摩擦角Φ/(°)渗透系数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 计算模型

根据树坪滑坡的地质条件和地形地貌特征，选取树坪滑坡三维数值计算模型的范围：沿水流方向为725 m，垂直水流方向为840 m，模型底面高程为25 m(计算模型在滑坡平面图上的位置如图2所示)。计算域包含滑体、滑带和基岩，整个计算域剖分了14 478个六面体单元，共计17 542个节点，三维计算模型与网格见图4。

图4 计算模型网格

4.3 计算初始条件

树坪滑坡坡内地下水位线往往随外界因素(库水、降雨)的变化而变化。根据地质模型确定的纵剖面地下水分布及在滑坡体前缘182 m处钻孔地下水位监测数据来初步计算确定初始地下水位(图5所示)，初始状态孔压分布见图6，前缘地下水位基本与实测水位接近。计算域的底部应力边界为法向约束，前后两侧及左右两侧均采用法向约束，坡体表面为自由边界。

图5 树坪滑坡前缘地下水位(库水稳定时)

图6 树坪滑坡三维计算模型初始状态孔隙水压力分布

4.4 计算工况

水库水位变化和降雨作用是诱发库区滑坡的最主要因素[10-11]，根据三峡库区库水调度及本项研究目标制定计算工况如下。

工况1：库水在175 m水位稳定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年一遇降雨(由文献[6]及不同周期降雨强度计算得100年一遇降雨强度为180.24 mm/d)，历时161 d；

工况4：库水在175 m水位约以0.186 m/d的速度下降至145 m水位的最后2 d时叠加100年一遇降雨，历时161 d。

5 计算结果及分析

为了分析树坪滑坡在库水下降及降雨作用下的稳定性，通过能进行岩土体饱和-非饱和流固耦合分析的有限元软件(ABAQUS)对树坪滑坡在4种工况下的滑坡变形破坏过程进行模拟，获得了该滑坡孔隙水压力、位移场、应力场和塑性区的变化过程、变化特征和变化规律。

图7 滑坡体在各工况下稳定后孔隙水压力分布

5.1 孔隙水压力结果

图7为滑坡体在各工况下稳定后孔隙水压力分布图。如图7所示，工况2、工况3前缘后缘孔压基本一致，说明滑坡体孔压场受降雨影响较小，主要原因是对于整个滑坡体来说，总降雨量较小，不足以改变稳定后的渗流场；对比图7(a)与图7(b)、图7(c)，可见，库水下降后的后缘孔压由175 m稳态时的1.957 MPa增大至2.521 MPa，说明库水下降对滑坡地下水渗流场有较为明显的影响，降雨的影响主要集中在降雨结束时滑坡体前缘的浅表层。

图8 各工况下水平方向位移

5.2 位移场变化规律

图8为滑坡体在各工况条件下水平方向位移图，对比图8中的(a)，(b)，(c)可见，库水下降初期叠加降雨较库水下降对滑坡位移的影响大，工况2和3的最大位移分别为25.44 cm和18.67 cm，由图8中(c)，(d)可知，库水下降后期时叠加降雨对位移的影响大于初期时叠加降雨，工况4的水平位移达到了44.05 cm，分析其原因为库水下降后坡体内应力场及岩土体力学特性发生很大变化，此时遇到强降雨时可能导致滑坡局部失稳。由上述4种工况计算结果可得到如下结论：强降雨对树坪滑坡位移的影响较库水作用大。

图9 各工况下水平位移和库水位变化趋势图

图10 滑坡GPS变形图

为进一步分析库水和降雨对树坪滑坡变形的趋势，由各工况计算结果绘制如图9所示水平位移和库水位变化关系图。图9(a)趋势图表明库水下降初期滑坡水平位移突然增大，后逐渐减缓，但仍保持继续增大的趋势；图9(b)趋势图表明库水下降初期时叠加降雨对滑坡水位位移影响很大，后水平位移有减小的趋势，但总体位移变形仍旧很大；图9(c)趋势图表明工况4前期同工况2趋势相同，仅在库水下降末期时叠加降雨后有突然增大的位移趋势，工况4为最危险工况。为比较数值模拟结果与实际位移趋势，选取树坪滑坡滑体东部GPS监测点ZG86，ZG87数据与滑体西区中部点ZG88数据(如图10)同计算结果进行比较，发现滑坡体模拟位移变形趋势同GPS基本吻合，如图10所示，2012年6月初滑坡区受到强降雨突然发生大变形。工况3和工况4位移模拟结果也表明坡体发生大变形的部位与实际情况基本类似。

图11 各工况下滑体第一主应力图

5.3 应力场及塑性区

为分析库水及降雨对树坪滑坡应力场的影响，选取各工况下滑体第一主应力图进行分析(如图11)，各工况下滑体第一主应力均有明显变化。图11(a)表明滑体仅在库水下降作用下只有滑体后部应力发生较大变化，最大主应力为0.436 9 MPa，滑体前缘应力变化较小。图11(b)、图11(c)表明在库水和降雨共同作用下滑体应力值有较大变化，工况3最大值为0.389 6 MPa，主要集中在中后部，工况4最大主应力为2.565 MPa，主要集中在滑体中部。由工况计算结果可知，影响滑坡应力场的主要因素为降雨作用，次要因素为坡体结构及库水作用。

图12所示为各工况下塑性区计算结果，分别比较各图发现，各工况下滑体上均没有出现明显塑性区，仅在滑体后缘有局部塑性区，不影响滑坡整体稳定，没有出现塑性贯通现象。

图12 各工况下塑性区分布

以上分析了树坪滑坡在库水下降联合降雨作用下渗流场、位移场及应力场变化规律，计算结果虽与实际情况接近，但三维流固耦合分析仍需取得新的发展。数值模拟结果的真实性取决于滑坡的初始状态(初始应力分布、初始孔压分布等)，本文利用勘察资料及地下水位监测资料初步获得滑坡体初始孔压分布，接近滑坡体真实状态，但初始应力场的准确获得仍较困难，尤其是库水的周期性涨落作用，必将引起滑坡土体的蠕变和松弛，初始应力状态应该结合滑坡土体的蠕变和松弛共同获得。

6 结 论

树坪滑坡属于多期型巨型滑坡，一旦失稳滑坡将会对长江通航及三峡大坝产生重大影响，本文通过研究树坪滑坡监测数据及应用ABAQUS软件对树坪滑坡进行数值模拟，并得到以下结论：

(1) 库水和降雨作用对树坪滑坡的变形影响较大，尤其是强降雨作用可能导致树坪滑坡发生滑动。

(2) 库水下降对树坪滑坡位移变形的趋势有先陡增后缓慢增加的特点。

(3) 库水下降导致树坪滑坡滑坡体应力场及土体物理力学特性发生较大的变化，特别是库水下降后期，此时若叠加强降雨将导致滑坡体的大变形，这也与CPS监测所显示树坪滑坡常在库区枯水期时经历强降雨发生大位移的结果相同。

参考文献：

[1] 王 缓,徐志英,速宝玉.复杂裂隙岩体渗流与应力弹塑性全耦合分析[J].岩石力学与工程学报,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] 高海鹰,夏颂佑.三维裂隙岩体渗流场与应力场耦合模型研究[J].岩土工程学报,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] 仵彦卿,柴军瑞.裂隙网络岩体三维渗流场与应力场祸合分析[J].西安理工大学学报,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] 杨林德,杨志锡.各向异性饱和土体的渗流耦合分析和数值模拟[J].岩石力学与工程学报,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] 黄润秋，许 强，戚国庆.降雨及水库诱发滑坡的评价与预测[M].北京：科学出版社，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] 王世梅，郑 宏，童富果.水库蓄水后三峡库区重大水库复活型滑坡空间预测评价研究[R].宜昌：三峡大学三峡库区地质灾害教育部重点实验室,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] 易 武，尚 敏，黄海峰,等. 秭归谭家河滑坡地质模型研究报告[R]. 宜昌：三峡大学三峡库区地质灾害教育部重点实验室，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流固耦合理论的库岸滑坡稳定性分析[D].济南：山东大学,2005. (ZHANG Xin. Stability Analysis of Reservoir Landslide Based on the ABAQUS Seepage-Stress Coupling Theory[D]. Jinan: Shandong University, 2005. (in Chinese))

[9] 费 康，张建伟.ABAQUS在岩土工程中的应用[M].北京：中国水利水电出版社，2010.(FEI Kang, ZHANG Jian-wei. Application of ABAQUS in Geotechnical Engineering[M]. Beijing: China Water Power Press, 2010.(in Chinese))

[10] 张学年,盛祝平.长江三峡工程库区顺层岸坡研究[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] 杜荣桓，刘新民.长江三峡工程库区滑坡与泥石流研究[M].成都：四川科学技术出版社，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))