强降雨作用下昔格达边坡渗流特性及稳定性分析

杨世豪，苏立君，张崇磊，李丞，胡兵立

(1.中国科学院 山地灾害与地表过程重点实验室；中国科学院、水利部 成都山地灾害与环境研究所，成都 610041；2.中国科学院 青藏高原地球科学卓越创新中心，北京 100101；3.中国科学院大学，北京 100049)

昔格达地层为广泛分布于中国西南地区的半成岩河湖相沉积地层，形成于晚更新世与上第三系之间，主要由粉砂土和黏土组成[1-2]。随着工程活动在昔格达地层分布区广泛开展，其工程性质引起学者关注。周平等[3]、王志杰等[4]研究发现，昔格达地层岩土体水稳性差，极易崩解。杜翔宇等[5]研究昔格达地层岩土体的微观结构特性，并对其强度指标与含水率关系进行阐述。因力学性质特殊，昔格达地层分布区滑坡频发，相关学者对其进行了研究。黄绍槟等[6]将昔格达滑坡归为覆盖层滑坡、昔格达基底滑坡和昔格达组层面滑坡3类。丁文富等[7]提出，昔格达地层与上覆松散堆积层的接触界面、与下伏地层的界面以及层间软弱沉积结构面是该地层中主要的不良结构面。综上，昔格达地层边坡多为层状结构且滑面多位于不同土层接触面上。

降雨是常见滑坡诱发因素[8]，对昔格达地层岩土体遇水软化的特殊土体而言，这种作用尤为显著。降雨对边坡稳定性的主要影响在于土体孔隙水压力增加造成抗剪强度降低[9]。此外，雨水入渗将在边坡表层岩土体形成暂态饱和区，使该区域土体自重增加，下滑力增大，边坡稳定性降低[10]。White等[11]的研究采用降雨条件下土体体积含水率、孔隙水压力的变化规律来表征边坡的渗流特性。鉴于降雨入渗会对边坡稳定性造成严重影响，相关研究建立了边坡渗流特性与稳定性的关系[12]。张建等[13]对降雨触发浅层坡体失稳的迟滞现象及其与土质参数的关联性进行研究。层状边坡渗流规律较均质边坡更复杂。韩同春等[14]认为，湿润锋至接触面时，引起接触面孔隙水压力上升是滑坡主要诱因。马吉倩等[15]通过数值模拟方法，得到降雨条件及坡积土层厚度对含水率、孔隙水压力沿高程分布的影响规律。石振明等[16]通过改进Green-Ampt模型，提出考虑降雨入渗的多层非饱和土边坡稳定性分析方法。上述研究阐述了层状边坡降雨渗流特征与稳定性变化规律。

昔格达地层边坡多为层状边坡，但目前针对降雨入渗作用下昔格达地层边坡内部渗流场分布规律研究较少。以四川省雅安市石棉县莫家岗滑坡为例，利用有限元方法模拟经历降雨边坡内渗流过程，并以Morgenstern-Price方法计算边坡各时刻稳定性系数。通过研究坡内渗流特点及边坡在降雨作用下稳定性系数的变化规律，为昔格达地层滑坡治理提供参考。

1 计算理论

降雨入渗实质为坡内土体含水率上升，由非饱和状态向饱和状态发展的过程。入渗过程中，含水率随深度分布可分为4部分，即：浅层饱和带、含水率变化较大的过渡带、含水率均匀分布的传导带和湿润度随深度减小的湿润区，湿润区前缘为湿润锋[17]。

坡体中，水的运动基本方程为根据达西定律及质量守恒定律推导得到的RICHARDS方程[18]。

式中：x为水平方向；z为竖直方向；hm为基质吸力水头；k(hm)为非饱和土渗透系数函数；C(hm)为比水容量及水土特征曲线斜率；t为时间。其中，非饱和土渗透系数k(hm)有别于饱和土，达西定律中，k(hm)等于常数渗透系数，为与体积含水率相关的函数，反映非饱和土中的水力关系对渗流的影响。

Morgenstern-Price法可对任意形状滑动面进行求解，同时，满足竖向力、水平力及力矩的平衡，是一种常用极限平衡求解方法。基本计算原理为，假定两相邻土条的法向条间力和切向条间力之间存在一对水平方向坐标的函数关系，根据整个滑动土体边界条件进行迭代，求出问题的解[7]。作为典型的复杂条分法，其计算过程比一般条分法复杂，但计算结果更为准确，因而，适用于利用数值模拟软件进行求解。

2 计算模型及工况设计

2.1 工程概况

莫家岗滑坡位于石棉县美罗乡狮子村4组，后缘地理坐标：N29°17′07.9″，E102°26′52.8″，主滑方向为184°,剖面呈上陡下缓形态，主滑方向坡高80 m，总体坡度16°，为中型浅层土质滑坡。坡体为典型层状结构，上部为第四系全新统残坡积物(Q4el+dl)，黄褐色粉质粘土，局部夹碎块石，厚度1.0～4.0 m；下部为第三系中统昔格达组粉砂土，强度较差；基岩为第三系中统昔格达组粉砂岩。

2.2 计算参数及边界条件

由试验测得,滑坡中第四纪残坡积物及昔格达组粉砂土相关参数见表1。残坡积物水土特征曲线由软件内置样本函数拟合；昔格达粉砂土水土特征曲线由粒径级配曲线拟合。两类土渗透系数函数根据相应土壤水分特征曲线，由V-G模型拟合[19]，昔格达组粉砂土粒径级配曲线如图1所示，两类土水土特征曲线及渗透系数函数曲线如图2所示。

表1 物理力学参数Table 1 Physical and mechanical parameters

图1 昔格达组砂土粒径级配Fig.1 Particle size distribution of Xigeda Formation soil

图2 两种土水土特征曲线及渗透系数曲线Fig.2 SWCC & permeability coefficient curve

采用Geo-studio中的seep/w模块计算边坡渗流规律，geostudio是适用于岩土工程数值模拟的计算软件，seep/w为软件内专门用于非饱和土渗流计算的分析模块。seep/w可进行稳态以及瞬态渗流计算分析，可以获得土体内部孔隙水压力、体积含水率等指标的空间及时间变化规律。此外，seep/w模块的计算结果可与其他模块进行耦合计算，如与SLPOE/W耦合，获得考虑边坡内部基质吸力空间以及时间分布条件下边坡的稳定性系数。

莫家岗滑坡平面图利用无人机航拍图像获得，如图3所示，剖面图如图4所示，其中，剖面图参考文献[20]。坡体地下水位以上网格尺寸为1 m，以下网格尺寸为5 m。坡顶、坡中及坡脚设置A、B、C三处观测面。稳态分析中，模型上表面为自由边界，两侧及底部为不透水边界。瞬态分析包括降雨和不降雨两段，其中，降雨期间，坡体表面设为单位流量边界，流量等于降雨量，雨停后为自由边界。

结合石棉县当地情况及中国气象局关于降雨强度的划分标准，设定模拟共历时7 d，前4天为降雨过程，降雨强度分别设为工况1大暴雨(0.006 m/h)和工况2大雨(0.001 6 m/h)两种，后3天无降雨作用。

图3 莫家岗滑坡平面图Fig.3 Plane map of Mojiagang landslide

图4 莫家岗滑坡剖面图(单位：m)Fig.4 Profile of Mojiagang landslide(unit:m)

3 计算结果分析

3.1 降雨作用下边坡含水率变化情况

3.1.1 坡顶监测面 图5、图6为两种工况下模拟过程监测面A的含水率分布。在大暴雨工况下，t=48 h，坡体出现饱和区，随降雨持续，饱和区扩大；t=96 h，最大入渗深度约3 m，饱和区深度达2 m；雨停后，饱和区消散，雨水继续入渗。t=168 h，入渗深度达3.5 m。A处渗流仅发生在覆盖层中，为均质土体中渗流。降雨时，含水率沿深度分布曲线含饱和区、过渡区、传导区、湿润区及湿润锋，与张建等[13]的结论一致。雨停后，高程87.5 m以上同深度含水率，随时间延长而减少；高程87.5 m以下同深度含水率，随时间延长而增大，湿润锋深度增加。在大雨工况条件下，降雨入渗总深度明显减小，t=96 h时，降雨入渗深度仅为2 m，但与大暴雨工况下的渗流场特征相对比，可以看出，在降雨结束后，大雨工况条件下，雨水的入渗深度及同深度处的含水率都在继续增长。

图5 大暴雨工况下坡顶监测面A体积含水率分布图Fig.5 The water content distribution with elevation on section A in situation 1

图6 大雨工况下坡顶监测面A体积含水率分布图Fig.6 The water content distribution with elevation on section A in situation 2

3.1.2 坡中监测面含水率 图7、图8为两种工况下模拟过程监测面B的含水率分布。在大暴雨工况下，0

图7 大暴雨工况下坡中监测面B体积含水率分布图Fig.7 The water content distribution with elevation on section B in situation 1

图8 大雨工况下坡中监测面B体积含水率分布图Fig.8 The water content distribution with elevation on section B in situation 2

3.1.3 坡脚监测面含水率 图9、图10为两种工况下模拟过程监测面C的含水率分布。坡脚覆盖层厚度1.04 m，t=48 h时，覆盖层已饱和，且雨水渗入昔格达地层中；t=96 h时，雨水最终入渗至粉砂土层内约3 m。因入渗深度范围内两种土体渗透性质不同，含水率分布曲线分层明显。降雨过程中，覆盖层迅速饱和，因而，从t=48 h开始，该层土体处于饱和状态；雨停后，层面附近，残坡积物含水率约0.35，昔格达组粉砂土含水率约0.17，表明，雨停后层面附近土体含水率较高。在大雨工况下，由于坡脚覆盖层厚度较薄，降雨期间，雨水已经渗透至昔格达地层内部，但由于入渗量较小，未形成明确分区。与大暴雨工况下的渗流场特征进行对比可发现，在降雨结束后，水分在坡体内的最大入渗深度仍在不断增加，且同深度处的土体含水率不断增加。

图9 大暴雨工况小坡脚监测面C体积含水率分布图Fig.9 The water content distribution with elevation on section C in situation 1

图10 大雨工况下坡脚监测面C体积含水率分布图Fig.10 The water content distribution with elevation on section C in situation 2

3.2 降雨作用下边坡孔隙水压力变化情况

3.2.1 坡顶监测面孔隙水压力变化情况 图11、图12为两种工况下模拟过程中监测面A孔隙水压力分布情况。两种工况下，A处渗流均发生在覆盖层中，为均质土体中渗流。在大暴雨工况下，雨水入渗深度范围内，孔隙水压力持续增长，土体发生孔隙水压力增长范围匀速拓展；由图11、图12可知，对于均质土体，雨停后，孔隙水压力随高程分布曲线可分为两个段，高程87.5 m以上土体中孔隙水压力消散；87.5 m以下土体中湿润锋深度增加，同深度孔隙水压力增加。在大雨工况下，水分渗透渗透深度更小，因此，全过程仅发生在上层覆盖物中。可以看出，模拟全过程土体中，孔隙水压力始终在逐渐增大，但未达到正值。

图11 大暴雨工况小坡顶监测面A孔隙水压力分布图Fig.11 The pore water pressure distribution with elevation on section A in situation 1

图12 大雨工况下坡顶监测面A孔隙水压力分布图Fig.12 The pore water pressure distribution with elevation on section A in situation 2

3.2.2 坡中监测面孔隙水压力变化情况 图13、图14为两种工况下模拟过程中监测面B孔隙水压力分布。在大暴雨工况下，0

图13 大暴雨工况下坡中监测面B孔隙水压力分布图Fig.13 The pore water pressure distribution with elevation on section B in situation 1

图14 大雨工况下坡中监测面B孔隙水压力分布图Fig.14 The pore water pressure distribution with elevation on section B in situation 2

3.2.3 坡脚监测面孔隙水压力变化情况 图15、图16为两种工况下模拟过程中监测面C孔隙水压力分布。在大暴雨工况下，因覆盖层厚度仅1.04 m，t=48 h时，昔格达组粉砂土中出现饱和区，孔隙水压力增长至正值，层面附近孔隙水压力增长范围加速拓展；由图15、图16可知，因两种土质渗流特性差异，雨停后，坡内孔隙水压力随高程分布曲线存在分区，分区特征与坡中监测面B类似，覆盖层中，近坡体表面孔隙水压力减少明显；昔格达粉砂土地层中，t=96 h时，昔格达粉砂土中含水率曲线饱和区与过渡区交点为同时刻坡内孔隙水压力最大部位，该深度以下同深度孔隙水压力随时间增长而增长，湿润锋深度加大。在大雨工况下，孔隙水压力随降雨过程的发展而不断增长，但在t=96 h时，土体内部最大孔孔隙水压力仍未达到正值，降雨结束后，由于土体中水分继续入渗，孔隙水压力逐渐增长为正值。

图15 大暴雨工况下坡脚监测面C孔隙水压力分布图Fig.15 The pore water pressure distribution with elevation on section C in situation 1

图16 大雨工况下坡脚监测面C孔隙水压力分布图Fig.16 The pore water pressure distribution with elevation on section C in situation 2

3.3 昔格达地层层状结构边坡渗流特征讨论

参考对层状结构边坡渗流特性的相关研究成果，并结合《工程地质手册》等资料，对普通岩土体相关参数的总结可以发现，与普通层状结构土质边坡相比，昔格达地层边坡中，上层松散堆积物覆盖层的渗透系数及饱和含水率的参数与下层昔格达地层岩土体的渗透系数及饱和含水率值之间的差距更大。这种层状结构边坡中，上下层土体之间在渗透系数以及饱和含水率等参数相差悬殊的特点，使得昔格达地层层状结构边坡在经历降雨过程时，在土层交界面处土体含水率以及孔隙水压力的突变较普通层状结构土质边坡更为剧烈，突变量值更大。

昔格达地层层状结构边坡所具有的相较于普通层状结构土质边坡更为剧烈的土层交界面处，土体含水率以及孔隙水压力的突变表明，该类边坡更易在土层交界面附近形成高含水率区域，且该区域内含水率更高，持续时间更久，消散更缓慢，这一特点对边坡稳定性造成较大影响。首先，由于接触面附近持续处于高含水率状态，降雨结束后，边坡内部仍会持续发生雨水渗流过程，使得边坡安全系数持续下降，其滞后性较普通层状结构土质边坡更为明显。其次，由于地层接触面附近的覆盖层土体含水率较高，且不能较快消散，相较于普通层状结构土质边坡，昔格达地层层状结构边坡土体内部水分更易发生沿土层接触面的流动，对边坡内部土颗粒产生较大渗流力，甚至对土体颗粒产生淘蚀。

3.4 边坡稳定性系数变化情况

图17为两种工况条件下边坡稳定性系数随降雨入渗的变化曲线。由图17可知，两种工况下，边坡安全系数的变化规律基本一致，降雨期间，稳定性系数持续下降，下降速率先增大后减小，反映出雨水入渗至坡体不同部位时，稳定性系数变化规律不同；24 h40 h时，水分渗入昔格达组粉砂土后，由于渗透深度持续增加导致边坡稳定性系数减小，但由于粉砂土渗透系数小于上层第四纪残坡积物，雨水下渗速率减慢，导致稳定性系数下降速率减慢。降雨结束后，由于坡体中水分并不能快速排出而是继续下渗，导致稳定性系数在雨停后继续下降；t=120 h时，稳定性系数降至最小；t=144 h为降雨结束后48 h，土体中部分水分排出，稳定性系数上升，后基本保持稳定，表明，边坡安全系数的变化趋势与降雨过程并不同步，反映出降雨对层状边坡稳定性影响的滞后性。但由于大雨工况下坡体雨水入渗深度整体较小，因此，大雨工况下边坡的安全系数高于大暴雨工况下边坡的安全系数。

3.5 边坡破坏机理分析

根据相关研究结论，莫家岗滑坡为受历时长的强降雨诱发的昔格达地层牵引式滑坡，滑面位于基覆界面部位，其变形过程与降雨过程密切相关。研究区滑坡最早发生于2008年汛期，降雨导致斜坡后缘岩土体发生强烈变形；2012年雨季斜坡后缘再次产生张拉裂缝；受暴雨影响，2013年7月斜坡发生进一步强烈变形，形成局部错台，斜坡中下部变形十分明显，附近乡村公路受路基出现不均匀沉降；2014年7月，研究区突降暴雨，致使斜坡局部强烈变形；2014年8月，研究区持续强降雨，斜坡体多处发生变形，斜坡后缘产生张拉裂缝[20]。

图17 安全系数变化规律Fig.17 Safety factor versus time

结合渗流场数值模拟结果，可以将莫家岗滑坡破坏类型概述为两种，即局部浅层滑塌以及深层牵引型滑坡。莫家岗滑坡为典型的层状边坡，雨水在边坡中逐层渗入，降雨入渗，首先对边坡上层第四纪残坡积物覆盖层产生影响，造成覆盖层土体处于饱和状态，加之表层土体受人类活动改造剧烈，易引起局部浅层滑塌。结合渗流场模拟结果以及边坡结构特点，可以对边坡深层牵引型滑坡的破坏机理做出解释：边坡中下部第四纪残坡积物覆盖层厚度较薄，降雨可渗流至下层昔格达组地层中，且由于两种土体的水土性质差异悬殊，造成界面处长期处于高含水率状态，由于昔格达地层土体力学性质差，在长期受雨水作用下，易造成坡体下部岩土体发生滑面位于土层交界面的变形破坏，因而，边坡强烈变形多发生于坡体中下部。坡体上部覆盖层土体厚度较大，雨水不易渗入下层昔格达地层中，因此，上部坡体能保持较好完整性，但由于下部坡体力学性质不断减弱，上部坡体抗力减少，整体受牵引下滑，后缘产生拉裂缝。

4 结论

1)昔格达地层边坡多为层状结构，该工程为上覆第四纪坡残积物的昔格达组粉砂土边坡；降雨过程中，雨水在坡体中逐层入渗；雨停之后，入渗过程仍将持续，但残坡积物与昔格达组粉砂土中含水率及孔隙水压力分布特征有明显差异。

2)土层交界面处含水率、孔隙水压力变化明显；界面上部土体长期处于高含水率状态，使昔格达地层靠近界面处土体长期受水分软化作用，易形成软弱夹层，诱发滑面位于界面的滑坡；土体中孔隙水压力发生增长的影响范围在界面附近会加速拓展，使边坡抗滑力减弱，导致稳定性系数加速降低。

3)降雨入渗导致边坡稳定性降低，当雨水入渗至坡体不同部位时，稳定性系数表现出不同的变化规律；降雨结束后，稳定性系数不会立即回升，表明边坡对降雨入渗的响应存在滞后性。