强降雨条件下碎屑岩滑坡远程运动模拟分析——以牛儿湾滑坡为例

吴伟乐 贺 凯高 杨李 滨刘朋飞

1.长安大学地质工程与测绘学院,陕西 西安 710054;

2.中国地质科学院地质力学研究所,北京 100081;

3.自然资源部活动构造与地质安全重点实验室,北京 100081;

4.中国地质调查局新构造与地壳稳定性研究中心,北京 100081;

5.重庆市地质环境监测总站,重庆 401122

0 引言

全球气候变化对极端天气事件(极端降雨、气温升高、强风和洪水灾害)的影响尤为强烈,增加了地质灾害发生的风险(高杨等,2017),尤其是强降雨条件下的高位远程滑坡灾害频发,给人类的生命和财产安全带来了极大威胁。高位远程滑坡是从高陡斜坡凌空剪出、加速下坠的一类滑坡,常常具有较大的能量,在运动过程中易转化为泥石流或碎屑流(许强等,2009;殷跃平,2010;殷跃平等,2017)。中国西南砂泥岩地层山区由暴雨灾害引发了多起高位远程地质灾害事件,如2014年9月初渝东北地区发生强降雨事件,导致多地多处发生远程滑坡(Li et al., 2022),主要有重庆奉节无山坪滑坡、咸池水库滑坡、白果寨滑坡。国内学者对强降雨型高位远程岩质滑坡的成灾模式和远程运动机理开展了相关研究(许强,2010;高杨等,2020,2022;高浩源等,2020;李壮等,2020),认为强降雨是导致滑坡远程运动的关键因素。也有学者对顺层滑坡开展了系统的研究,余飞等(2005)采用数值模拟方法,展现了顺层滑坡的渐进破坏过程;邹宗兴等(2012)将顺层岩质滑坡的变形破坏地质力学模式分为了两类,并给出了相应的稳定性计算方法;龙建辉等(2019)探讨了含有软弱夹层的顺层岩质滑坡的变形破坏规律;胡乐等(2021)通过离心机模型实验,探索了硬土软岩型滑坡沿顺层破坏时的运动特点。滑坡运动过程中的动能保持和基底摩擦阻力是影响滑坡远距离运动本质条件,不少学者提出了诸如气垫效应、流化效应和液化效应等观点,主要从摩擦阻力下降的角度解释了其远程致灾的现象(Buss et al., 1881;Kent,1966;Sassa, 1989;Hungr,1995;Davies et al.,1999;邢爱国等,2004;邢爱国和殷跃平,2009;De Blasio and Crosta, 2015;Lai et al., 2017;Gao et al, 2022)。

随着计算机技术的发展,数值模拟技术因其成本低、操作方便、可重复试验等优点,成为了反演分析滑坡运动过程和动力学机制的最为高效的手段之一。基于等效流体理论以及非连续力学方法,对滑坡运动的数值模拟主要包括连续介质算法以及离散元算法。①连续介质算法:如DANW、SPH2D、FlAC3D等基于剖面的数值模拟技术,把复杂的滑坡体看成等效流体,通过不同的流变模型确定基底的摩擦阻力,可以较好地模拟滑坡的运动速度和堆积厚度(Hungr,1995;Evans,2001;Sassa et al.,2004,2010;张远娇等,2012;Gao et al.,2017,2019,2020;张艳玲等,2021);②离 散 元 算 法:如PFC、EDEM和MatDEM等计算方法,通过设置不同的接触类型和摩擦参数,可以有效地模拟崩塌、碎屑流等大变形问题。如张龙等(2012)利用 PFC3D软件分析了摩擦因数与黏结强度的变化对鸡尾山滑坡的运动、堆积特征的影响。孟桓羽等(2022)利用离散元软件PFC3D,对山阳滑坡运动全过程进行模拟,再现了滑坡体失稳后,在重力作用下沿视倾向起动、碰撞和堆积的全过程。

牛儿湾滑坡位于重庆市武隆区白马镇鱼光村,在持续降雨影响下,滑体于2020年7月13日沿层面发生整体滑移失稳,造成8栋房屋倒塌、S529省道及页岩气输气管道损毁。目前仍有约273.5×104m3滑体处于基本稳定状态,若遭遇强降雨天气,极有可能再次失稳下滑,严重威胁居民房屋、S529省道、村道、高压线铁塔等设施的安全。文章基于对牛儿湾滑坡的详细野外工程地质调查与分析,系统阐述了滑坡体运动堆积特征,并利用离散元软件PFC3D还原了牛儿湾滑坡运动过程及堆积特征,揭示了砂泥岩地层山区滑坡整体失稳和远程流化运动的成灾模式,认为在砂泥岩地层山区的远程流化滑坡风险调查与预测过程中,应当充分基于滑体远程流化运动的成灾特点进行调查与评价,为防灾减灾工作提供定量化科学依据。

1 滑坡地质环境条件

1.1 自然地理

研究区位于重庆市武隆区白马镇,经纬度为107°30′41.82″E,29°23′59.47″N,属于低山地貌区(图1)。研究区西侧和西南侧为山顶,总体为西高东低,北部地形转折,形成小山脊,继续往北转折为冲沟,向东延伸至研究区东部(滑坡前缘)(图2)。目前已滑移区地势较周边低,形成凹槽。除已滑移区外,研究区总体为上缓下陡,上部地形坡角一般为10°～15°,下部地形坡角一般为20°～28°,整体坡向约115°;最高点位于西侧,高程约为625～666 m,最低点位于下部前缘堆积体所处冲沟,高程约为334～339 m,相对高差约291～327 m。

图2 牛儿湾滑坡遥感影像图及现场照片Fig.2 Remote sensing images and pictures of the Niuerwan landslide

1.2 地质构造及地层岩性

研究区位于扬子准地台四川台坳川东南坳褶带南端与川黔经向构造带交汇部位,主控构造为白马向斜,轴线位于白马镇以西。向斜南段轴向为南西向;中段轴向由南西向转为北北东向;北段轴向为北北东向。由于研究区处于复杂的断裂和褶皱的交汇地带,强烈构造活动使得区内岩体较为破碎(图3)。

图3 研究区区域地质图Fig.3 Regional geological map of the study area

图4 滑坡岩体赤平投影及出露基岩Fig.4 Stereographic projection of the rock mass and outcropping of the bedrock

1.3 气象水文

武隆区属典型的亚热带湿润季风气候,降雨充沛,多集中在5～10月,年平均降水量1100.1 mm,月均降水量见图5。研究区为斜坡地形,总体西高东低,有利于地表水的排泄,滑坡区域内无常年地表水,主要通过冲沟及公路排水沟进行地表水的排泄。

图5 研究区多年月均降雨量图Fig.5 Monthly rainfall in the study area

地表水易向冲沟汇集,且区内岩土体裂隙发育,故地表水易向岩土体中下渗;同时,由于部分路段的公路排水沟直接排向路边岩土体中,也极易造成岩土体饱水。在这样不良的排水条件下,持续性降雨成为该滑坡失稳的重要诱因,松散的残坡积层及发育的基岩裂隙,导致雨水大量入渗,不但加大了滑体容重,还削弱了岩土体的抗剪强度,极大地降低了滑坡稳定性。

2 滑坡运动特征

根据滑坡发生前后的遥感图对比和现场调查分析显示,在持续的强降雨条件下,滑坡发生整体失稳,随后滑体前部开始下滑。失去前部的阻挡,滑体中后部随之下滑,在运动一段距离后,受山脊地形影响,滑体分流为N45°E和N81°E两部分,最大运动距离约1.4 km,前缘呈尖灭状。滑坡面积约为13.4×104m2,堆积体平均厚度约为11.08 m,体积约为147.5×104m3(图6,图7),造成8栋房屋倒塌、S529省道及页岩气输气管道损毁。根据滑坡的运动堆积特征,将滑坡分为了滑源区、主堆积区和碎屑流堆积区(图6)。

图6 牛儿湾滑坡平面图Fig.6 Plan map of the Niuerwan landslide

图7 牛儿湾滑坡剖面图Fig.7 Profile map of the Niuerwan landslide

2.1 滑源区

滑源区地形坡度约16°,滑体基岩为三叠系巴东组砂质泥岩,表面覆盖厚度约0～9.30 m的第四系残坡积物。滑源区平面形态呈长条带状,宽约150～200 m,长约600 m。滑坡后缘高程约620 m,剪出口高程约450 m,高差170 m,面积约6×104m2;滑体厚度10～30 m,体积约147.5×104m3。滑体从剪出口剪出,整体下错约10～20 m,滑动方向约为N115°E,与岩层产状接近,后因地形影响滑动方向偏转为N81°E,下部基岩发生破碎但仍保持一定的结构性,上部松散且饱水的残坡积层迅速解体形成碎屑流向下流动。

2.2 主堆积区

主堆积区为滑坡堆积体的主要场所之一,范围约7×104m2,滑坡失稳后,约有70×104m3的滑体下滑形成远距离运动,最远运动距离约1400 m,仍有约70×104m3滑体留在滑源区范围内形成堆积。堆积体中可见被破坏的房屋和公路(图1),成分为碎石土,母岩主要为砂质泥岩。该区域已经发生滑动,土体松散,在暴雨情况下,堆积体有发生进一步滑动和形成泥石流的可能。

2.3 碎屑流堆积区

碎屑流沿主滑方向(N81°E)运动过程中,受凸起的山脊地形影响,碎屑流分流成两部分,分别朝N45°E和N81°E方向沿沟谷运动。其中大部分滑体沿N81°E方向运动,受山谷地形变窄影响,滑体动能降低、流速变慢并逐渐堆积,前缘呈现尖灭状;小部分滑体沿N45°E方向运动,越过了凸起的山脊,抛洒于此处形成堆积。整个碎屑流堆积区平面形态近似钳状,沟谷后缘高程约450 m,前缘高程约335 m,高差115 m,长约750 m,面积约6×104m2,最大堆积厚度约20 m。该区堆积物主要为第四系残坡积土以及砂质泥岩的松散混合物,碎石粒径明显比主堆积区中堆积物的粒径小。

3 牛儿湾滑坡运动过程数值模拟

3.1 滑坡模型建立

本次数值模拟采用PFC3D软件,重点在于重现滑体失稳后的运动堆积过程,因此采用PFC3D中的ball-wall模型进行建模,以wall刚性墙体作为滑面,ball颗粒作为滑体,提高了数值模拟的计算效率。

首先,根据牛儿湾滑坡滑后1∶1000地形图生成wall刚性墙体作为计算模型的边界;然后采用generate颗粒生成法生成ball颗粒模拟滑坡岩土体。在构建滑体的过程中,将滑体概化为上下两层的特殊二元结构。滑体总厚度约20 m,上下层颗粒厚度均按10 m生成两层滑体,上层颗粒采用线性接触模型模拟松散土体,下层颗粒则采用平行黏结模型模拟胶结的岩石(图8)。

图8 牛儿湾滑坡PFC3D模型Fig.8 PFC3D model of the Niuerwan landslide

3.2 模型微观参数标定

文中采用虚拟三轴实验进行参数标定,应力-应变曲线如图9所示,为了与模拟时较大的颗粒半径保持一致,将三轴实验试样设置为直径50 m、高度100 m的圆柱形,颗粒数目约12000个,最小半径Rmin为1.2 m,最大半径Rmax为2 m。

图9 PFC3D三轴压缩试验 Fig.9 Virtual triaxial compression test using the PFC3Dsoftware

由图9分析知,对于上层土体试样,由于没有设置粘结强度,三种围压下的峰值强度基本一致;对于下层基岩试样,在低围压下,试样的抗压强度较低,随着围压的增加,试样的抗压强度也随之增大,与实际情况相吻合。

将模拟结果与实际滑坡的运动特征和堆积形态进行反演对比,经过不断的优化和调节微观参数,最终实现模拟滑坡运动与实际运动特征基本一致。由此最终确定牛儿湾滑坡数值模拟的微观参数(表1)。

表1 PFC3D模型微观参数表 Table 1 Micro-parameters for the PFC3Dmodel

3.3 滑坡运动过程及堆积结果

模拟结果显示,滑坡运动总时间约360 s,堆积状态与野外调查结果基本一致,同时得益于数值模拟的便利性,将上下两层滑体分别建模,便于更好地分析上下层滑体的运动堆积特征(图10)。

图10 牛儿湾滑坡滑体分组运动情况图Fig.10 Diagrams showing the movement of the upper and lower layers of the sliding body in different time periods

由图10对上下层滑体进行分析可知,上层滑体为松散的残坡积物,受降雨影响滑体接近饱和,当滑坡失稳后,上层滑体迅速解体形成碎屑流,呈流态化运动(图10a—10d),大部分上层滑体都堆积到了滑坡前缘,运动距离超1 km;而下层滑体为强风化且节理发育的砂质泥岩,同样发生了解体的现象,并与上层滑体混合向下运动(图10e、10f),但下层滑体仍保持有一定的结构性(图10g—10i),并且运动距离较短,仍有大量滑体留在滑源区内。不同时间段内滑坡的运动情况如下。

0～45 s:滑坡发生整体失稳后,滑源区前部滑体从剪出口剪出,开始下滑,主滑方向为N81°E,部分滑体最大运动速度超过15 m/s(图11a、11b),上层滑体发生流态化滑动,下层滑体也逐渐发生解体,两层滑体混合,形成碎屑流向下流动(图10b);由于运动路径上山脊的阻挡,滑体开始产生分流,大部分滑体沿N81°E方向继续运动,小部分滑体则沿N45°E方向运动。

45～135 s:滑体持续运动,分流现象更加明显,同时滑源区中后部滑体失去了前部滑体的阻挡,随之发生更大规模的下滑,但受地形坡度制约,下滑速度相对较小,大部分滑体的运动速度处于4～6 m/s(图11c、11d)。

135～225 s:该阶段滑体的运动更加明显,中部滑体继续以4～10 m/s的速度向下运动,而前部滑体受到逐渐变窄的山谷地形影响,已开始停止运动(图11e、11f);同时流态化特征更加明显,上层滑体持续下滑并裹挟下层滑体不断汇入到堆积区内(图10e、10f)。

225～360 s:前部滑体已停止运动形成堆积体,后方滑体在前方堆积体的阻挡以及阻力作用下也逐渐停止运动(图11g—11i),滑坡的最终堆积形态与实际情况较为一致。

图11 牛儿湾滑坡滑体速度分布图Fig.11 Velocity of the sliding body

3.4 速度曲线及位移轨迹分析

如图8所示,为了监测滑坡的速度和运动轨迹,在滑体前、中、后部的上层滑体、下层滑体分别设置了共26个监测颗粒,对其中19—26号颗粒进行了运动轨迹的监测(图12),同时对滑体前部1—6号、中部7—12号和后部13—18号颗粒进行速度的监测(图13)。

前部滑体的上、下两层滑体都发生了解体,混合形成碎屑流共同运动,1—6号颗粒的运动距离基本一样,约为435～472 m(图12),同时运动速度变化也保持一致,最大速度约18 m/s(图13b);中部滑体的7—12号颗粒的速度变化也基本一致,最大速度约5 m/s(图13c),但上层颗粒(21号)的运动距离约403 m,而下层颗粒(22号)的运动距离仅155 m;后部滑体只有13号颗粒的最大速度达到约5 m/s,14—18号颗粒几乎没有发生运动(图13d),上层颗粒(19号)发生下滑,运动距离约200 m,下层颗粒(20号)运动距离仅72 m。

图12 监测颗粒运动轨迹图Fig.12 Trajectory of specific particles

图13 监测颗粒运动速度曲线图Fig.13 Velocity curves of the monitoring particles

经过以上的分析,发现上下两层滑体的运动特征存在明显的差异,在滑体前部,上下两层滑体的解体破碎现象较为彻底,形成碎屑流,运动特征保持一致;而滑体中后部,上层滑体解体更加彻底,运动距离可以达到下层滑体的2倍,具有明显的流态化运动特征,下层滑体则保持有一定的结构性,呈现整体运动(图12)。滑体前中后部的运动方向也展示出了其特殊的失稳模式,后部滑体沿N112°E～N122°E方向发生滑动,与岩层倾向(N115°E)基本一致,认为后部块体沿真倾向方向滑动;而前部和中部滑体受到滑坡右侧边界山体的影响,滑动方向则偏转为N70°E～N82°E,表现出了沿视倾向方向滑动的特点(图12)。

4 成灾模式探索讨论

牛儿湾滑坡为顺层滑坡,地层倾角较为平缓,滑体运动速度较典型的高位远程滑坡慢,但该滑坡仍产生了远程运动的现象,其远程致灾模式值得进一步的探索与讨论。文章根据滑坡的运动、堆积特征,提出了其远程运动的概化模式,可分为3个阶段(图14)。

图14 牛儿湾滑坡远程运动模式概化图Fig.14 Generalized model for the long-runout movement of the Niuerwan landslide

(1)整体失稳阶段:滑坡沿地层层面发生整体失稳,滑源区范围较大,势动能转换后滑体不同部位的运动速度差异较大,前部滑体动能较大,中后部滑体动能则相对较小。

(2)混合加速阶段:滑体开始解体形成碎屑流,但上下两层性质迥异的滑体展现了不同的运动特征。下层滑体受到基底层摩擦阻力以及上覆土体的影响,加上地形较为平缓,没有强烈的碰撞,保持了一定的结构性;上层滑体为饱和的残坡积土,受到扰动后彻底破碎形成碎屑流,呈现出明显的流态化运动特征,同时,上层滑体裹挟下层滑体顶部较为破碎的块体,逐渐形成混合碎屑流开始加速下滑。

(3)运动流化堆积阶段:随着滑体的进一步运动,混合碎屑流形成远距离运动,持续汇入下方堆积区。

5 结论

文章以牛儿湾滑坡为研究对象,以详细的野外调查资料为基础,对牛儿湾滑坡开展了系统的分析以及数值模拟工作,取得以下结果。

(1)牛儿湾滑坡运动历时360 s,最大速度约18 m/s,最大平均速度约8 m/s,滑坡最大运动距离约1250 m。值得注意的是,远距离运动的物源大部分来自上层滑体,上层滑体的流态化运动特征大大增加了滑坡的致灾范围。

(2)上下两层滑体的运动有明显差异。滑体前部,上下两层滑体的解体破碎现象较为彻底,运动特征基本一致;而滑体中后部,上层滑体解体相对更加彻底,运动距离可以达到下层滑体的2倍,具有明显的流态化运动特征,下层滑体则仍保持有一定的结构性,呈现整体运动。滑体前中后部的运动方向也展示出了其特殊的失稳模式,后部滑体沿真倾向方向滑动,而前部和中部滑体受到地形影响,表现出了沿视倾向方向滑动的特点。

(3)将牛儿湾滑坡的远程运动过程模式分为整体失稳、混合加速和运动流化堆积三个阶段,为研究该类滑坡的远程运动模式做出一定的探索。

(4)强降雨条件和“上土下岩”的二元结构,造成滑坡软弱夹层地下水富集,表层残破积土层饱水,不但降低了岩体的抗剪强度,同时增加了下滑力,是导致滑坡深层失稳整体下滑,表层流化远程运动的关键因素。