圆弧滑动有限元土坡稳定分析

文章编号：1000－7598－(2005) 10―1525―05

圆弧滑动有限元土坡稳定分析

殷宗泽1，吕擎峰1, 2

（1.河海大学 岩土工程研究所，南京 210098；2.兰州大学 土木工程与力学学院，兰州 730000）

摘 要：提出一种基于有限元应力变形计算的边坡稳定分析方法，仍假定滑动面形状为圆弧形，有限元网格由一组同心圆作为纬线，一组竖向线为经线构成。对两相邻圆弧线所夹的弧形带分析滑动力和抗滑力，建立平衡方程，确定安全系数，其中小值安全系数对应的弧形带为可能的滑动带。变化圆心位置用优选方法寻找最小安全系数对应的圆心，从而得出真正的滑动面。算例分析表明，该方法计算所得安全系数与Bishop法接近，是合理的，其突出优点是由有限元计算直接得出滑面上的应力，而不须作近似的插值处理，因而应力更准确。该方法可以考虑土的非线性变形特性，也更符合土的实际情况。此外，用有限元计算得出位移，亦可将稳定分析和变形联系起来，为现场通过位移监测来估计边坡的稳定性提供了可能性，同时也为膨胀土边坡的稳定分析中考虑膨胀性的影响提供了可能性。 关 键 词：边坡稳定；有限元法；圆弧法；非线性 中图分类号：TU 457 文献标识码：A

Finite element analysis of soil slope based on

circular slip surface assumption

YIN Zong-ze, LV Qing-fen

(1.Geotechnical Research Institute, Hohai University, Nanjing 210098, China;

2.School of Resources and Environmental Sciences, Lanzhou University, Lanzhou 730000, China)

Abstrct: A slope stability analysis method based on finite element analysis is developed. The slip surface is still assumed to be circular arc. The mesh of finite element consists of a family of circular arc lines and a family of vertical lines. The slip force and the resistance are analyzed for the arc belt which is between two neighbouring arc lines, and the equilibrium equation is set up. Then the factor of safety is determined. The arc belt corresponding to the smallest factor of safety is the possible slip belt. The position of circle center is changed and the circle center of real slip surface is seeked by the optimum seek method. For the real slip surface, the factor of safety is minimum. The computation examples show that the calculated factor of safety is quite close to that of Bishop’s method. And this shows that the proposed method is reasonable. The outstanding merit of the proposed method is that the stress state on slip surface can be obtained directly from finite element computation, and is more accurate. The method can reflect nonlinear deformation behavior of soil and is more close to the real situation. Besides, as one of the methods of slope stability analysis by finite element, the displacements are calculated. This offers the possibility to estimate the slope stability during construction by displacement monitor in situ. And it can also be used to analyze the slope of expansive soil in considering the influence of expansive deformation.

Key words: slope stability; finite element method; circular arc method; nonlinear deformation

1 引 言

边坡稳定分析方法一直是岩土工程中的重要研究课题。采用条分法进行土坡稳定分析已几十年，到摩根斯坦法的出现已发展到相当完善的程度。因此继续在条分法的框架下作出新的发展已经很难了。20世纪60年代以来，有限元法开始应用于土

坡的稳定性分析，为土坡稳定分析提供了新的思路。

有限元法用于土坡的稳定分析有其独特优点：（1）滑动面上的计算应力比较真实，可考虑应力-应变的非线性和弹塑性。（2）将稳定和位移的发展联系起来。现场观测只能测出位移的发展，无法测出安全系数，将稳定性和位移的发展联系起来就可以依据位移估计安全性的变化，这对施工中检测和

收稿日期：2004-04-28 修改稿收到日期：2004-09-13

作者简介：殷宗泽，男，1937年生，教授，博士生导师，主要从事土体本构关系及土工数值分析等研究。E-mail：yinzz@jlonline.com

1526 岩 土 力 学 2005年

控制边坡的稳定性是十分重要的。（3）用有限元可计算膨胀土边坡和黄土边坡在浸水条件下的变形，进而反映对稳定性的影响。

目前基于有限元方法分析边坡稳定主要有3种类型：

（1）对边坡作非线性有限元分析，得出边坡完整的应力变形成果，预测边坡区域由弹性变为塑性的完整演进过程，为边坡的治理、施工方法提供依据[1]。但它对稳定性的分析是定性的，不能给出安全系数的大小。

（2）用有限元计算的单元应力内插确定滑动

插计算带来的误差。由于要假定许多滑动面，每一个滑面对应一种有限元网格，网格必须能由程序自动生成，网格的剖分就要尽可能规律化、简单化。 2.1 计算方法的思路和步骤

（1）仍然假定滑动面为圆弧，不过不是一条弧线，而是有一定厚度的弧形带。

（2）在绘制有限元网格时，用一族同心圆弧线作为横向线，另外再用一族竖向线与其相交织成 网，编制程序自动生成网格。

（3）作第一次有限元计算得到各单元的应力分量。

（4）将相邻两同心圆弧所夹的弧形带作为可能的滑动带。对滑动带i，计算带内各单元形心处与弧面方向平行的面上的法向应力σnj和剪应力τj。过滑动带各单元形心作一圆弧，取该圆弧线以上的土体为脱离体，并令发挥的强度为τf/Fs，这里Fs为边坡整体滑动安全系数，建立力矩平衡方程：

nnσtanϕ+c

nj

R∑τjlj=R∑lj （1）

Fsj=1j=1

面上的应力，再根据力的平衡确定安全系数。条分

法要作出土条是刚性的假定，还要假定条间力的方向和分布。有限元法确定应力避免了这些假定，可使滑面上的应力更符合实际。加拿大的Geo-slope软件、Donald和Giam[2]的模式搜索法、邵龙潭[3]、Zou和Takuo[4] 等的方法都属于这一类，但其中有些方法比较复杂，难以应用。

（3）强度折减法。该方法将土的抗剪强度除以安全系数，用于有限元计算，如果算得的土坡恰恰达到破坏，则所选的安全系数就是实际安全系数。为此须要将安全系数从1.0逐步增大进行试算，直

～

到土体破坏[57]。该方法的关键是确定临界破坏状态的判断标准，目前判断标准主要有两种：一是最大节点水平位移达到某值；二是作迭代计算时达到不收敛。前者是用最大节点水平位移与折减系数的关系曲线的激烈转折点作为临界破坏状态的判断标准，但有时曲线没有明显的转折点。而对于后者，影响计算不收敛的因素很多，其安全系数难以较好确定。

进一步研究利用有限元法计算边坡稳定安全系数是十分必要的。本文提出了一种以有限元应力分析为基础的边坡稳定分析方法，仍假定滑动面是圆弧，有限元网格与滑面相应，可直接获得滑面上的应力，通过自适应网格剖分和优选方法搜索最危险的滑动圆弧，从而确定安全系数。

∑(σ则 Fs=

j=1

n

nj

tanϕ+c)lj

n

jj

∑τl

j=1

式中 j为滑动带上第j个单元；n为滑动带内单元数；lj为滑动带内第j个单元的弧长(单元内上下两圆弧长的平均值)。对各滑动带均作上述计算，求最小Fs，相应的滑动带就是局部危险滑动带。

（5）将找到的滑动带用等间距的同心圆弧细 分成10条弧形带，相应单元加密，竖向也适当加 密，其他滑动带不加密，重新形成网格体系后作第二次有限元计算，可得到较窄的危险滑动带，相应的Fs也更符合实际。

（6）以上仅仅是对一个给定的圆心确定了危 险滑动带和相应的安全系数。还要在平面内变动圆心的位置寻找最危险的滑动带。对于一个选择的圆心重复上述计算，计算危险滑动带所对应的Fs，用

2 圆弧滑动有限元分析方法

有限元法分析边坡稳定性，主要是利用有限元计算得出的应力场，推求滑动面上的法向应力和切向应力来计算安全系数。不少学者只通过一次有限元计算确定应力场，假定许多可能的滑动面，用某种内插方法求滑动面上的应力分量。本文所提出的方法是使网格的划分迁就滑动面，用有限元计算直接获得滑动面上的应力，省去内插计算，也避免内

0.618优选法在平面内优选，最后确定全局最小安 全系数Fsmin。 2.2 网格的自动剖分

网格生成必须由计算机程序自动实现，这是本文所提方法的关键之一。自动剖分程序可按以下步骤进行：

（1）对给定的圆心，选定一个有限元计算范 围，用最大半径Rm表示。可取Rm=a+2h，h为坡高；a为圆心到坡顶的垂直距离。对于有较深软土

第10期 殷宗泽等：圆弧滑动有限元土坡稳定分析

1527

的情况，可取Rm=a+3h。以Rm为半径画圆弧与坡面交于A，B两点（如图1所示），有限元计算区域就限制在弧AB以上。

（2）以圆心O到坡面上转折点C，D，…的距离为半径画圆弧，得若干圆弧形土带，土带不宜太宽，可规定粗网格的最小控制尺寸d=h/10，当相邻两圆弧的的间距，即弧形土带的的宽度∆R>d时，在其间等间距内插圆弧（如图1中细线所示）。内插圆弧数n由下式计算：

应力τ，根据单元所在位置的滑弧的倾角，由单元

x，y方向的应力分量求得，即

σn=σxsin2θ+σycos2θ−2τxysinθcosθ ;⎫⎪τ=τxy(sin2θ−cos2θ)−(σy−σx)sinθcosθ ⎪⎭

为单元处滑面与水平面的夹角。 2.4 圆心的优选

⎬ （3）

式中 σx，σy，τxy为单元形心处的应力状态；θ估计危险滑动带圆心的可能变化范围，即xo的上下限xo1，xo2和yo的上下限yo1，yo2。在y向用

n=int(∆R/d) （2）

则加密后圆弧间距为∆R1=∆R/(n+1)，如图1中圆弧CC1和DD1的半径差大于d，将该圆弧带加密为

0.618法在yo1和yo2选择yo。对每个yo在x向用0.618法在xo1和xo2之间优选xo。对每组圆心作有限元计算确定滑动带和Fs，通过上述优选确定最小安全系数及相应的滑动带。

与条分法相比，有限元法要复杂的多，但在优选寻找危险滑弧方面简化了。条分法确定危险滑弧一般要确定3个参数，圆心坐标(x, y)和滑弧半径。本文方法危险滑弧半径是在有限元计算时确定的，故优选只要确定圆心坐标(x, y)两个因素。

3条。加密粗圆弧带从而得到一系列圆弧。

（3）从所有圆弧与坡面的交点处作竖向线与 各弧线相交。

（4）若相邻两垂线的间距大于d，则用上述同样的方法加密。

这样就构成了初步计算网格。在第二次计算要加密滑动带时，可取∆R1=∆R/10，而竖向线间距

∆X1=∆X/5，且仅在初步确定的滑动带范围内加密（图2为算例1的最后网格）。 O 3 算例及分析

3.1 均质土坡

一均质土坡，坡高10 m，坡度1:1，土层的非线性弹性模型参数见表1。

a Rm h A h C D1 D C1 B

表1 土层的双曲线模型参数

Table 1 Parameters of hyperbolic model

γ / kN·m-3ϕ / (º)

c/ kPa

Rf

K n G

F D

1

22 28 30 0.8

150 0.40 0.35 0.01

图1 有限元网格的剖分

Fig.1 Mesh of finite element method

用本文方法计算的最小安全系数为1.86，危险滑动带的位置及相应的单元剖分见图2。此例

2.3 有限元计算

（1）边界线取最大半径圆弧，其上各节点的x和y方向位移为零，边界线离滑面要远些。

（2）本构模型可取邓肯模型或其它模型，包括弹塑性模型，这些模型一般都是增量形式。它本可以反映荷载的逐步施加，还可以模拟填筑和开挖过程，但这会使网格的自动生成增加困难。本方法中的网格一次生成，荷载(各单元自重)一次施加，用中点增量法进行一次重复计算作为结果应力状态。

（3）土坡和地基中可能存在若干土层，在网格划分时不考虑地层分界线，在有限元计算中单元参数取其形心所在的土层参数。

（4）滑面方向的应力分量，法向应力σn和剪

Bishop法的最小安全系数为1.85，滑弧的位置如图2中虚线所示。本文方法和Bishop法的计算结果接近，计算所得的最大水平位移为23.1 cm，位于坡面的中部，最大沉降为56.6 cm，位于坡顶滑面附近。

Bishop法滑弧 10 m y1 1 x

图2 有限元网格及危险滑动带的位置

Fig.2 Mesh of finite element method and the position of

critical slip surface

1528 岩 土 力 学 2005年

图3为有限元法和毕肖普法得到的滑动带上的剪应力和抗剪强度分布。可以看出，两种方法的计算结果基本规律一致，仅在坡脚部位有一些差异。毕肖普法假定土体为刚性体，还对土条间的作用力及分布作了假定，而有限元法用非线性应力-应变关系，差别是必然的，应该说有限元计算的应力更接近实际。

为18.1 cm，位于坡脚部位滑面的附近。最大沉降为25.9 cm，位于坡顶。有限元法和毕肖普法得到的滑动带上的剪应力和抗剪强度分布如图6所示。由图可见，尽管土层分布不均，两种方法所得结果是接近的，这也证明了本文方法的合理性。 Bishop法滑弧

图5 潜在危险滑动带的位置及相应的单元剖分 Fig.5 Meshing corresponding to the minimum factor and the critical slip

图3 滑动带上剪应力和抗剪强度分布

Fig.3 Distributing of shear stress and shear strength along

the slip surface

3.2 复杂土层边坡

一土坡的断面和土层分布如图4所示，土层的邓肯模型参数见表2。

2 10 m 1 土层② 土层① 图6 滑动带上剪应力和抗剪强度分布 Fig.6 Distributing of shear stress and shear strength along the slip 土层③ 4 结 语 本文方法得出的安全系数与条分法中的Bishop法大体相当，这表明本方法有足够的精度，与本文方法中滑动面上的应力由有限元直接算得，避免了内插产生的误差是有关的。此外，作为一种用有限元法分析边坡稳定性的方法，可以同时解出位 移，从而可与填筑或开挖过程中的位移观测值相联系评Rf图4 土坡的断面和土层分布 Fig.4 Section of slope and soil layers 表2 土层的非线性弹性参数 Table 2 Nonlinear hyperbolic model parameters 土层 ① ② ③ γ 18 φ c / kPa / kN·m-3 /(º) K n G F D价边坡的安全性。虽然该法要重复进行有限元计算，工作量较大，但现在微机计算速度很快，重复作二维的有限元计算实际上是很快完成的，用这种方法使滑面上应力更符合实际。 本文仅仅研究了圆弧滑动面的情况，对于非圆弧滑动面还要进一步研究。此外，对于填筑和开挖的模拟、膨胀土和黄土浸水变形的模拟虽然在道理上是可行的，但也还有待深入研究。 38.1 17 450 0.33 0.39 0.393.530.82730.72 0.43 0.142.920.6917.1 30.7 18 131 0.45 0.33 0.131.700.6616.4 27.7 20 本文方法计算的安全系数为2.40，Bishop法计算的安全系数为2.45。危险滑动带位置和Bishop法的危险滑弧如图5所示。计算所得的最大水平位移第10期 殷宗泽等：圆弧滑动有限元土坡稳定分析 1529

参 考 文 献

[1] Duncan J M. State of the art: limit equilibrium and finite

element anlysis of slopes[J]. Journal of Geotechnical Engineering, American Society of Civil Engineers, 1996, 122(7): 477－596.

[2] Giam S K, Donald I B. Determination of Critical Slip

Surfaces for Slopes Via Stress-strain Calculations[A]. 5th Australia New Zealand Conference Geomechanics[C]. Sydney: [s.n.]. 1988. 461－464.

[3] 邵龙潭, 唐洪祥, 韩国城. 有限元边坡稳定分析方法及

其应用[J]. 计算力学学报, 2002, 18(1): 81－87.

SHAO Long-tan, TANG Hong-xiang, HAN Guo-cheng. finite element method stress analysis slope stability pattern search method[J]. Chinese Journal of

[4]

[5] [6]

[7]

Computational Mechanics, 2002, 18(1): 81－87

Zou J Z, Williams D J, Xiong W L. Search for critical slip surfaces based on finite element method[J]. Canadian Geotechnical Journal, 1995, 32: 233－246.

Griffiths D V, Lane P A. Slope stability analysis by finite elements[J]. Geotechnique, 1999, 49(3): 387－403.

Dawson E M, Roth W H. Slope stability analysis by strength reduction[J]. Geotechnique, 1999, 49(6): 835－840.

张鲁渝, 时卫民, 郑颖人. 平面应变条件下土坡稳定有限元分析[J]. 岩土工程学报, 2002, 24(4): 487－490. ZHANG Lu-yu, SHI Wei-ming, ZHEN Yin-ren. The slope stability analysis by FEM under the plane strain condition[J]. Chinese Journal of Geotechnical Enginee- ring, 2002, 24(4): 487－490.

2005年《岩土力学》第6期被EI收录论文（38篇，收录率100 %）

No. 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38

论 文 题 名

瓦斯煤层裸露面蠕变失稳的时间预测研究 隧道掘进拱部运动力学模型研究

隧道及地下工程围岩稳定性及可靠性分析的极限位移判别 非均匀外载作用下油井套管强度特性的实验研究 一次动压煤矿巷道预应力锚索支护设计与参数优化 深部采空区尾砂胶结充填体强度特性试验研究 基坑土体侧向卸荷真三轴试验研究

准各向同性裂隙岩体中有效动弹性参数与弹性波速关系的研究 深埋巷道围岩稳定性分析与控制技术研究 浅埋采场初次来压顶板砂土层载荷传递研究

不排水循环荷载作用过程中累积孔压对细砂弹性剪切模量 Gmax的影响

隧道与地表构筑物相互影响的研究

坚硬顶煤弱化爆破的宏观损伤破坏程度研究 盾构法隧道施工变形智能模糊控制方法研究 万州侏罗纪红层软弱夹层特征 公路隧道全断面揭煤防突技术

公路隧道穿越采空区治理技术研究

粘土的粘塑性特性试验与三要素本构模型 通渝隧道涌突泥成因分析

泥岩蠕变行为的实验研究及其描述

石油开采引起的油藏压实与地面沉陷预测 残积土特殊应力路径的三轴试验研究

深埋隧道软硬交替复合顶板岩体变形破坏分析 单一岩石变形特性及本构关系的研究

岩石声发射Kaiser效应应用于地应力测试的研究 超前预测预报在笔架山隧道施工中的应用 堆载预压过程中的地基稳定性计算方法 隧道地应力测试及岩爆预测研究

隧道开挖施工的爆破振动监测与控制技术

某水利枢纽厂房开挖三维可视化与数值模拟耦合研究 米箭沟尾矿坝加高方案渗流场数值分析 高压水射流破岩规律的数值模拟研究 深部软岩巷道锚注支护数值模拟研究 岩石损伤断裂的沙堆元胞自动机模拟

盾构法隧道衬砌施工阶段受力特性的三维有限元模拟 块裂结构岩质边坡水力学模型及数值模拟

考虑饱和-非饱和渗流、温度和应力耦合的三维有限元程序研制 岩石力学参数数据库系统的开发和研究

作 者 页 码

鲜学福，李晓红姜德义，等 841－844 尹光志，岳 顺，代高飞 845－849 李晓红，王宏图，贾剑青，等 850－854 邓金根，田 红，王治中，等 855－858 高 谦，刘福军，赵 静 859－864 李一帆，张建明，邓 飞，等 865－868 何世秀，朱志政，杨雪强 869－872 王宏图，李晓红，杨春和，等 873－876 靖洪文，李元海，许国安 877－880 黄庆享 881－883 姬美秀，陈云敏 884－888 许 江，顾义磊，康骥鸣 889－892 索永录 893－895 朱忠隆，孙 钧，张庆贺 896－900 简文星，殷坤龙，马昌前，等 901－905 姜德义，刘 春，张广洋，等 906－909 李晓红，姜德义，刘 春，等 910－914 李建中，彭芳乐，龙冈文夫 915－919 顾义磊，李晓红，赵 瑜，等 920－923 万 玲，彭向和，杨春和，等 924－928 田 红，邓金根，周建良，等 929－931 方祥位，陈正汉，申春妮，等 932－936 贾剑青，王宏图，李晓红，等 937－940 姜永东，鲜学福，粟 健 941－945 姜永东，鲜学福，许 江 946－950 王青海，李晓红，夏彬伟，等 951－954 徐书平，刘祖德，鲍 华，等 955－958 康 勇，李晓红，王青海，等 959－963 陈 庆，王宏图，胡国忠，等 964－967 夏艳华，白世伟，倪才胜 968－972 柴军瑞，李守义，李康宏，等 973－977 孙清德，汪志明，于军泉，等 978－982 王连国，缪协兴，董健涛，等 983－985 尹光志，张东明，王 浩 986－989 张海波，殷宗泽，朱俊高，等 990－994 赵 瑜，李晓红，卢义玉，等 995－999

韦立德，杨春和 1 000－1 004徐高巍，白世伟，贺怀建 1 005－1 008