基于现场监测数据的某滑坡稳定性有限元分析

韩彪 李凯

摘要：基于某滑坡现场监测数据，文章采用FLAC3D数值分析软件对该滑坡稳定性进行数值模拟，分析了天然工况下滑坡所处状况，确定了滑坡稳定系数。计算结果表明：反演分析得到的计算参数具有实际应用意义，监测资料揭示的滑面与数值分析滑面基本吻合，滑坡稳定系数计算值为1.05，数值模拟值为1.02，均显示滑坡处于欠稳定状态;暴雨工况下，稳定系数为0.96，滑坡由欠稳定状况转变为不稳定状态。

关键词：监测数据;反演分析;数值模拟;稳定性

0 引言

修建高速公路不可避免地要开挖原有地貌，在山体原有平衡被打破后，公路沿线形成大量的临空面，在不良地质条件（如软弱夹层、断层发育）和外在影响因素（如工程开挖、降雨等）的作用下极易形成滑坡灾害。无论滑坡出现在公路沿线还是居民区，一旦发生将严重影响人类生命财产安全。因此，对滑坡稳定性的研究分析一直是工程地质界的热门话题，而对滑坡治理技术的研究则具有更加重要的现实意义。

在滑坡工程领域，计算机数值分析是一种重要的研究和辅助手段，常使用的方法有有限单元法（FEM）、离散单元法（DEM）、边界单元法（BEM）、拉格朗日单元法以及块体理论等[1-6]。滑坡稳定性数值分析方法以有限元法为代表。有限元法能够适用于各种工况，满足较高的精度要求。FLAC3D在岩土工程领域使用最早也是最广泛的。例如，2003年，张雪东等分析某滑坡在天然和蓄水的两种工况下利用FLAC3D进行模拟，指出两种状态下不平衡力曲线的不同[7]。2006年，孙书勤等利用FLAC3D对天台乡基岩顺层滑坡进行三维数值模拟，天然状况下滑坡是稳定的，但是遇到暴雨时会诱发牵引式滑坡[8]。2007年，高圣益等对万州安乐寺滑坡在各种工况下进行流固耦合模拟，得到塑性区和剪应变增量图并对其进行评价，为以后滑坡的分析提供量化依据[9]。张志沛等[10-11]利用数值模拟对比分析了有限元法和极限平衡法的异同，并指出滑坡变形破坏与降雨存在一定的滞后现象。2013年，张建勋对新岩滑坡稳定性数值分析搜索到两个危险滑动面的位置，同时也对位移突变带、塑性区图和剪切应变量分布图进行了分析，为工程的勘察设计提供了参考[12]。基于现场监测数据的某滑坡稳定性有限元分析/韩 彪，李 凯[=JP2]1 工程概况

1.1 工程背景

某滑坡位置所处里程段范围内路线主要以路堑形式通过，该段路堑开挖最高形成6级边坡，高度为57.5 m，其余处路堑边坡为2～5级。2017-11-27开挖至第4级边坡时，坡面出现多条纵向裂缝，2018-01-13开挖至第3级边坡时，距堑顶约30 m的后缘出现大范围的弧形错台，第3级边坡坡面有明显剪出迹象。前缘反压之后，通过监测发现2018-01-25之后边坡变形速率明显减缓并逐渐趋于稳定，堑顶的变形量得到控制。

1.2 监测剖面及监测数据分析

本文采用反演法分析滑坡稳定性，所用监测孔ZK2、ZK3、ZK4位于4-4剖面，监测孔ZK6、ZK7、ZK8位于3-3剖面，3-3剖面和4-4剖面布置如图1所示（由于篇幅有限仅给出4-4剖面监测数据，3-3剖面数据不做展示）。

3-3剖面监测数据揭示出该剖面处滑坡中部和后缘处于蠕动变形阶段，累计变形量较小，变形范围为±10 mm。滑坡前缘累计监测变形持续增加，最大监测位移出现在3.5 m深处，为30 mm，前缘先于后缘变形，说明该滑坡为牵引式滑坡。

由图2～4可知：4-4剖面处滑坡活动明显较为剧烈，剖面前缘处于稳定状态，后缘12.5 m深处出现最大80 mm的垂直路线方向位移，且剖面中部一直处于变形挤压阶段，累计变形量仍有增长趋势。监测数据确定滑坡后缘经一次大的变形后，下滑趋势逐渐减弱。

结合两断面监测数据可知：4-4剖面滑坡活动较为剧烈，滑面深度为9.0～12.5 m，滑坡性质属牵引式滑坡。由监测数据推断出的滑面位置如图5所示。

2 数值模拟

2.1模型建立

计算模型选取滑坡4-4剖面，采用ANSYS软件生成FLAC3D识别的模型单元文件和节点文件。为精确计算，对研究滑坡滑动带进行网格细分，模型长285 m，宽30 m，高度为140 m。边界条件为下部固定约束，左右两侧法向约束，上部为自由边界。所建模型如图6所示。

2.2 参数选取

结合原位试验数据以及位移监测数据，获得模型岩土体参数，如表1所示。

2.3 计算结果分析

数值模拟位移计算值与监测位移结果对比如表2所示。

对比4-4剖面滑带处实际监测位移与数值模拟计算值，ZK2处实际监测位移值为20 mm，数值计算值为27 mm;ZK3处实际监测值为8 mm，数值计算值为22 mm;ZK4处实际监测值为32 mm，数值计算值为37 mm。其中ZK3处计算值误差较大，为14 mm，但并未影响滑坡整体趋势。数值计算得到的稳定系数为1.02，勘查资料给出的该滑坡天然工况下稳定系数为1.05。由《滑坡防治工程勘察规范》（DZT0218-2006）判断滑坡处于欠稳定状态。

数值计算结果与实际监测值基本吻合，即数值计算结果可反映边坡实际情况。

2.4 天然工况模拟结果分析

通过数值模拟，得到该滑坡天然工况下稳定性相关结果，如图7所示。

图7为天然工况下滑坡稳定性分析结果图，计算所得滑坡稳定系数为1.02，滑坡处于欠稳定状况，分析X方向位移云图（即垂直路线方向），滑带处位移量与三个监测孔处位移量大致吻合，但位移云图并未显示出明显的滑带位置，说明天然工况下滑坡变形并不剧烈，符合滑面尚处于欠稳定状态的实际状况。

由剪应变增量图可以看出：研究滑坡以监测孔ZK2所处的滑坡中部为分界，最大剪应变增量分布并不连续，自滑坡前缘至后缘呈现出先降低后增大的变化规律（整体剪应变增量最大值并未贯穿整个坡体）。剪应变增量最大值分布于滑坡前緣，说明该区域变形速度最大，坡体易沿此区域发生破坏，由模拟结果判断坡体存在三种破坏形式：（1）滑体前-中-后三个区域滑带出现连通之后发生整体滑动;（2）前部滑体以ZK2处堑顶为滑坡后壁，形成次级滑坡;（3）上部滑坡整体从ZK3监测孔上部平台剪出，形成上部滑坡。并未贯通整个坡体的最大剪应变增量分布表明，坡体在该工况下尚未发生整体下滑，对应滑坡处于欠稳定状态。

塑性区模拟结果可知：滑坡前、后部下部滑带为正在发生剪切部分，上部滑体多为过去剪切部分。滑坡后缘为过去发生剪切滑体，该部分位于正在发生剪切变形的滑带上，有向下滑动的先决条件。受后缘滑体挤压，中部滑体可分为上下两层，其中下层部分（长约90 m）为过去剪切-过去张拉带，其将正在剪切-正在张拉带分隔为前后两段，上部为正在发生剪切变形的坡体，占据了滑体1/4体积，说明中部滑体正在发生剪切变形的部分逐步上移，中部滑带部分不再发生剪切破坏，由剪切带向抗滑带演变。

由分析结果可知：天然工况下研究滑坡有三种破坏形式，即整体下滑、上部剪出下滑和次级前缘滑坡下滑。塑性区显示虽然整体下滑趋势并不明显，但滑坡中部抗滑段滑体受到明显挤压，滑坡稳定系数为1.02，滑坡处于欠稳定状态。

2.5 暴雨工况结果分析

该滑坡滑面由泥质粉砂岩与泥岩交界处的软弱带组成，滑床部分为风化程度较低的泥质粉砂岩，经暴雨浸湿后滑带强度降低，依据实际情况对泥岩和含砾粉质黏土参数做出调整。暴雨工况下滑坡数值计算结果如图8所示。

通过暴雨工况下的计算结果可知，垂直路线方向位移最大值出现在坡体前缘，整个坡体被滑带贯穿，滑体内部整体位移变化趋势呈现出前缘变形大于后缘变形的情况，前缘位移增长速率大于后缘增长速率。滑带处平均位移量为4～5 cm，相比天然工况下坡体前缘监测孔ZK4位移处（变形量为26 mm），其量值增长了一倍。说明受雨水的影响，滑带处变形量显著上升，滑坡整体位移变形由天然工况的前缘大-中部小-后缘大转变为前缘最大-中部小-后缘最小的变化规律。位移量与位移分布规律的变化均显示滑带连续，滑坡由不稳定状态向欠稳定状态过度。

对比图8（c）与图7（c），由于滑带土受雨水影响，滑体内正在剪切-过去剪切部位向下延伸至滑带区域处，与滑床相接，逐渐取代了原有的过去发生剪切-过去发生张拉部位，原本占滑带部位约90 m的过去剪切-过去张拉段缩减为20 m。相比之下正在剪切-过去剪切和正在剪切-过去剪切-过去张拉两部分总共占据了整个坡体的3/5，且从滑坡中部向两端不断延伸。

受降雨影响，滑坡前缘位移值增长了一倍，原本并不连贯的滑带发生连贯，最易发生剪切破坏的区域已贯穿整个坡体。整个坡体有3/5处于正在剪切变形阶段，该部分集中分布于滑带中部，仍不断向滑坡前缘延伸，表明坡体中部应力环境已经由天然工况下的挤压剪切阶段逐步过渡到纯剪切阶段，即实际中由蠕滑阶段向剧滑阶段过渡。经计算滑坡稳定系数为0.96，相比天然工况稳定系数降低5.8%，由欠稳定状况转变为不稳定状态。

3 结语

本文借助FLAC3D软件，反演分析并对比滑坡天然工况与数值模拟的结果，得到了暴雨工况下滑坡所处的状况，通过分析垂直路线方向位移云图、最大剪应变增量云图和塑性区分布，得到各工况下位移变化规律和剪应变增量变化规律，并计算不同工况下稳定系数。通过上述研究工作，可以得到以下结论：

（1）通过分析监测资料，4-4剖面后缘经过一次大的变形后，滑坡中部处于挤压变形阶段。监测数据揭露出的滑面位置分布于泥岩与粉质泥砂岩交界的软弱夹层处，该软弱夹层强度低，浸水易软化，滑体极易沿着该交界面发生滑动，确定出滑面深度为9.0～12.5 m。

（2）天然工况下，滑体内部滑带并未贯穿整个滑体，滑坡稳定系数为1.02，存在三种破坏形式：①滑带贯通后，滑体整体滑动;②前部滑体以ZK2处堑顶为滑坡后壁，形成次级滑坡;③上部滑坡整体从ZK2监测孔上部平台剪出，形成上部滑坡。

（3）暴雨工况下，稳定系数为0.96，相比天然工况稳定系数降低5.8%，滑坡由欠稳定状况转变为不稳定状态。监测数据与数值计算值对比结果显示，滑坡模型的建立是成功的，模拟计算得到的结果可以反映出滑坡所处的实际状况。

参考文献：

[1]Sarada Prasad Pradhan，Tariq Siddique.Stability assessment of landslide-prone road cut rock slopes in Himalayan terrain：A finite element method based approach[J].Journal of Rock Mechanics and Geotechnical Engineering，2020，12（1）：59-73.

[2]VikramGupta，Rajinder K.Bhasin，Amir M.Kaynia，et al.Finite element analysis of failed slope by shear strength reduction technique： a case study for Surabhi Resort Landslide，Mussoorie township，Garhwal Himalaya[J]. Taylor & Francis，2016，7（5）：1 677-1 690.

[3]ZHUYao，ZHU Chengguang，LIU Chun，et al.A Discrete Element Analysis of the Sliding Friction Heat in High-Speed and Long-Runout Landslides[J].Earthquake Research in China，2020（1）：103-116.

[4]L.Najarzadeh，B.Movahedian，M.Azhari.Stability analysis of the thin plates with arbitrary shapes subjected to non-uniform stress fields using boundary element and radial integration methods[J].Engineering Analysis with Boundary Elements，2018（87）：111-121.

[5]Mehdi Panji，Hamid Koohsari，Mohammad Adampira，et al.Stability analysis of shallow tunnels subjected to eccentric loads by a boundary element method[J].Journal of Rock Mechanics and Geotechnical Engineering，2016，8（4）：480-488.

[6]Rcyes S F，Deere D U.Elastic-plastic Analysis of Underground Openings by the Finite Element Method.[D].USA：University of Illinois at Urbana-Champaign，1966.

[7]张雪东，陈剑平，黄润秋，等.呷爬滑坡稳定性的FLAC3D数值模拟分析[J].岩土力学，2003（24）：13-116.

[8]孙书勤，黄润秋，丁秀美.天台乡滑坡特征及稳定性的FLAC3D分析[J].水土保持研究，2006，13（5）：30-32.

[9]高圣益，魏学勇，周 晃.基于FLAC3D的滑坡稳定性数值模拟分析[J].长江工程职业技术学院学报，2007（4）：1-4.

[10]張志沛，高旭和.基于MIDAS/GTS的黄土边坡稳定性分析方法及应用[J].水利与建筑工程学报，2016，14（2）：182-184，217.

[11]张志沛，武素琴，彭 惠.芜铜高速公路上马村滑坡稳定性分析与研究[J].路基工程，2011（3）：227-230.

[12]张建勋，陈福全，简洪钮.被动桩中土拱效应问题的数值分析[J].岩土力学，2004，25（2）：173-178.