基于能量突变的强度折减法边坡失稳判据

刘新荣，涂义亮，钟祖良，刘永权



基于能量突变的强度折减法边坡失稳判据

刘新荣1, 2，涂义亮1, 2，钟祖良1, 2，刘永权1, 2

(1. 重庆大学土木工程学院，重庆，400045；2. 重庆大学山地城镇建设与新技术教育部重点实验室，重庆，400045)

为了建立有限元强度折减法(SSR−FEM)中边坡失稳的能量突变判据，基于能量原理，推导强度折减过程中边坡能量反应的计算理论，二次开发折减过程中能量反应的FLAC3D计算程序。获得二维边坡能量反应随折减系数的变化规律，与常见3类失稳判据判定的安全系数对比。并利用三维边坡算例进一步验证能量突变判据的适用性。研究结果表明：当边坡达到临界失稳状态时，重力势能降、动能增量和耗散能增量均突然增大，而弹性应变能增量突然减小，表明能量突变可作为一类边坡失稳判据。能量突变判据具有物理意义明确、整体性强和易于判断等优点，尤其适用于三维边坡。

有限元强度折减法；边坡；失稳判据；能量突变

有限元强度折减法具有无需假定滑动面位置和形状、考虑了岩土体本构关系和能够得到边坡变形破坏过程等优点，因而在边坡工程中得到了较广泛的应 用[1−3]。但该方法仍存在一些限制其应用的缺点，如边坡失稳判据没有达成统一。依据不同的失稳判据判定的安全系数可能不同，甚至差异较大，因此，研究有限元强度折减法中边坡的失稳判据具有重要的理论意义和工程价值。目前，国内外许多专家学者对此展开了深入的研究，根据研究内容，他们主要从2方面展开。一方面，许多学者对当前最常见的数值计算收敛性判据[4−6]、塑性区贯通判据[7−9]和特征点位移突变判据[10−11]3类判据的影响因素、计算精度、适用条件及统一性等进行了研究。如，赵尚毅等[6]认为塑性区贯通是边坡破坏的必要不充分条件，可把有限元计算是否收敛作为边坡破坏的依据。栾茂田等[7]认为广义塑性应变分布及其发展状态可以作为边坡失稳判据。刘祚秋等[9]认为采用一定幅值的等效塑性应变从坡脚到坡顶上下贯通作为边坡破坏的标准是适宜的。裴利剑等[12]认为常见3类失稳判据本质上具有一致性和统一性，人为误判和数值计算误差是导致各判据判定结果差异的原因。另一方面，一些专家学者在常见3类失稳判据的基础上，提出了新的失稳判据。如刘晓宇 等[13]提出局部化带贯通判据；郑颖人等[14]通过算例分析，提出了剪切应变增量判据；吴春秋等[15]提出了动力学失稳判据，并证明了该判据的优点；柴红保等[16]认为弹性应变能随折减系数的变化规律可以作为边坡失稳判据。以往的失稳判据主要是从应力−应变和强度变形等角度进行研究，而边坡岩土体作为一种复杂的非线性材料，注定其变形失稳过程的应力−应变状态非常复杂，从而，有可能给这些失稳判据的研究带来困难[16]。事实上，边坡在变形失稳过程中，始终与外部环境交换着能量，且遵循能量守恒定律，边坡发生失稳的根本原因是能量的驱动作用。本文作者基于能量原理，推导有限元强度折减过程中边坡能量反应的计算理论，分析折减过程中各能量的反应变化，建立基于能量突变的失稳判据，并通过算例分析验证该失稳判据的适用性。

1 能量突变失稳判据的建立

1.1 边坡变形失稳过程的能量反应分析

利用能量原理研究边坡的变形失稳过程，首先需要确定1个研究体系。通常可以将承受力、热、光、电、生物和化学等作用的边坡作为研究对象，而边坡存在的空间称为外部环境[17]。

边坡体系与外部环境处于动态平衡的能量场中，外部环境通过力、位移和温度等对边坡的约束作用本质上是能量的输入、传递、转化和耗散过程。而边坡岩土体既是能量传递的介质，又是能量转化和耗散的载体。外荷载、自重应力、温度、电磁辐射等作用在边坡上，以机械能、内能和辐射能等能量形式输入到边坡内。输入的能量一部分使边坡产生弹性变形，并以弹性应变能的形式存储于岩土体内，另一部分被边坡岩土体的塑性变形、内部损伤、各种声热及电磁辐射等耗散。但边坡承受输入能量的能力有限，即边坡岩土体所能存储的弹性应变能和所能耗散的能量有限[17]。当输入的总能量大于所能承受能量时，多余的能量将转化为边坡动能，而此时边坡表现为滑动破坏。边坡能量的输入、传递、转化和耗散过程始终处于动态平衡状态，并且遵循能量守恒定律，即外界输入的总能量等于边坡体系存储的能量、耗散的能量和输出的能量总和，可表示为

式中：Δ为外荷载、自重应力等对边坡所做的机械功；Δt为环境温度、电磁辐射等输入到边坡内的热能；Δe为边坡弹性应变能增量；Δd为边坡耗散能增量；Δk为边坡动能增量。

1.2 边坡强度折减过程的能量计算

在有限元强度折减过程中，边坡体系的能量反应满足能量守恒定律。通常折减过程只考虑自重应力的作用，忽略外部荷载、环境温度、声辐射和电磁辐射等的影响，即Δt=0 J。因此，由式(1)可得，对于任意的折减系数，边坡体系的能量平衡方程可简化为

以边坡体系中任意单位体积的单元为研究对象，单元处于3向应力状态，若3个主应力分别为1，2和3，主应变分别为1，2和3，单元重心的速度和在铅垂方向的高度分别为和，材料弹性模量和泊松比分别为和，密度为，则单元所存储的重力势能g、弹性应变能e和动能k为：

边坡中的应力−应变、速度−位移及材料参数等并非均匀分布，而是坐标的函数。因此，由式(3)~(5)可求整个边坡体系存储的重力势能g、弹性应变能e和动能k，即

由式(6)~(8)可求边坡从初始状态到折减系数为过程中，体系的重力势能降Δg、弹性应变能增量Δe和动能增量Δk，如式(9)~(11)所示。

将式(9)~(11)代入式(2)，得边坡从初始状态到折减系数为过程中，边坡的耗散能增量为

1.3 边坡强度折减过程的能量反应分析

为了研究边坡强度折减过程的能量反应，采用文献[6−8]中的边坡算例进行分析，如图1所示。土体为服从摩尔−库仑屈服准则与非关联流动法则的理想弹塑性材料，土的重度=20 kN/m3，黏聚力=42 kPa，内摩擦角=17°，弹性模量=100 MPa，泊松比=0.3。采用FLAC3D软件分析，建立平面应变模型，厚度为1 m，共包含1 842个节点和850个单元。边界条件为边坡两侧水平约束，底部固定约束。特征点1~3分别位于边坡坡顶、坡中和坡脚。收敛条件为最大不平衡力与节点力平均值之比小于10−5。

单位：m

利用强度折减法折减材料强度参数，采用FLAC3D内嵌的FISH语言，根据式(9)~(12)编制折减系数下边坡重力势能降Δg、弹性应变能增量Δe、动能增量Δk和耗散能增量Δd的计算程序，并提取特征点1的竖直方向位移Δ1z和特征点2的水平方向位移Δ2x。表1所示为程序的部分计算结果。

表1 部分折减系数下边坡的特征点位移和各能量增量

利用数值计算收敛性判据、特征点位移突变判据和塑性区贯通判据判定边坡安全系数。

由表1可知：当折减系数为1.23时，计算收敛；当折减系数为1.24 时，计算将不收敛，因此，采用数值计算收敛性判据判定边坡的安全系数等于1.24。

作出折减系数−特征点位移Δ曲线，如图2所示，位移以向右和向上为正。由图2可知：当折减系数从1.22增至1.23时，特征点位移发生突变。因此，采用特征点位移突变判据判定该边坡的安全系数应等于1.23。

1—特征点1竖直方向位移；2—特征点2竖直方向位移；3—特征点2水平方向位移。

观察各折减系数下边坡的等效塑性应变云图，当折减系数为1.21时，塑性区尚未贯通；当折减系数为1.22时，塑性区完全贯通；当折减系数为1.23时，边坡沿着塑性区产生较大滑动，如图3所示。因此，采用塑性区贯通判据判定边坡的安全系数应等于1.22。

F：(a) 1.21；(b) 1.22；(c) 1.23

为了分析边坡能量反应随折减系数的变化规律，作折减系数−边坡能量增量曲线，如图4所示。

(a) 所有能量增量；(b) 弹性应变能增量；(c) 动能增量

由图4和表1可发现：1) 重力势能降随着折减系数增大而增大，当折减系数小于1.24时，其增长速度比较平稳，当折减系数大于1.24时，增长速度突然增大；2) 重力势能降绝大部分转化为耗散能，且耗散能随折减系数的变化规律近似于重力势能降的变化规律，当折减系数为1.24时，曲线出现拐点；3) 重力势能降仅有一小部分转化为弹性应变能，且弹性应变能增量随折减系数的变化规律表现为先增大后减小，当折减系数为1.24时，弹性应变能增量达到峰值；4) 当折减系数小于1.24时，边坡的动能极小，考虑FLAC3D计算的误差，可以认为此时边坡的动能为 0 J，当折减系数大于1.24时，动能突然增大到不可忽略的水平。综上所述，当折减系数为1.24时，4种能量增量都出现了明显的突变。

事实上，在强度折减过程中，边坡的能量反应与稳定性存在密切的关系。折减过程中输入边坡的能量全部来自重力势能降，其一部分因材料的塑性变形而耗散，另一部分以弹性应变能的形式存储于材料内部。重力势能降随着折减系数的增大而平稳增大，转化成的弹性应变能和耗散能也将平稳增大。但边坡存储和耗散能量的能力是有限的，当重力势能降增大到这一极限能量时，边坡达到临界失稳状态；当重力势能降进一步增大时，多余的能量将转化为动能，滑体由静止变成运动，边坡发生失稳破坏。重力势能降因滑体的突然滑动而突然增大；边坡滑动将耗散大量的能量，耗散能也将突然增大；而边坡所存储的弹性应变能将得到释放，弹性应变能突然减小。

因此，可以利用边坡的重力势能降、耗散能增量、弹性应变能增量和动能增量的突变作为强度折减法的失稳判据。采用能量突变判据可判定上述边坡的安全系数应等于1.24，与常见3类失稳判据判定的安全系数非常接近。

2 有限元算例验证

前文通过分析一个简单的算例并不足以说明能量判据的适用性，为了不失一般性，防止分析计算的偶然因素影响，下面将利用前面的方法进一步分析2个算例，以验证能量突变判据的适用性。

2.1 算例一

本算例采用图1所示的计算模型，材料参数和边界条件完全不变，几何尺寸除坡角改变外，其他均不变。分别计算坡角为 30°，35°，40°，45°和 50°时，边坡的安全系数。

通过计算得到不同坡角边坡的折减系数−能量增量关系曲线如图5所示。根据不同的失稳判据判定的安全系数如表2所示。由表2可知：基于能量突变判据判定的安全系数与常见的3类失稳判据判定的安全系数近似相等。

(a) 重力势能降；(b) 耗散能增量；(c) 弹性应变能增量；(d) 动能增量坡角/(°)：1—30；2—35；3—40；4—45；5—50。

表2 不同失稳判据下各坡度边坡的安全系数

2.2 算例二

三维边坡模型如图6所示，边坡材料为服从摩 尔−库仑屈服准则与非关联流动法则的理想弹塑性材料，重度=25 kN/m3，黏聚力=100 kPa，内摩擦角=45°，剪切模量=200 MPa，体积模量=100 MPa。采用FLAC3D软件分析，建立模型共包含5 394个节点和4 510个单元。边界条件为边坡前后左右4个面水平约束，底部固定约束，特征点1~3如图6所示。收敛条件为最大不平衡力与节点力平均值小于10−5。

单位：m

采用强度折减法计算该边坡安全系数，根据数值计算收敛性判据可判定边坡安全系数为2.25。折减系数−特征点位移曲线如图7所示。由图7可知：当折减系数从2.23增为2.24时，特征点1和特征点2均发生了较大的位移突变，因此采用特征点位移突变判据判定边坡的安全系数应等于2.24。图8所示为边坡的等效塑性应变云图。由图8可知：当折减系数为2.23时，塑性区尚未贯通；当折减系数为2.24时，塑性区完全贯通，因此采用塑性区贯通判据判定该边坡的安全系数应等于2.24。

1—特征点1X方向位移；2—特征点1Z方向位移；3—特征点2X方向位移；4—特征点2Z方向位移。

单位：m

折减系数−能量增量曲线如图9所示。从图9可发现：当折减系数从2.23增为2.24时，重力势能降、耗散能增量、动能增量突然增大，因此，采用重力势能降、耗散能增量、动能增量突变判据判定该边坡的安全系数应约等于2.24；当折减系数增为2.25时，弹性应变能增量突然减小，因此采用弹性应变能增量突变判据判定该边坡的安全系数应约等于2.25。

(a) 所有能量增量；(b) 弹性应变能增量；(c) 动能增量

通过前面的分析可得三维边坡的安全系数如表3所示。由表3可得：基于能量突变判据判定的安全系数与常见的3类失稳判据判定的安全系数近似相等，证明能量突变判据具有一定的适用性。

表3 不同判据判定的三维边坡安全系数

与常见失稳判据相比，能量突变判据以边坡整体作为研究对象，无需考虑边坡复杂的应力应变、强度位移和滑动面的位置和形状，只需定量计算折减过程中的能量变化，通过各能量突变来判定安全系数。该判据具有物理意义较明确、整体性较强和判断简易等优点，尤其适用于三维边坡。

3 结论

1) 边坡变形失稳的根本原因是能量的驱动作用。在强度折减过程中，输入边坡的重力势能降大部分转化为耗散能，少部分存储为弹性应变能。边坡耗散和存储能量的能力是有限的，当重力势能降大于耗散能加弹性应变能之和时，多余的能量转化为动能，表现为边坡出现滑动破坏。

2) 边坡达到临界失稳状态之前，随着折减系数增大，重力势能降、弹性应变能增量和耗散能增量均平稳增大，而动能保持为0 J。边坡达到临界失稳状态时，各能量均出现突变。其中，重力势能降和耗散能增量增长速度均突然增大，动能突然从0 J增大到一定数值，而弹性应变能增量则突然减小。

3) 通过二维边坡和三维边坡算例分析，发现利用能量突变来判定的安全系数与常见3类失稳判据得到的安全系数近似相等，表明能量突变可以作为边坡失稳判据，具有一定的适用性。

4) 能量突变判据以边坡整体作为研究对象，无需考虑边坡复杂的应力应变、强度位移和滑动面位置与形状，整体性较强，物理意义较明确，尤其适用于三维边坡。

[1] 李春光, 朱宇飞, 刘丰，等. 基于下限原理有限元的强度折减法[J]. 岩土力学, 2012, 33(6): 1816−1821.LI Chunguang, ZHU Yufei, LIU Feng, et al. Evaluation of strength reduction factor by lower bound limit analysis using finite element method[J]. Rock and Soil Mechanics, 2012, 33(6): 1816−1821.

[2] 陈国庆, 黄润秋, 石豫川, 等. 基于动态和整体强度折减法的边坡稳定性分析[J]. 岩石力学与工程学报, 2014, 33(2): 253−256.CHEN Guoqing, HUANG Runqiu, SHI Yuchuan, et al. Stability analysis of slope based on dynamic and whole strength reduction methods[J]. Chinese Journal of Rock Mechanics and Engineering, 2014, 33(2): 253−256.

[3] 陈力华, 靳晓光. 有限元强度折减法中边坡三种失效判据的适用性研究[J]. 土木工程学报, 2012, 45(9): 136−146.CHEN Lihua, JIN Xiaoguang. Study on the applicability of three criteria for slope instability using finite element strength reduction method[J]. China Civil Engineering Journal, 2012, 45(9): 136−146.

[4] 周元辅, 邓建辉, 崔玉龙, 等. 基于强度折减法的三维边坡失稳判据[J]. 岩土力学, 2014, 35(5): 1429−1437.ZHOU Yuanfu, DENG Jianhui, CUI Yulong, et al. Instability criterion of three-dimensional slope based on strength reduction method[J]. Rock and Soil Mechanics, 2014, 35(5): 1429−1437.

[5] SLOAN S W. Geotechnical stability analysis[J]. Geotechnique, 2013, 63(7): 531−572.

[6] 赵尚毅, 郑颖人, 张玉芳. 有限元强度折减法中边坡失稳的判据探讨[J]. 岩土力学, 2005, 26(2): 332−336. ZHAO Shangyi, ZHENG Yingren, ZHANG Yufang. Study on slope failure criterion in strength reduction finite element method[J]. Rock and Soil Mechanics, 2005, 26(2): 332−336.

[7] 栾茂田, 武亚军, 年廷凯. 强度折减有限元法中边坡失稳的塑性区判据及其应用[J]. 防灾减灾工程学报, 2003, 23(3): 1−8. LUAN Maotian, WU Yajun, NIAN Tingkai. A criterion for evaluating slope stability based on development of plastic zone by shear strength reduction FEM[J]. Journal of Disaster Prevention and Mitigation Engineering, 2003, 23(3): 1−8.

[8] 连镇营, 韩国城, 孔宪京. 强度折减有限元法研究开挖边坡的稳定性[J]. 岩土工程学报, 2001, 23(4): 407−411. LIAN Zhenying, HAN Guocheng, KONG Xianjing. Stability analysis of excavation by strength reduction FEM[J]. Chinese Journal of Geotechnical Engineering, 2001, 23(4): 407−411.

[9] 刘祚秋, 周翠英, 董立国, 等. 边坡稳定及加固分析的有限元强度折减法[J]. 岩土力学, 2005, 26(4): 558−561. LIU Zuoqiu, ZHOU Cuiying, DONG Liguo, et al. Slope stability and strengthening analysis by strength reduction FEM[J]. Rock and Soil Mechanics, 2005, 26(4): 558−561.

[10] MANZARI M T, NOUR M A. Significance of soil dilatancy in slope stability analysis[J]. Journal of Geotechnical and Geoenvironmental Engineering, ACSE, 2000, 123(1): 75−80.

[11] 陈国庆, 黄润秋, 周辉, 等. 边坡渐进破坏的动态强度折减法研究[J]. 岩土力学, 2013, 34(4): 1140−1146. CHEN Guoqing, HUANG Runqiu, Zhou Hui, et al. Research on progressive failure for slope using dynamic strength reduction method[J]. Rock and Soil Mechanics, 2013, 34(4): 1140−1146.

[12] 裴利剑, 屈本宁, 钱闪光. 有限元强度折减法边坡失稳判据的统一性[J]. 岩土力学, 2010, 31(10): 3337−3341. PEI Lijian, QU Benning, Qian Shanguang. Uniformity of slope instability criteria of strength reduction with FEM[J]. Rock and Soil Mechanics, 2010, 31(10): 3337−3341.

[13] 刘晓宇, 赵颖, 刘洋, 等. 土质边坡极限平衡状态及临界滑动面的判定方法[J]. 岩石力学与工程学报, 2012, 31(7): 1369−1378. LIU Xiaoyu, ZHAO Ying, LIU Yang, et al. Determination method of limit equilibrium state and critical slip surface of soil slope[J]. Chinese Journal of Rock Mechanics and Engineering, 2012, 31(7): 1369−1378.

[14] 郑颖人. 岩土数值极限分析方法的发展与应用[J]. 岩石力学与工程学报, 2012, 31(7): 1297−1316. ZHENG Yingren. Development and application of numerical limit analysis for geological materials[J]. Chinese Journal of Rock Mechanics and Engineering, 2012, 31(7): 1297−1316.

[15] 吴春秋, 朱以文, 蔡元奇. 边坡稳定临界破坏状态的动力学评判方法[J]. 岩土力学, 2005, 26(5): 784−788. WU Chunqiu, ZHU Yiwen, CAI Yuanqi. Dynamic method to assess critical state of slope stability[J]. Rock and Soil Mechanics, 2005, 26(5): 784−788.

[16] 柴红保, 曹平, 林杭, 等. 采用边坡稳定性强度折减法分析弹性应变能突变判据[J]. 中南大学学报(自然科学版), 2009, 40(4): 1054−1058. CHAI Hongbao, CAO Ping, LIN Hang, et al. Criteria of elastic strain energy in slope stability analysis using strength reduction method[J]. Journal of Central South University (Science and Technology), 2009, 40(4): 1054−1058.

[17] 谢和平, 彭瑞东. 岩石破坏的能量分析初探[J]. 岩石力学与工程学报, 2005, 24(15): 2603−2608. XIE Heping, PENG Ruidong. On energy analysis of rock failure[J]. Chinese Journal of Rock Mechanics and Engineering, 2005, 24(15): 2603−2608.

(编辑 刘锦伟)

Slope’s failure criterion based on energy catastrophe in shear strength reduction method

LIU Xinrong1, 2, TU Yiliang1, 2, ZHONG Zuliang1, 2, LIU Yongquan1, 2

(1. College of Civil Engineering, Chongqing University, Chongqing 400045, China;2. Key Laboratory of New Technology for Construction of Cities in Mountain Area(Chongqing University),Ministry of Education, Chongqing 400045, China)

In order to establish the slope’s failure criterion of energy catastrophe in finite element method with shear strength reduction (SSR−FEM), the theory for calculating the slope’s energy response was deduced based on the energy principle, and the procedure for calculating the energy response was re-developed with FLAC3D in the process of SSR−FEM. Then, the variation law of the energy response of the two-dimensional slope with the reduction factor was obtained, which was compared with the safety factors judged by the three common failure criteria. Finally, the applicability of energy catastrophe criterion was verified by a three-dimensional slope example. The results show that the gravitational potential energy decrement, kinetic energy increment and dissipated energy increment all increase suddenly when the slopes reach the state of critical failure, but the elastic strain energy increment contrariwise decreases. Thus, the energy catastrophe can be used as a slope failure criterion. The energy catastrophe criterion possesses some advantages, such as definite physical meaning, strong integrity and easy to judge. Moreover, it is especially suitable for three- dimensional slopes.

strength reduction finite element method; slope; failure criterion; energy catastrophe

10.11817/j.issn.1672-7207.2016.06.034

TU4

A

1672−7207(2016)06−2065−08

2015−06−09；

2015−08−03

国家自然科学基金资助项目(51108485，41372356)；高等学校博士学科点专项科研基金资助项目(20110191120033)；中央高校基本科研业务费科研专项基金资助项目(CDJZR12200012)(Projects(51108485, 41372356) supported by the National Natural Science Foundation of China; Project(20110191120033) supported by the Specialized Research Fund for the Doctoral Program of Higher Education of China; Project(CDJZR12200012) supported by the Fundamental Research Funds for the Central Universities)

涂义亮，博士，从事边坡工程理论、试验和数值计算研究；E-mail：tyl_ok@126.com