基于FLAC^3D 的滑坡稳定性数值模拟分析

飞鱼游鹰

第24卷第4期
2007年12月
长江工程职业技术学院学报
Journal of Changjiang Engine~ng Vocational College
Vo1．24 NO．4
Dec．2OO7
基于FIAC3。的滑坡稳定性数值模拟分析
高圣益 ，魏学勇2，周晃3
(1．长江空间信息技术工程有限公司，武汉430019；2．中国地震局地壳应力研究所，北京100085；
3．重庆市巫山县国土资源局，重庆404700)
摘要：在介绍R AC， 基本特点的基础上，利用R A 软件对万州安乐寺滑坡在水库不同蓄水位条件下的滑坡
稳定性进行了流—伺耦合模拟分析研究。通过建立数值计算模型和选择合理的参数，模拟得出滑坡体的塑性区分
布特征和剪应变增量集中带的发展趋势。运用三维耦合分析，对水库不同蓄水水位条件下的滑坡稳定性做出评
价，为滑坡的治理提供了必要的理论依据。
关键词：R AC ；滑坡稳定性；流一固耦合分析
中图分类号：TU457 文献标识码：A 文章编号：1673—0496(2007)04 —0001—04
Simulation Analysis of Landslide Stability Based on FLAC～
' '
GAO Sheng—yi ，W EI Xue—yong ，ZH0U Huang
(1，Yangtze Space Information Technology Engineering Ltd．，Wuhan 430019，China；
2，Crusml Stress Institute of China Earthquake Bureau，Beijing 100085，China；
3，Chbngqing Wushan Land Resources Bureau，Chongqing 404700，China)
Abstract：Some fluid—solid coupling simulation an alysis was made on the landslide stability of Anle Temple in
W an zhou under vadable water levels of the resen，oir with the software of FLAC3D
． trhe distribufion characteristics
of plastic zone an d trend of shear strain increment concentration of sliding land were simulated out by setting up
data model an d selected parameters．In the end，the landslide stability Was apprmsed according to the reservoir
waterlevelanditwilldo somejobinthelandprotectionafterwards．
Key words：vI．~C ；stability of landslide；fluid—solid coupling analysis
近年来，各种数值模拟技术在岩土力学中有了很大的
发展和广泛的应用。然而，这些数值分析方法其理论本身
以及采用的算法都有各自的局限性。例如有限元和边界元
都有小变形的假设，且需要大量的内存。近年来发展起来
的快速拉格朗日分析(Fast Lagrangian Analysis ofContinua，
简称FLAC)是在较好吸取上述方法的优点和克服其缺点基
础上形成的一种新型数值分析方法。FLAC~ 是美国Itasca
Consulting Group，Inc．公司为地质工程应用而开发的基于拉
格朗日差分法的一种三维显式有限差分程序，它不仅适宜
于处理大尺度、大变形工程和地质问题，而且可以在初始
模型中加入诸如断裂、节理构造等地质因素。目前，FLAC~
软件已经广泛应用于工程地质、岩土力学以及构造地质学
收稿日期：20o7__08___23
作者简介：高圣益(1964一 )，男，湖南人，高级工程
师，大学，主要从事地质勘探研究。
和矿山工程领域。本文以三峡工程万州库区安乐寺滑坡为
例， 阐述了FLAc3D在滑坡稳定性流一固耦合分析中的应
用。
1 FLAC 的解析特点和计算流程
1．1 解析特点
无论是静力还是动力问题，三维快速拉格朗日分析都
利用动态的运动方程进行求解，这使得三维快速拉格朗日
分析很容易模拟动态问题，如振动、失稳、大变形等。同
有限元相比，FLAC~ 具有以下几个特点：
(1)FLA 基本原理类同于离散单元法，但它能像有
限元法那样适用于多种材料模式与边界条件的非规则区域
的连续问题求解。
(2)该程序采用了离散元的动态松弛法，不需要求解
大型联立方程组(刚度矩阵)，便于在微机上实现。
(3)基于拉格朗日算法，适合模拟大变形，且能模拟
维普资讯 http://www.cqvip.com
高圣益等 基Y-FLAC~。的滑坡稳定性数值模拟分析
岩体沿某-软弱面产生的滑动变形。
(4)能针对不同材料特性，使用相应的本构方程来比
较真实地反映实际材料的动态行为。
当然， 由于算法自身的原因也造成了三维快速拉格朗
日分析的一些固有缺陷，与其他方法相比，快速拉格朗日
分析的计算效率较低，但随着计算机运算速度的加快，这
种缺憾已经微不足道。
l-2 计算流程
FLAC~ 的计算过程如图l所示。 ·
图l FLAC3D分析方法的计算循环
在每个时步内按下列顺序来进行计算：
(1)由节点速度求出新的应变速度。
(2)运用结构方程式由应变速度和应力计算出各时步
内的新的应力。
(3)运用运动方程式由应力和失衡力算出新的质点速
度和位移。
在各个时步内进行卜述循环计算。当计算达到平衡状
态时，失衡力趋近于0。如果由塑性流动产生，则失衡力趋
近于某一数值。
2 FLAC' 内置基本单元体和本构模型
2．1 网格生成器
FLAC～中用于计算的有限差分网格限定了待研究分析
问题的实体域。最小的可能的网格可以只有一个单元构
成，但对于大多数问题，用于计算的网格需要由成千上万
个单元组成。FLAC" 具有一个功能强大的网格生成器，内
置有l2种基本形状的单元体可供选择。利用这些基本单元
体，几乎可以构成任何形状的空间 体模型。
在FLAC" 中，这l2种内置的基本形状单元体的关键宁
包括brick，wedge，tetra，pyramid， cylinder， dbrick，
radbrick，radtunnel，radcylinder，cshell，cyllnt，tunlnt。这
些基本单元体是构成评价或预测模型的基础， 岂们可以单
独用于一个模型，也可以联合用于一个更加复杂的模型。．
2．2 本构模型
FLAC～主要是为地质工程应用而开发的岩土体力学数
值评价计算程序， 自身设计有九种材料本构模型： 空模
型、弹性各向同性材料模型、弹性各向异性材料模型、德
拉克一普拉格弹性材料模型、莫尔一库伦弹塑性材料模型
等。本构模型反映了材料变形的变形特征和破坏机制，在
计算过程中，它对模拟计算的精度和准确件都有较大的影
响。因此，在实际应用中应根据研究对象和实际材料的变
一2一
形特性选用相应的本构模型。
3 流一固耦合分析方法
FLAC～基于多孔连续介质的渗流理论，模拟流体存孔
隙介质中的流动及其相互的耦合作用，采用显式有限差分
法来求解渗流运动方程。其基本假定为：(1)渗流运动服从
达西定律；(2)流体和固体均被当作是可变形体；(3)固体应
力状态遵循太沙基有效应力原理；(4)流一固耦合过程满足
Biot方程。当流体在多孔介质中流动时，主要引起孔隙压
力、饱和状态和渗透流量等变量的改变，这些变量通过流
体质点平衡方程与Darcy定律间的相互关系来描述流体的流
动。
在任何计算中，边界条件是很重要的，它的好坏直接
影响到计算的结果。在FLAC~ 中，有4种类型的边界条
件，它们分别是：①给定孔隙水压力；②给定边界外法线
方向流量分量；③不透水边界(不透水边界程序中默认)。④
透水边界，透水边界采用如下形式给出：
： 一 ) (1)
．
～
式中， 是边界外法线方向流量分量；k是渗透系数
(m／s)； 为水的容重；P是边界处的孔隙水压力； 是渗
流出口处的孔隙水压力。其他边界条件，如部分边界上的
常压力和常流量，则需通过程序内嵌的FISH语言来实现。
4 滑坡稳定陛流一固耦合分析
4．1 滑坡概况
滑坡区位于苎溪河北岸，为典型的中低山切割丘陵河
谷地貌。地形起伏大，北高南低，东西向地形起伏不大。
地面高程在140~220m之间，滑坡区由二级台阶组成：第一
级台阶断续分布， 高程在l77～l85m之 ， 台阶最宽
130m； 第二级台阶发育完整， 分布连续， 高程205～
220m，宽120~250m，地形平缓，坡角2。～5。，后壁较
陡。一、二级台阶之间的斜坡倾角为20。～40。。滑坡区内
覆盖层厚1．65~38．00m，由人工填土、粉质粘土夹块石、淤
泥质粉质粘土、粉砂土及下伏基岩泥岩央砂岩组成。人T
填土厚0．50～3．10m，粉砂土厚1．70～ 13．70m，均属软弱
土；粉质粘土灾块石(厚1．65~34．40m)属巾软土、基岩为坚
硬土。由于滑坡体巨大，加之三峡水库蓄水，使滑坡地质
条件进一步恶化，安全系数进一步降低。冈此，对该滑坡
体的变形情况和稳定性做出评价具有重要意义。
4．2 计算模型和边界条件
在滑坡勘察数据资料的基础上，利用FLAC~ 提供的网
格生成器建立滑坡稳定性流一固耦合分析的计算模型，如
图2所示。计算模型范围，垂直河流方向(X轴方向，指向
河流为正)宽360m：顺河方向(Y轴方向，指向河流上游
为正)取l0米；铅直方向(Z轴方向，向上为正)由地表取
一8Om高程。
应力边界条件为左、右两侧采用水平约束，底面采用
双向约束：渗流边界条件为模型底面假定为不透水边界，
维普资讯 http://www.cqvip.com
2007年12月 长江工程职业技术学院学报 第24卷第4期
图2 评价模型网格剖分图
左侧和右侧边界为给定孔隙水压边界，水位以上为自由透
水面。初始应力场为自重应力场，初始渗流场根据钻孔地
下水位资料选取，其中，左侧边界水位为210m。考虑到后
缘边界距库区较远，因而计算中左侧边界的水位保持210m
不变，右侧边界水位由水库不同蓄水位确定。
4．3 力学模型与力学#,-It的确定
岩石当荷载达到屈服强度后将发生破坏、弱化，应属
于弹塑性体。在FtACe 中，对于弹塑性材料，其破坏判据
准则有两个：德拉克— 普拉格准则和莫尔一库仑准则。本
项研究选择莫尔—库仑准则。
计算模型中采用的各岩层力学参数包括弹模、泊松
比、内聚力、内摩擦角、抗拉强度和密度等参数，见表1。
表1 计算模型中各种介质力学参数
粤 体积模量剪切模量内聚力内摩擦角密度抗拉强度渗透系数
(MPa) ( a) (MPa) (o ) (kg／m') ( a) 10"
4．4 计算工况
根据三峡水库不同蓄水高程，确定以下几种计算工
况。工况1：135m~位；工况2：145mzk位；工况3：156m
水位；工况4：175rI冰位。
(a)135mTge．
(C)156m~位
模拟计算思路是：首先模拟滑坡的自重应力场和初始
渗流场，将历史位移初始化为零，再将水位逐步抬升到
135、145、156和175m。通过耦合分析可以得到每一计算工
况的塑性区和剪应变增量集中带分布。最后对滑坡的变形
和稳定性做出综合评价和预测。
5 模拟结果分析
利mFUAC" 程序计算了不同工况条件下的滑坡体的塑
性区、剪应变增量集中带分布以及它们的发展趋势，通过
对塑性区、剪应变增量分布进行综合分析，从而对不同工
况条件下的滑坡稳定性做出评价。
5．1 塑性区分析
图3给出了不同水位高程下滑坡体的塑性区分布。由图
可知，当水库蓄水至135m时，滑坡体的后缘和坡脚处发生
了部分单元的塑性破坏。其中，滑体后缘塑性区以拉剪破
坏为主，坡脚处主要为剪切破坏。滑坡体的后缘和坡脚的
塑性破坏区没有贯通，因此，在135m蓄水条件下，滑坡处
于稳定状态。随着水库水位的抬升，滑坡体的塑性区也随
之扩大。从图3可以看出，当蓄水至145m1~，塑性区分布在
滑坡前缘被水体淹没的部位有较大幅度的增加，滑体塑性
区面积大致占整个滑坡体面积的30％，以拉伸破坏为主。
此时，由于水位较高以及渗流的影响，滑坡前缘可能会形
成局部滑移。
当三峡水库蓄水至156m时，滑坡体的变形趋势发生了
明显的变化，塑性区范围也有大幅度的增加。此时，滑坡
体塑性区面积大约是整个滑坡体面积的50％，较145m7]~位
时塑性区面积增加了近25％。从塑性区连续性看，滑坡体
塑性破坏基本贯通，滑坡处于极限平衡状态。
水库水位抬升到175mN，滑坡体的塑性区已经完全贯
通，其面积大致占整个滑坡体面积的75％。此时，滑坡在
受到外力(如暴雨)的作用下，将发生滑动破坏。此时，
如果用一条线将破坏区连接起来，就形成了滑坡滑动时的
(b)145m7]~位
图3 不同水位高程下滑坡体的塑性区分布图
(d)175m7]~位
一3一
维普资讯 http://www.cqvip.com
高圣益等 基-T-Ft,hC~。的滑坡稳定性数值模拟分析
滑移 。
5．2 剪应变增量
判断滑坡体的(潜在)滑动面(带)，可根据其剪应变增量
的大小来判断：剪应变增量较大(绝对值)的部位，则为其
(潜在)滑动面(带)，变形破坏也多沿此处发生；剪应变增量
较小或基本上没有发生变化的部位，一般不会有潜在滑动
面的产生，因此，这些部位也不会发生较大的变形或破
坏。不同蓄水位条件下的剪应变增量集中带分布见图4。
由图4可知，水库水位上升到135m)=~，坡体剪应变集中
带没有贯通，滑坡体处于稳定状态；当水位上升到145m
时，剪应变增量的量值、范围都有所增大，但仍然没有贯
通整个滑坡体，滑坡体处于基本稳定状态；当水位升至
156m时，剪应变增量集中带基本贯通滑坡体，此时，滑坡
体处于极限平衡状态；水位上升到175m后，剪应变增量带
已经贯穿整个滑坡体，此时，如果滑坡体受到外力(地震
力)影响或暴雨条件下，将可能失稳破坏。剪应变增量集
中带的发展趋势反映了滑坡体在不同蓄水位条件下滑坡体
破坏的渐进破坏过程，剪应变增量集中带的分布位置显示
了滑坡可能存在的潜在滑动面。通过以上分析，滑坡的潜
在滑面基本受基岩面控制，这与通过勘察获得的滑坡潜在
滑动面位置基本—致。
一4—
6 结论
(1)通过对安乐寺滑坡稳定性的耦合分析表明，该滑
坡在135m水位下处于稳定状态。随着水库蓄水位的抬升，
滑坡的塑性区逐渐增大，剪应变增量集中带范围也进一步扩
大，但仍处于基本稳定状态。当蓄水至156m．~，滑坡的塑
性区剪应变集中带基本贯通整个滑坡体，滑坡处于极限平
衡状态。水库继续蓄水至175m，滑坡的塑性区已经贯通整
个滑体，此时，滑坡在受到外力(如地震、暴雨)影响而
可能失稳破坏。因此，该滑坡应该采取工程措施加固。
(2)水库蓄水是诱发库岸滑坡的重要因素之一，库水
作用下滑坡的变形及稳定性分析属于固一液两场耦合问
题。FLAC~ 作为—个大型的有限差分计算软件，为滑坡稳
定耦合分析提供了强有力的工具。
参考文献：
【1l】丁秀丽，付敬，张奇华。三峡水库水位涨落条件下奉节
南桥头滑坡稳定性分析【J】．岩石力学与工程学报，2004，
23(17)．
【2】FLAC manual。Itasca：Consulting Group,Inc．1997．5。
【3】许强，严明，黄润秋．某水电站左岸深裂缝对工程荷载
下边坡稳定性影响的 分析叨．地质灾害与环境保
护，2002， (3)．
维普资讯 http://www.cqvip.com