有限差分法与极限平衡法在某滑坡稳定性研究中的应用

高 凯,王林维,范新宇,康文军,贾志献

(中国电建集团西北勘测设计研究院有限公司，西安 710065)

有限差分法与极限平衡法在某滑坡稳定性研究中的应用

高 凯,王林维,范新宇,康文军,贾志献

(中国电建集团西北勘测设计研究院有限公司，西安 710065)

基于云南某水电站滑坡的地质勘察及设计资料，通过有限差分和极限平衡2种方法对滑坡进行对比分析，达到2种方法的优势互补，客观模拟了该滑坡变形特征及潜在破坏情况，为边坡开挖设计及滑坡治理提供了科学依据与技术指导。关键词：滑坡稳定性；有限差分法；极限平衡法；FLAC3D

0 前 言

中国幅员辽阔，地质条件复杂多样，滑坡、泥石流等灾害频发。结合不同区域的工程地质条件，采用一定的方法和技术手段分析工程滑坡稳定性，并据此采取相应的工程对策，对安全、经济、高效的滑坡治理工作具有重要意义。

目前滑坡稳定性分析的主要方法可分为基于强度折减技术的数值分析方法和基于极限平衡理论的各类条分法2类。为了准确真实地反映某电站滑坡目前的实际情况及稳定性，本文在地质调查的基础上，采用有限差分[1]及极限平衡[2-3]2种理论方法对其进行模拟计算分析。通过数值模拟计算，预测滑坡区的稳定性及潜在破坏情况，为边坡开挖设计及滑坡治理提供科学依据。

1 工程概况

本文滑坡位于澜沧江某水电站右岸，以土质为主，主要为含砾黏土和坡积块碎石土，其中滑带土厚12～20 m。目前滑坡上部路面错断并下沉3～4.2 m，高程1 904.00～1 855.00 m之间的边坡拉裂并向外挤出，整个滑坡体长约250～300 m，宽约130 m，高差约50 m，初步估算该滑坡体体积约60万m2，目前1 855.00～1 904.00 m之间的边坡已全部解体(见图1)。

该滑坡区下部为当地居民出入及工程施工的必经之路，因此有必要准确查明滑坡的发育规模、坡体结构、形成机制，研究其滑带土的物理力学特性，计算各种工况下的稳定性，以便对其进行合理的预防和治理。

图1 电站滑坡全貌

2 岩土体物理力学参数取值

通过对本次物理力学试验数据的整理分析，结合有关资料确定该滑坡岩土体的抗剪强度指标，其物理参数建议取值见表1。

表1 天然状态下计算参数表

饱和状态下C、φ折减系数为0.8，滑坡所属区域地震基本烈度为Ⅶ度，相应的地震动峰值加速度为0.10g。

3 有限差分法数值分析

3.1 计算原理

FLAC3D是美国Itasca公司开发的三维显式有限差分法程序，它可以模拟岩土或其它材料的三维力学行为。FLAC3D软件的基本原理[4-5]是拉格朗日差分法，它是一种利用拖带坐标系分析大变形问题的数值方法，并利用差分格式按时步积分求解，适合于解决岩土工程中经常遇到的大变形问题，是一种理想的岩土工程计算软件。本文将通过该软件对滑坡整体计算分析[6-13]，重点研究滑坡区所在的位移场、应力场和塑性区特征，评价滑坡在自重条件下的稳定性，探讨其潜在破坏情况。

3.2 模型建立

(1) 实体模型

本文计算区域边界按照滑坡区的地形，并满足一定的边界效应要求，在平面上(X、Y)取586 m×537 m的区域，高程(Z)范围为1 600.00～2 200.00 m。

地质勘察分析表明，最危险的部位位于滑坡中的滑带土，考虑到计算条件因素并结合实际情况，将该滑坡模型简化，模型上部为坡积黏土层，下伏冲洪积碎石土层为绢云母石英千枚岩，其中滑带土位于模型中上部，黏土层与碎石土层之间。

(2) 计算模型

根据实际情况对重点研究区域进行网格加密处理，采用四面体网格离散计算模型，经过接口程序将模型转换后，将实体模型离散为94 559个单元和20 930个节点；至此，实体模型已成功导入FLAC3D，并可以进行三维有限差分计算。

3.3 边界条件及物理力学参数

将外部山体上表面设为自由边界，模型底部设为固定约束边界，四周为单向边界。计算使用的物理力学参数见表1。

3.4 计算方法及步骤

计算时，按下述步骤进行：首先，选择弹性本构模型，按前述约束条件，只考虑重力，进行弹性求解，计算至平衡后对位移场和速度场清零，生成初始应力场；在此基础上，对滑坡区各地层赋表1参数值，选择Mohr-Coulomb模型进行弹塑性求解，直至系统达到平衡。

3.5 计算结果分析

将电站滑坡模型计算结果进行后处理，绘制滑坡区的整体和典型剖面(Y=260 m)位移及主应力等值线云图，综合分析滑坡的位移场、应力场以及塑性区的分布规律及特征。典型剖面见图2。

图2 典型剖面位置图

(1) 位移场分析

由滑坡区整体位移云图(见图3)可知位移变形较大区域集中于滑坡中后部，最大位移可达700 mm，其余部位位移主要在0～50 mm之间。表明在有限范围内位移增长极为迅速，为滑坡位移的转折部位，同时滑坡的潜在破坏以浅表层圆弧形剪切破坏为主。另外，滑坡体后缘上部岩土体位移值较大，说明滑坡之后对上部土层产生较大影响，使其产生较大变形。

图3 整体位移云图

(2) 应力场分析

从典型剖面最大主应力云图(见图4)分析可知，滑坡区所受应力均为负值，主应力等值线基本顺着坡面方向，并一直延伸到坡脚，但在上表黏土层底部发生跳跃，说明该部位已发生明显的应力集中，这对滑坡区稳定性不利。而滑坡区内部，最大主应力方向与水平轴的夹角逐步变大，表明滑坡区深部岩土体主要受垂直方向的压应力作用，体现为受压屈服，且在上部黏土层内发生应力集中，对滑坡区稳定性不利。

图4 整体最大主应力云图

(3) 塑性区分析

为定量分析电站滑坡区的稳定性，采用FLAC3D内置的solve fos命令计算滑坡的安全系数Fs。由于滑坡区内部表现为受压屈服，其破坏形式为剪切破坏。经计算，其安全系数Fs=0.95，表明其在天然状态欠稳定。由图6可知滑坡受剪切作用较大的区域位于上表黏土层及碎块石土层；该部位等值线较为密集，局部土体塑性区已贯通破坏，对滑坡的稳定性极为不利。

图5 整体最小主应力云图

图6 整体剪应力云图

4 极限平衡法分析

本文继续对滑坡区采用二维极限平衡法计算分析，以便更好地找出滑面所在位置及分析验证该滑坡稳定性。土层结构的物理力学参数取值见表1，其破坏模式遵循摩尔-库仑强度准则。

4.1 计算原理

采用瑞典条分法原理，即：假定滑裂面是圆弧，将滑裂体划分为若干刚性竖直土条，且不考虑土条两侧的作用力，通过滑动土体的整体力矩平衡来计算稳定系数。其中土条底部抗剪力τR,i满足极限平衡条件，可根据莫尔-库仑准则和边坡稳定系数Fs确定，即：

(1)

式中：Ci和φi分别为第i土条底部的抗剪强度指标；σi是该土条底部的法向应力。因此，如不考虑孔隙水压力，瑞典条分法的稳定系数Fs可写作：

(2)

式中：Wi、li和αi分别是第i土条的总重量、底部长度和底部倾角。

4.2 计算剖面的选择

计算使用的断面为1-1、2-2及3-3作为典型剖面(见图7)赋值计算，阴影部分为滑坡范围。

图7 滑坡区平面地质简图

4.3 计算结果及分析

计算方法采用简化Bishop法，分别计算1-1～3-3剖面在天然、地震及暴雨工况下的稳定性，由于计算图片较多，本文只附2-2剖面(图8)相关图片。

图8 2-2剖面图

通过所建立的剖面模型分别对模型不同材料赋值，采用圆弧法对潜在最易失稳部分进行搜索，并计算对应不同工况下的最小稳定系数。经计算，电站滑坡区2-2剖面计算成果见图9～11，3条剖面的稳定系数结果见表2。

表2 危险面稳定计算成果表

图9 2-2剖面天然工况图

图10 2-2剖面暴雨工况图

图11 2-2剖面地震工况图

从计算结果看，天然工况下滑坡区3条剖面稳定系数介于0.9～1.1之间，整体处于极限平衡状态，2-2剖面位置已经出现局部破坏；暴雨及地震工况稳定系数均小于1.0，结合图片表明此2种工况下滑坡区已发生浅层破坏；且在滑坡后缘仍然存在不稳定的圆弧破坏面，在恶劣条件及人为因素干扰下，该部位极有可能再次发生拉裂、变形、塌滑现象。

5 2种方法比较分析

从上述2种方法的计算过程及结果可知：

(1) 有限差分法和极限平衡法计算出的安全系数较为接近，滑动面位置也近乎相同，证明这2种方法在安全系数的定义方面具有一致性，在工程实际中可以相互比较及验证。

(2) 有限差分法是通过建立三维模型，对不同地层区分赋参数计算，不仅可以得出某个剖面的安全系数，还可以得出滑坡区整体的稳定性，结果更符合实际情况。

(3) 有限差分法得出的滑动面是一条滑动带，而极限平衡法得出是一个圆弧，滑动带是坡体中塑性破坏区域，并不是严格的滑动面，但这符合实际情况，极限平衡法得出的滑动面并不是真实的滑动面。

(4) 有限差分法不仅能得出合理的滑动面和安全系数，而且能得出滑坡区的位移及应力应变图，可以清楚地看到坡体的下沉和坡脚的隆起，而极限平衡法不能做此分析。

(5) 有限差分法建立模型等前处理工作较为复杂，计算分析速度较慢，而极限平衡法建模相对简单，分析速度较快，且可以同时计算多条滑面的安全系数。

综上所述，有限差分法与极限平衡法在计算滑坡稳定性方面各有利弊，但整体来说前者更具有优越性，是当今岩土体稳定性计算的趋势。

6 结 语

(1) 滑坡是影响工程建设的重要因素之一，及时发现及治理是工程顺利进行的必要条件。

(2) 通过三维模型的数值分析可知，电站滑坡区位移变形较大区域位于上表黏土层及碎块石土层；该部位等值线较为密集，局部土体塑性区已贯通破坏；发生滑坡之后，由于土体下错、塌滑，导致滑坡体后缘过陡，进而影响到上部土体受力状态，计算结果显示1 930 .00 m高程以上位移变形也较大，有应力集中现象；在极限平衡法验证过程中，此次滑坡范围土体由于应力释放暂时趋于稳定状态，反而1 930.00 m高程以上后缘土体出现危险滑面，稳定系数较低，极有可能再次发生拉裂、变形、塌滑现象。

(3) FLAC3D在计算边坡的稳定安全系数方面的精确性与极限平衡法相当。与极限平衡法相比,FLAC3D考虑了岩体内部的应力应变关系,可以直观地得出岩体内各点的应力应变、位移曲线图和滑动面位置,分析边坡破坏的发生和发展过程,说明FLAC3D在边坡稳定分析中比极限平衡法更具有优越性。

(4) 2种方法计算结果基本一致，反映了滑坡变形破坏的特点，计算结果可以作为滑坡治理的依据。建议对该滑坡进行治理，以保证当地居民及作业人员生命财产安全。

[1] 黄润秋,许强.显式拉格朗日差分分析在岩石边坡工程中的应用[J].岩石力学与工程学报，1995(04)：346-354.

[2] 梁素韬,洪海,张红梅.基于极限平衡法的折线形滑坡体稳定性分析[J]. 水利水电技术, 2013(07)：48-50.

[3] 丁美林.极限平衡法在公路边坡稳定性分析中的应用研究[J].青年科学,2010.

[4] 孙书伟,林杭,任连伟,等. FLAC3D 在岩土工程中的应用[M].北京：中国水利水电出版社，2011.

[5] 陈育民.FLAC /FLAC3D 基础与工程实例[M].北京：中国水利水电出版社，2009.

[6] 李良平,胡伏生.鄂尔多斯盆地白垩系三维地质建模研究[J].西北地质,2007,40(02)：109-112.

[7] 许国,李敦仁.Gocad地质三维建模技术及其在水电工程中的应用[J].红水河,2007(26)：113-116.

[8] 梁晓波,张为法.厄瓜多尔德尔西水电站右坝肩边坡稳定性模拟分析[J].西北水电,2014(04):19-22.

[9] ansys建模分网指南[M].北京：中国水利水电出版社.

[10] 王向东,文江泉.用FLac3d进行土质高边坡稳定性分析[J].西华大学学报,2005,24(03)：87-89.

[11] 李旭东.FLac3d在边坡稳定性分析中的应用[J].中国水运,2008,8(04):77-79.

[12] 乔国文,王运生.西南某电站右岸开挖边坡稳定性的FLac3d分析[J].工程地质学报,2014,12(03)280-284.

[13] 陈祖煜.土质边坡稳定分析-原理、方法、程序[M].北京;中国水利水电出版社,2003.

Application of Limit Equilibrium Method and Finite Differential Method in Study of Landslide Stability

GAO Kai, WANG Linwei, FAN Xinyu, KANG Wenjun, JIA Zhixian

(Northwest Engineering Corporation Limited, Xi'an 710065, China)

In accordance with the geological investigation data on landslide and the design data, the deformation characteristics and potential failure of the landslide are objectively simulated through analysis and comparison on the landslide by application of limit equilibrium method and finite differential method. This provide the slope excavation design and landslide treatment with scientific basis and technical guidance.Key words:landslide stability; limit equilibrium method; finite differential method; FLAC3D

1006—2610(2016)06—0022—05

2016-09-29

高凯(1987- )，男，山东省滨州市人，助理工程师，从事水利水电工程设计工作.

TV698.232

A

10.3969/j.issn.1006-2610.2016.06.006