三维数值模拟在某高速公路边坡稳定性分析中的应用

宋民崇，余云燕

(1.中国市政工程华北设计研究总院有限公司，天津 300074;2.兰州交通大学 土木工程学院,甘肃 兰州 730070)



三维数值模拟在某高速公路边坡稳定性分析中的应用

宋民崇1，余云燕2

(1.中国市政工程华北设计研究总院有限公司，天津 300074;2.兰州交通大学 土木工程学院,甘肃 兰州 730070)

摘要：近年来三维数值模拟在边坡稳定性分析中被广泛应用，文中以某滑坡为工程背景，运用FLAC3D软件分析滑坡的稳定状态，对滑坡体的位移、应力、应变、最大不平衡力和特定记录点的位移与时步关系进行分析，模拟出滑动面的位置，综合评价滑坡的稳定状态，从而为其它同类工程提供有意义的经验。

关键词：滑坡；稳定性分析；FlAC3D；数值模拟

随着国家经济的发展，对基础设施的投资日益加大，特别是西部大开发的开展，铁路、公路建设量逐年增长，而且不断向山区延伸，因地形、地貌等条件的约束，不可避免地遇到大量的边坡工程，有些地区甚至出现高陡边坡，不少公路、铁路线段常有滑坡等地质灾害发生，造成运输中断，大量人员伤亡和严重的经济损失。所以进行边坡稳定性研究非常重要。

1工程概况

1.1地形地貌

坡体位于西南山区，地貌形态属于侵蚀构造中山地貌，山体主要由板岩、千枚岩等组成。境内山峦重叠，山坡多在35°以上，悬崖峭壁随处可见，山脉为石质山地，北侧为下古生界碧口群石英岩夹变质砂岩组成的岩质岸坡，上部覆盖第四系薄层坡积碎石土，多为岩石风化残积物。

1.2地层岩性

滑坡区地层岩性主要为第四系坡积碎石、滑坡堆积块石、碎石和土、下元古界碧口群，岩性主要有板岩、砂质板岩、变质砂岩、石英片岩。从新到老依次如下：

第四系坡积碎石层：分布于滑坡北侧平缓山坡上。灰黄色，主要由碎石和土组成，碎石成分为变砂岩、板岩等，分级差，成棱角状，结构松散，局部具有架空现象。

滑坡堆积块石、碎石层：分布于南侧，主要由块、碎石和土组成，地表植被发育，富含腐殖质。前缘岩体破碎，多弯曲变形。

下元古界碧口群板岩：岩石的力学性能较好，具有抗风化能力强、岩体致密而坚硬的特点。变质砂岩夹石英片岩，新鲜面呈青灰色，微具片理，变余细粒结构，厚层状结构，单层厚0.6～0.8 m，岩性致密坚硬。

1.3水文条件

所处区域降水量较大，对坡体的影响主要为两种途径：①基岩裂隙水。基岩裂隙水主要受大气降水补给，沿层理、节理、裂隙下渗向河谷排泄。②孔隙潜水。赋存于第四系松散覆盖层中，受大气降水补给，排泄以地下径流和蒸发为主。

地下水水质及水化学类型：经地表水取样并进行水质分析，不含侵蚀性CO2、水化学类型基本为HCO3-Ca-Mg、SO4-Ga型水。滑坡区地下水对混凝土无腐蚀性，对钢筋混凝土结构中的钢筋无腐蚀性。

滑坡体目前主要受自然条件的影响，包括地形地貌、地层岩性、水文地质条件等，没有人类活动如施工扰动的影响。某滑坡工程地质断面图如图1所示。

图1　某滑坡工程地质断面图

2数值计算

本文采用FLAC3D对滑坡的稳定安全系数进行计算，基于强度折减法的数值计算相对于传统的极限平衡法所得到的安全系数有很多优点：可以模拟出复杂的边界条件，地质结构特征，各种作用荷载和岩土体的各种力学特性，并且可以考虑岩土体的应力与应变的关系。基于FLAC3D的强度折减法不需要预先假定滑面，在求解安全系数时就能够自动搜索出滑移面的位置和形状，可以对滑坡的渐近破坏过程进行一定的分析。

强度折减方法应用到有限差分或有限元分析中可以表述为:保持岩土体的重力加速度为常数，通过逐步减小抗剪强度指标，将c，φ值同时除以折减系数Fs，得到一组新的强度指标ci，φi，然后进行有限元分析，反复计算直至边坡达到临界破坏状态，此时采用的强度指标与岩土体原具有的强度指标之比即为该边坡的安全系数Fs。公式如下:

2.1计算模型的建立

根据地质资料建立三维计算模型，FLAC3D建模工程浩大、繁琐，本文采用ANSYS建立模型，然后将模型导入FLAC3D进行数值计算，边坡前缘的脚点为坐标原点(0,0,0) ，X方向沿着坡脚方向向右为正，Y方向以滑动向前为负，Z向上为正，向下为负，高度为150 m，边坡长150 m，宽150 m，初始地应力为自重地应力场。边界条件为：X,Y,Z方向为固定约束，滑坡体主要由碎石土组成，滑床主要由板岩构成。计算参数选取根据勘测资料，按照表1取值，计算工况为天然状态。

表1　岩土体力学参数

在计算岩土体的剪切模量和体积模量时，应根据以下公式计算：

式中：E为岩土体弹性模量，B为岩土体体积模量，S为岩土体剪切模量，μ为材料泊松比。

2.2数值计算结果

应用FLAC3D软件进行数值模拟，计算得到滑坡安全系数为0.97，如图2所示，并对位移、应力、应变、最大不平衡力和特定记录点的位移与时步关系进行分析。图3为ANSYS网格划分示意图，图4为滑坡体示意图。

图2　安全系数计算结果

图3　ANSYS网格划分示意图

图4　滑坡体示意图

由于滑坡体主要成分为碎石土，力学性质较低，滑床为板岩，较为坚固，滑坡内部应力场由于外界条件的影响产生变化，本文应用TECPLOT图形处理软件出图，下面对计算结果进行分析。

图5为滑坡体总位移云图，图6为滑坡体总体位移中线位置剖面云图，在总体位移云图中可以看到，滑坡体最大位移为15 cm，滑坡体位移从后缘中心位置向外逐渐减小，最大位移产生在滑坡体的后缘中心位置处，从中线位置剖面云图中可以看出最大位移是从滑坡体后缘中心位置向下逐渐减小。

图5滑坡体总位移云图

图6　滑坡体总位移中线位置剖面云图

图7为滑坡体X方向位移云图，图8为滑坡体X方向中线位置位移剖面云图。由图7、图8可以看出，滑坡体在X方向的最大位移为3 cm，在X方向上产生位移很小，产生的位移是由于滑坡体一侧略微高于另一侧，在总位移云图中可以看出蓝色区域位移为0，其它区域位移为毫米级。

图7　滑坡体X方向位移云图

图8　滑坡体X方向中线位置位移剖面云图

图9为滑坡体Y方向位移云图，图10为滑坡体Y方向中线位置位移剖面云图。由图9、图10可以看出，滑坡体在Y方向产生位移，即为滑坡滑动方向，最大位移为11 cm，在Y方向产生位移主要是滑坡向下滑动导致的，在中线剖面云图中也看出位移是由上到下逐渐减小。

图9　滑坡体Y方向位移云图

图10　滑坡体Y方向中线位置位移剖面云图

图11为滑坡体Z方向位移云图，图12为滑坡体Z方向中线位置位移剖面云图，从图11可以看出，滑坡体最大位移为11 cm，即在上下方向上产生竖向位移，从图12可以看出，Z方向位移也是由滑坡体后缘中心位置由上到下逐渐减小的。

图11　滑坡体Z方向位移云图

图12　滑坡体Z方向中线位置位移剖面云图

图13为主剖面最大主应力云图，图14为主剖面最小主应力云图，由图13、图14可以看出，滑坡体最大主应力在分布上基本平行于斜坡临空面，最大主应力约为0.2～2.8 MPa，最小主应力约为0.05～0.55 MPa，都为压应力。

图13　主剖面最大主应力云图

图14　主剖面最小主应力云图

图15为滑坡体X方向应力云图，图16为滑坡体Y方向应力云图，图17为滑坡体Z方向应力云图，由图15～17可以看出，X方向应力较小，这与X方向上产生位移较小相符合，Y方向的应力最大值出现在坡脚，Z方向应力较大，即竖直方向上产生的应力值较大，并且随着深度的增加而增加，这符合地应力赋存规律。

图15　滑坡体X方向应力云图

图16　滑坡体Y方向应力云图

图17　滑坡体Z方向应力云图

图18、图19分别为滑坡体剪应变增量云图和剪应变增变率云图，由图18、图19可以看出，最大值出现在滑坡体的后缘，因此，滑坡体的破坏应该是后缘先开裂。在剪应变增量云图中同时可以看出滑坡体滑动面的位置。

图18　滑坡体剪应变增量云图

图19　滑坡体剪应变增变率云图

图20　滑坡体内点(75,95,120)X方向位移变化历史图

在FLAC3D计算时，可在滑坡体内设置一点(75,95,120)，记录这一点随时步的变化，同时也可以记录滑坡体最大不平衡力的变化历史。图20～22分别为滑坡体内某一点X方向、Y方向和Z方向的位移变化历史，图23为最大不平衡力变化历史。由图20～22可以看出，X方向产生微小的位移，到计算结束时，记录点有2.2 cm的位移，随着时步的增加，位移是从0逐渐增加的，记录点Y方向产生8.1 cm的位移，即指向临空面方向，也是滑坡滑动方向，记录点在Z方向上产生7 cm的位移，即向下有7 cm的位移，此记录点在滑坡体内中心位置。对记录点的选取有一定的代表性，记录点要能反映滑坡体的稳定状态，通过位移的变化来说明滑坡体所处状态。由图23可以看出，最大不平衡力逐渐趋于0，计算是收敛的。

图21　滑坡体内点(75,95,120)Y方向位移变化历史

图22　滑坡体内点(75,95,120)Z方向位移变化历史

图23　最大不平衡力变化历史

图24、图25分别为滑坡体位移矢量图和速度矢量图。由图24、图25可以看出，这两个矢量云图的矢量方向基本是重合的，显示出了滑坡体的运动趋势，即向着主滑方向滑动，主要是Z方向和Y方向，这两张云图同时也显示了滑坡体的范围。

图24　滑坡体位移矢量图

图25　滑坡体速度矢量图

3结论

本文主要是对某高速公路滑坡进行稳定性分析，应用FLAC3D软件对其进行了定量的计算，所做的主要工作和得到的相关结论可以概括为以下几点：

1) 根据地质资料运用ANSYS建立三维计算模型，然后将其导入FLAC3D软件进行计算，计算得到滑坡体安全系数为0.97，并对位移、应力、应变、最大不平衡力和特定记录点的位移与时步关系进行定量的分析。

2)通过三维数值模拟，对滑坡体运动趋势进行分析，在其主滑方向即Y方向最大位移是11 cm，竖直向下即Z方向的位移为11 cm，X方向上位移较小，最大位移仅为3 cm。

3)对滑坡体应力的分析得出Y方向的应力最大值出现在坡脚，Z方向的应力值较大，并且随着深度的增加而增加，符合地应力赋存规律，对剪应变增量的分析可以得出滑坡的破坏将是后缘先开裂，同时也可以判定出滑坡体滑动面的位置。

滑坡体位于高速公路沿线，目前处于不稳定状态，所处区域多雨且地震活动频繁，为了保证日后公路运营的安全，建议采取边坡防护、治理措施。

参考文献

[1]徐邦栋.滑坡分析与防治[M].北京:中国铁道出版社,2001.

[2]郑颖人,赵尚毅.有限元强度折减法在土坡与岩坡中的应用[J].岩石力学与工程学报,2004,23(19):3381-3388.

[3]刘波,韩彦辉.FLAC原理、实例与应用指南[M].北京:人民交通出版社,2005.

[4]王恭先,徐峻龄,刘光代,等.滑坡学与滑坡防治技术[M].北京:中国铁道出版社，2007.

[5]宋民崇，余云燕，赵德安.FLAC软件在边坡稳定性分析中的应用[J].黑龙江工程学院学报(自然科学版)，2009，23(4)：24-27.

[6]陈育民,徐鼎平.FLAC/FLAC3D基础与工程实例[M].北京:中国水利水电出版社,2009.

[7]许强,黄润秋,程谦恭,等.三峡库区泄滩滑坡滑带土特征研究[J].工程地质学报，2003，11(4)：354-359.

[8]张永安.滇中红层路堑边坡与路堤稳定性研究[D].昆明：昆明理工大学，2010.

[9]莫素华.边坡稳定性研究方法[J].中国水运(学术版)，2007，7(2):100-102.

[10] 俸锦福;郝勇;颜丙山.边坡稳定性分析理论及方法综述[J].中国水运(下半月),2009,9(7):183-184.

[责任编辑：刘文霞]

Application of three-dimensional numerical simulation tothe stability analysis of a highway slope

SONG Minchong1，YU Yunyan2

(1.North China Municipal Engineering Design and Research Institute Co.,Ltd，Tianjin 300074，China；2.College of Civil Engineering，Lanzhou Jiaotong University，Lanzhou 730070，China)

Abstract：In recent years, three-dimensional numerical simulation has been widely used in landslide stability analysis.Taking a landslids as the engineering background,this paper makes use of FLAC3Dsoftware to analyze the stability of this landslide, such as landslide displacement, stress, strain, the maximum unbalanced force and the relationship between the displacement of a specific record point and the time step, and simulates the position of the sliding surface, and the comprehensive evaluation of the stability of the landslide,in order to provide a meaningful experience for othersimilar projects.

Key words：landslide; stability analysis; FlAC3D; numerical simulation

中图分类号：TU457

文献标识码：A

文章编号：1671-4679(2016)01-0004-06

作者简介：宋民崇(1983-)，男，工程师，研究方向：道路与桥梁.

收稿日期：2015-09-23