数值模拟与极限平衡法相结合的边坡稳定性分析

钟福生

(1.厦门紫金工程设计有限公司；2.紫金矿业集团股份有限公司)



数值模拟与极限平衡法相结合的边坡稳定性分析

钟福生1,2

(1.厦门紫金工程设计有限公司；2.紫金矿业集团股份有限公司)

摘要考虑到数值模拟法与极限平衡分析法在岩质边坡稳定性分析中各自的局限性，提出将数值模拟法与极限平衡分析法相结合对岩质边坡进行稳定性分析。以某矿山边坡为例，首先以离散单元法对边坡进行数值模拟，将边坡及其内部的滑动面划分成若干节理单元体，计算出节理单元上的法向和切向应力，根据极限平衡法计算出边坡稳定性安全系数，既考虑了滑动面内应力对岩质边坡稳定性的影响，又通过计算边坡的安全系数来评价边坡的稳定性，具有较好的推广前景。

关键词数值模拟极限平衡分析岩质边坡离散单元法安全系数

目前，边坡稳定性分析方法以极限平衡法与数值模拟法为主[1]。前者基于极限平衡理论,将滑动体视为刚体，可以算出边坡的稳定性系数并搜索出可能的破坏面，但没考虑岩土体内部的变形效应对边坡稳定的影响，同时还简化和假设条块间力的作用大小和方向，无法分析边坡的破坏机理[2-3]；数值模拟法虽然考虑了岩土体内部的应力应变效应，可分析边坡破坏的发展过程，却难以给出一个意义明确的安全系数，限制了在工程中的应用[4-6]。鉴于这两种方法在边坡稳定性分析中的局限性，将边坡稳定性分析的数值模拟法与极限平衡法相结合，对边坡进行数值模拟，分析出边坡的破坏机制和潜在滑动面，利用数值模拟的应力计算结果得出滑动面上各节点应力，通过应力变换，将其转化为与节理单元垂直和相切的应力，最后根据极限平衡法计算出边坡稳定性安全系数。

本文以数值模拟法中的离散单位法为例，应用上述方法对某岩质边坡进行了稳定性分析，并与单独采用数值模拟和极限平衡法的计算结果进行了比较。

1基本原理

1971年，Cundall提出了著名的离散单元法模型，经过几十年的不断完善，被广泛应用于不连续岩土介质的稳定性分析中[7-8]。同有限元法和有限差分法类似，离散单元法也需要将模型划分成单元。不同之处在于，在划分单元之前需要将模型根据实际存在的不连续面和假定的不连续面划分为块体，再划分成单元，与不连续面两侧相邻的单元节点既可以接触也可以分离，单元之间相互作用的力可以根据力和位移的关系求出[9]。

1.1物理方程

如图1所示，离散单元法将块体之间的作用分成法向力Fn和切向力Fs，用假想的法向弹簧和切向弹簧来描述块体之间的位移和力的关系[10]。其中法向弹簧只有在受压状态下起作用，切向弹簧只有在受拉状态下起作用。

图1　离散单元块体之间的法向力与切向力

Fn和Fs的计算公式分别为：

(1)

(2)

式中，kn和ks分别为法向刚度和切向刚度；un和us分别为法向位移和切向位移，un必须为压缩位移，即un<0，当un>0时，即块体间受到拉张超过接触面之间的抗拉强度而分离，上诉两式均不成立。

块体除了受拉张分离而导致岩体破坏的情况之外，接触面间的切向力过大，超过其剪切强度也会导致岩体破坏，即Fs≥tanφFn+c时，式(2)失效。

1.2动力学方程

单质体的运动决定于其所受的不平衡外力和力矩，其平动和转动需要用动力学方程来求解。假设一个单质体在某时刻所受外力是F(t)，由牛顿运动定律可得：

(3)

根据中心差分公式，加速度可表示为：

(4)

将式(3)代入，式(4)可改写为：

(5)

假设t=0时刻单元所受外力已知，则可计算出同时刻单元的加速度，可以根据式(5)计算Δt/2时刻的单元速度，根据式(6)可计算出Δt时刻的单元位移，这样又回到式(3)中。

(6)

1.3安全系数的定义

目前边坡安全系数的定义主要有3种：①抗剪力与剪应力之比；②岩土体强度的安全储备系数；③特定滑动面上抗滑力与致滑力之比。对于给定的滑动面，适合采用第3种定义作为安全系数。根据定义可得安全系数Fs：

(7)

式中，φ为潜在滑动面的内摩擦角；m为滑动面上测点总数；σni为t时刻滑动面上第i个测点的法向应力；τsi为t时刻滑动面上第i个测点的剪切应力；Δl为测点间距；l为滑动面总长度。

2工程实例

2.1工程概况

某大型露天铁矿是一座尚未开发的矿山，设计时下盘边坡布置在顺层岩体内，有发生滑动破坏的可能，需要对该区域的边坡进行稳定性分析，以确定是否有重新布置的必要，或采取加固措施。图2显示边坡内的岩层分布情况。从图2中可以看出，最靠近坡面的岩层面出露于临空面，形成的切割块体很可能沿岩层面发生平面滑动。

图2　边坡内岩层分布情况

2.2建立计算模型

表1为中区下盘边坡岩体力学参数，模型的X轴和Y轴范围分别为1 150和520 m。参照表1提供的岩体参数对模型进行初始地应力计算，获取地应力分布情况后对模型进行开挖解算，形成图3所示的边坡模型。计算采用的岩体节理参数见表2。

表1　中区下盘边坡岩体力学参数

注：①～⑤为图2中编号。

2.3数值模拟结果

开挖解算后获取的竖直地应力和剪切应力场如图4所示。开挖后的边坡在竖直方向的应力分布较为平滑，靠近节理处的应力梯度近似平行于节理。

图3　中区下盘边坡数值模型网格剖分图(开挖后)

编号抗拉强度/kPa抗剪强度ϕ/(°)c'/MPa法向刚度/GPa切向刚度/GPa①03015012060②0.013222012060③0.013427012060

注：①～③为图3中编号。

图4　开挖后的边坡竖直应力和剪切应力场

将图3层面1划分成34个节理单元，共35个节点；提取该层面上各节理单元的法向力Fni和切向力Fsi后，绘制图5所示的应力分布曲线。从图5中可以看出，坡体上部的Fni和Fsi变化较大，造成这种现象的原因可能是滑动体上部岩土体承受了较大的拉张力，岩体受拉张而屈服。靠近坡脚处节理单元上的Fni较小。

图5　节理单元的法向力和切向力沿层面1的分布

2.4边坡稳定性安全系数

计算节理单元局部安全系数：

(8)

根据计算结果绘制局部安全系数柱状图，如图6所示。随着接触单元接逐渐近坡体下部，其安全系数逐渐减小，直至小于1。说明坡体的下缘承受了较大的挤压力，有受挤压下滑的趋势。根据式(7)得出边坡安全系数为1.31。

2.5与极限平衡法比较

若层面1作为滑动面，按照传统的极限平衡法计算，潜在滑体安全系数为1.38。与式(7)计算结果相差不大。

图6　层面1上的安全系数分布

3结论

将数值模拟与极限平衡计算相结合，以岩质边坡为例进行了边坡稳定性分析，计算结果表明：①在滑动面上，可用数值模拟法计算出节理单元上的法向力和切向力，计算节理单元上的局部安全系数，有助于判断边坡的破坏模式;②靠近坡体滑动面上缘的岩土体易受拉张屈服破坏，从侧面解释了岩质边坡破坏时通常伴随着张裂缝的现象；靠近坡体滑动面下缘的岩土体易受压剪屈服破坏;③本文所介绍的方法用于岩质边坡稳定性分析是可行的，可解决滑动面复杂的岩质边坡安全系数计算问题。由于本文所提供的算例较为简单，该法尚需在工程实践中进一步验证。

参考文献

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

[2]曾亚武,田伟明.边坡稳定性分析的有限元法与极限平衡法的结合[J].岩石力学与工程学报,2005,24(S2):5355-5359.

[3]刘明维,郑颖人.基于有限元强度折减法确定滑坡多滑动面方法[J].岩石力学与工程学报,2006,25(8):1544-1549.

[4]赵尚毅,郑颖人,邓卫东.用有限元强度折减法进行节理岩质边坡稳定性分析[J].岩石力学与工程学报,2003,22(2):254-260.

[5]陈建宏,钟福生,杨珊.基于有限元折减强度法与极限平衡法结合的岩质边坡稳定性分析[J].科技导报,2012,30(11):38-42.

[6]郑宏,田斌,刘德富,等.关于有限元边坡稳定性分析中安全系数的定义问题[J].岩石力学与工程学报,2005,24(13):2225-2230.

[7]Itasca Consulting Group. Inc. UDEC (Universal Distinct Element Code) User's Manual, Version3.0[R]. MinneaPolis: Itasca Consulting Group,1996.

[8]赵红亮,朱焕春,朱永生,等.UDEC及其在岩质高边坡稳定性分析中的应用[J].四川大学学报:工程科学版,2007,39(S):192-195.

[9]王泳嘉,邢纪波.离散单元法同拉格朗日元法及其在岩土力学中的应用[J].岩土力学,1995,16(2):1-14.

[10]Hart R D. An introduction to distinct element modeling for rock engineering. Comprehensive Rock Engineering[M]. Sao Paulo: Pergamon Press, 1993.

(收稿日期2015-10-13)

Slope Stability Analysis by Combining Numerical Simulation with Limit Equilibrium Method

Zhong Fusheng1,2

(1. Xiamen Zijin Engineering Design Co., Ltd.; 2. Zijin Mining Group Co., Ltd.)

AbstractConsidering the limitations of numerical simulation and limit equilibrium method in rock stability analysis,combing numerical simulation with limit equilibrium method to analyze the slope stability problems.Taking a open-pit slope as an example,firstly,the discrete element method is adopted to conduct numerical simulation,the sliding surface of the open-pit slope is divided into a number of joint elements to calculate the normal and shear stresses;then,based on the limit equilibrium method to calculate the safety coefficient of the open-pit slope.The research method proposed in this paper not only consider the influence of the stress in sliding surface to the rock slope stability,but also evaluate the stability of the open-pit slope by calculating the safety coefficient,therefore,it has good promotion prospects.

KeywordsNumerical simulation, Limit equilibrium method, Rock slope, Discrete element method, Safety coefficient

钟福生(1987—)，男，工程师,硕士，361006 福建省厦门市湖里区翔云三路128号。