基于非线性Baker准则的深埋硐室拱顶稳定性分析

覃长兵; 李跃阳; 戴辰宇; 石雨沙; 仉文岗

doi:10.11779/CJGE20230662

基于非线性Baker准则的深埋硐室拱顶稳定性分析 English Version

1.
重庆大学土木工程学院, 重庆 400045
2.
东南大学土木工程学院, 江苏南京 211189

基金项目:

重庆市技术创新与应用发展专项重点项目 CSTB2022TIAD-KPX0101

详细信息

作者简介:
覃长兵(1989—)，男，博士，教授，主要从事边坡与隧道方面的教学与科研工作。E-mail：qinchangbing@cqu.edu.cn

中图分类号: TU411
计量
- 文章访问数: 238
- HTML全文浏览量: 27
- PDF下载量: 46
出版历程
- 收稿日期: 2023-07-05
- 网络出版日期: 2024-07-15
- 刊出日期: 2025-01-31

Roof stability analysis of deeply-buried cavities based on nonlinear Baker criterion

1.
School of Civil Engineering, Chongqing University, Chongqing 400045, China
2.
School of Civil Engineering, Southeast University, Nanjing 211189, China

摘要

摘要: 针对Hoek-Brown和Mohr-Coulomb准则主要分别适用于岩体和土体的情形，采用更具适用性的广义非线性Baker准则，并结合极限分析上限方法，旨在研究深埋硐室拱顶岩土体的动力稳定性。通过构建满足Baker准则下拱顶塌落的曲线型破坏机构，考虑竖向地震荷载的影响，计算内外力功率并建立功率平衡方程，基于变分原理推导在有无竖向地震力作用时拱顶破坏机制的上限表达式，从而求解得到破坏面形状、塌落拱高度和宽度的解析解。同时采用ABAQUS数值模拟验证解析解的有效性及合理性，研究分析不同参数下拱顶破坏机制的变化规律。结果表明除岩土体特性参数外，尤其是向上的竖向地震力对拱顶破坏范围有显著影响。
- 深埋硐室 /
- 拱顶稳定性 /
- 非线性Baker准则 /
- 竖向地震力 /
- 上限分析
Abstract: The dynamic roof stability analysis of deeply-buried cavities is investigated by using the upper bound limit analysis method adopting a more general nonlinear Baker criterion, in contrast to the Hoek-Brown and Mohr-Coulomb criteria which are mainly applicable to rock and soil, respectively. A curved failure mechanism for roof collapse is proposed in the realm of the Baker criterion. The vertical seismic loading is considered herein. The balance equation for work rate is then established after computing the external and internal rates of work. Based on the variational principle, the upper-bound formulation for roof collapse mechanism is derived with/without considerations of the vertical earthquake effects. Accordingly, the closed-form solutions for the failure surface, collapse height and width are explicitly obtained. At the same time, the ABAQUS modelling is used to verify the robustness and validity of closed-form solutions. The parametric studies are carried out to investigate the change laws of the roof collapse mechanism under different parameters. The results indicate that apart from rock/soil properties, the upward seismic force has a significant effect on the failure region above the cavity roof.
- deeply-buried cavity /
- roof stability /
- nonlinear Baker criterion /
- vertical seismic force /
- upper bound analysis

HTML全文

0. 引言

边坡是人类活动中极为重要的自然地质环境之一。如何保证边坡稳定成为了岩土工程一个重要的研究领域^[1]。在进行边坡稳定分析时，最主要的判断依据就是稳定安全系数。根据安全系数的大小来评判边坡是否稳定，这就涉及安全系数的定义问题。由于各种分析方法所依据的安全系数定义的不同^[2]，对于同一个计算对象，其计算出来的安全系数值一般也有所不同。但现有的一些安全系数的定义式中，很少考虑到滑动的矢量特征，且对边坡真实状态引入了过多的人为假定。为了克服这方面的局限性，葛修润院士从滑动的基本概念出发，指出滑动本身就是一个矢量变化的过程，并从矢量的角度入手，可求出滑裂面的下滑力与抗滑力，从而进一步求出边坡的安全系数。这一求解方法，称之为矢量和分析方法（vector sum method, VSM）^[3-4]。矢量和分析方法是基于力是矢量这一基本概念，安全系数定义为坡体在受各种荷载（包括自重）作用下，潜在滑面上提供的抵抗坡体滑动的合力矢与作用在潜在滑动面上贡献滑动力的合力矢在滑体的整体滑动趋势方向上投影的比值。该方法具有明确的物理意义，考虑了滑动方向，能较好地反映边坡的整体稳定性。

矢量和分析方法是根据力的矢量特性提出的一种边坡稳定性分析方法，该方法在求到了滑体的合力后，通过假定主滑方向，认为滑体沿主滑动方向平面滑动，目前常用如下4种方法来确定主滑方向：①基于潜在滑面各点真实剪应力合力矢单位方向；②基于极限下滑力的合力矢方向；③基于边坡潜在滑面剪出点指向剪入点的位置矢量方向；④基于最小势能原理。然而，在滑坡发生时，滑动面通常呈现弧形、平面、折线或曲面，这意味着坡体滑动不仅涉及平动问题，还需考虑坡体的转动，因此本文针对坡体转动提出了一种新的边坡稳定性分析方法：力矩和法（torque sum method, TSM）。此方法定义边坡稳定安全系数F为：坡体在受各种荷载（包括自重）作用条件下，潜在滑动面所能提供的抵抗坡体绕某点转动的极限强度所提供的力矩之和与滑动体所受各种荷载（包括自重）所引发坡体绕这点转动的力矩之和的比值。力矩和法可以直接采用显式计算，可以对任意滑动面上的安全系数进行求解计算，具有清晰的物理含义，概念清楚，可以直接根据滑动面的实际应力分布以及滑动面所能提供的抗滑能力进行计算，弥补了矢量和方法不能考虑坡体可能发生转动失稳的不足。本文给出了力矩和法的定义表达式，适用于一般的边坡岩土体力学分析。运用本方法计算了两个典型的边坡算例，并与传统分析方法的计算结果进行了对比分析。

1. 基于力矩和法的安全系数分析方法

1.1 安全系数的定义

对于可能发生坡体转动的情况，参照矢量和方法，边坡稳定安全系数F可以通过计算坡体应力状态获得。力矩和方法的安全系数F定义为潜在滑动面所能提供的绕某转动中心的极限抗转动力矩的总和∑M_r与作用在滑动面上计算应力绕转动中心的滑动力矩的总和∑M_s之比，即

$F=\frac{\sum \boldsymbol{M}_{\text{r}}}{\sum \boldsymbol{M}_{\text{s}}} 。$

(1)

而极限抗转动力矩是由滑动面上的极限应力强度提供。由于力是一个矢量，不但有作用的大小，而且有作用的方向，所以式(1)中求和是考虑了力的矢量性质，遵循矢量求和运算法则合成的一个力矩和。

1.2 基本假定

对于一种稳定性分析方法，若要将其所涉及到的所有因素都考虑进去，问题将会是复杂和繁琐的^[5]。如果通过适当的假设，可以让分析方法更加的有针对性，并把握住问题的主要矛盾，进而降低问题的复杂性。对力矩和法，将建立在如下假定的基础之上：

（1）对某一特定的边坡进行稳定性分析的先决条件，是通过地质勘察，确定边坡的基本物理力学参数和边坡的载荷以及边坡的力学边界条件，而坡体内应力分布状态可由其它应力分析方法（如有限元法、离散元法、条分法等）计算获得。

（2）在二维边坡问题中，确定边坡所构成的计算范围为 $\mathit{Ω}$ ，假设在该计算范围内已知边坡的潜在滑动面L，构成滑动体的区域为S。

（3）在边坡潜在滑动面上，假设岩土体的材料强度特性服从Mohr-Coulomb准则，抗剪强度为

${\tau _{\text{f}}} = c + \sigma \tan \varphi 。$

(2)

（4）对绕某点转动的力矩和法安全系数F(x_p, y_p)定义为：转动点p坐标为(x_p, y_p)，滑动面上所提供抵抗坡体转动的各力沿转动点的力矩代数和∑M_r(x_p, y_p)与滑动面上计算应力沿转动点的力矩代数和∑M_s(x_p, y_p)的比值：

$F({x_{\text{p}}},{y_{\text{p}}}) = \frac{{\sum {{M_{\text{r}}}({x_{\text{p}}},{y_{\text{p}}})} }}{{\sum {{M_{\text{s}}}({x_{\text{p}}},{y_{\text{p}}})} }} 。$

(3)

（5）滑动面是滑动体与母体脱离发生绕某点转动的面，可以是任意形状。对于滑动（转动）的滑体，滑动面是约束边界。做力矩和计算时，滑动体假定为刚体。

1.3 滑动体转动中心的确定

当分析某一具体的潜在滑动体时，也即确定了其对应的潜在滑动面。作为平面刚体，刚体最基本的运动形式为沿某直线的平动和绕某点的转动，如果不考虑其后的运动轨迹，从静止时的运动趋势考虑，可以认为刚体只有平动和转动，而且其启动加速度在理论上是可以计算的。

对于边坡滑体，以往的知识认为滑体是沿滑动面滑动（或转动），笔者认为滑动面是边坡滑体失稳后，产生运动形成的轨迹，与边坡初始失稳时形成的失稳面（或破坏面）是有区别的，如果将滑动面的形式与位移或运动的形式相关联，则无法解释现有的分析滑动面为非圆弧或非直线的方法中，滑动体作为刚体时如何发生。当边坡初始失稳时，无论滑动面形式如何，滑动面只是滑动体的一个边界，其相互作用力可以认为是边界上的力边界条件。

因为刚体的基本位移是平动和转动，分析滑动体时，应该按照其基本的位移特征进行分析。矢量和方法考虑了滑动体的平动位移，认为滑动体将沿某一直线运动，并通过力矢量求和等相应的假设，获得其安全系数。作为对矢量和方法的一个补充，本文将只考虑滑动体的转动。对于同时考虑平动和转动的复杂情况，将在往后的工作中开展。

值得强调的是，这里所指的滑动体的平动和转动，均与滑动面无关，仅与滑动体本身和其受力有关，平动和转动也只是指滑动体在初始失稳时的位移趋势，与滑动体失稳后的位移和滑动面的形成无直接关系。实质上，边坡的失稳是指的滑动体在初始位移趋势上的失稳，并在考虑初始位移趋势上计算其失稳安全系数。

如何获得滑动体的初始位移趋势或初始转动中心，可以按照边坡稳定性分析的普遍方法，也即求解边坡最小安全系数的方法来获得。定义其初始位移趋势（或转动中心）是使滑动体具有最小安全系数时对应的转动中心。

仅考虑滑动体转动时，如图 1为边坡稳定性分析时的安全系数求解示意图。中可以通过有限元法或其它方法，计算得到潜在滑动面上任何一点i的正应力 ${\sigma _i}$ 和切应力 ${\tau _i}$ 。设滑动面上点i处岩土体内聚力为 ${c_i}$ ，内摩擦角为 ${\varphi _i}$ ，中 ${\alpha _i}$ 为滑动面在i点上的切线相对于坐标系x轴的正向夹角。

图 1 转动点搜索示意图

Figure 1. Schematic diagram of rotating point searching

下载: 全尺寸图片幻灯片

设点p(x_p, y_p)为转动点，为了寻找所计算滑动体具有最小安全系数的转动点，本文采用了较直观的网格（Grid）搜索法，转动点假设位于设定的一矩形搜索区域，矩形区域均布网格节点，作为试算的转动点。每一节点只需计算一次安全系数，就可通过比较获得区域内具有安全系数最小的节点。如果获得的节点位于矩形区域边界上，则需调整区域。该方法与圆弧滑动面用Grid方法搜索滑动圆心类似（图 2）。

图 2 安全系数求解示意图

Figure 2. Solution of factor of safety

下载: 全尺寸图片幻灯片

1.4 力矩和法的安全系数计算

力矩和法安全系数的计算包括转动力矩和计算与抗转动力矩和计算两方面：

（1）转动力矩和计算

滑动体上任一微段 ${{\Delta }}{l_i}$ 的正应力 ${\sigma _i}$ 和切应力 ${\tau _i}$ 是由滑体自重和作用在滑体上的外载荷所产生的，理论上滑动面（也是边界）上的所有正应力 ${\sigma _i}$ 和切应力 ${\tau _i}$ 与滑动体的自重和外荷载是一个静力平衡系。正应力 ${\sigma _i}$ 和切应力 ${\tau _i}$ 对转动中心求力矩并沿滑动面求和，就是滑体的转动力矩和。同时，岩床反作用在滑体上的力为 ${\sigma '_i}$ 和 ${\tau '_i}$ ，它们数值相等，方向相反。如果能够求到滑动体的重心，转动力矩和的求解可以直接求重心对转动中心的力矩与滑体上的外荷载对转动中心的力矩的总和。这与沿滑动面由 ${\sigma _i}$ 和 ${\tau _i}$ 求力矩和是一致的。

将滑动面上的应力对转动点力矩沿滑动面积分（或求和）即为式（3）中的分母项。

滑动面上任何一点i的切应力 ${\tau _i}$ 的斜率为 ${k_{\tau i}}$ ，所在直线方程为

$y={k}_{\tau i}x+{b}_{i} 。$

(4)

式中： ${k_{\tau i}} = \tan {\alpha _i}$ ， ${b_i} = {y_i} - {k_{\tau i}}{x_i}$ ，(x_i, y_i)为i点坐标。

转动点(x_p, y_p)到 ${\tau _i}$ 所在直线的距离为

${d}_{\tau i}=\frac{\mid {k}_{\tau i}{x}_{\text{p}}+\text{(}-1\text{)}{y}_{\text{p}}+{b}_{i}\mid }{\sqrt{{k}_{\tau i}{}^{2}+{\text{(}-1\text{)}}^{2}}} 。$

(5)

切应力 ${\tau _i}$ 提供转动力矩，力矩和为 $\sum _{i=1}^{n}{\tau }_{i}\text{Δ}{l}_{i}{d}_{\tau i}$ 。

滑动面上任何一点i的正应力 ${\sigma _i}$ 的斜率为 ${k_{\sigma i}}$ ，所在直线方程为

$y={k}_{\sigma i}x+{B}_{i} 。$

(6)

其中， ${k_{\sigma i}} = - \frac{1}{{{k_{\tau i}}}} = - \frac{1}{{\tan {\alpha _i}}}$ ， ${B_i} = {y_i} - {k_{\sigma i}}{x_i}$ ，(x_i, y_i)为i点坐标。

转动点(x_p, y_p)到 ${\sigma _i}$ 所在直线的距离为

${d}_{\sigma i}=\frac{\mid {k}_{\sigma i}{x}_{\text{p}}+(-1){y}_{\text{p}}+{B}_{i}\mid }{\sqrt{{k}_{\sigma i}{}^{2}+{(-1)}^{2}}} 。$

(7)

正应力 ${\sigma _i}$ 提供转动力矩，力矩和为 $\sum\limits_{i=1}^{n}{\sigma }_{i}\text{Δ}{l}_{i}{d}_{\sigma i}$ 。

滑动体的总的转动力矩和为

$\sum\limits_{i=1}^{n}{M}_{\text{s}}({x}_{\text{p}},{y}_{\text{p}}) ={\displaystyle \sum _{i=1}^{n}{\tau }_{i}\text{Δ}{l}_{i}{d}_{\tau i}} +\sum\limits_{i=1}^{n}{\sigma }_{i}\text{Δ}{l}_{i}{d}_{\sigma i} 。$

(8)

由于转动点的位置不同，滑动面的形状可能不规则，某一弧段的 ${\sigma _i}$ 可能提供下滑力矩也可能提供抗滑力矩，体现其计算力矩时可正可负（如果定义沿坡面倾向方向转动为正），但作为力矩和方法，不需要区分此时的应力是否是产生下滑力矩或抗滑力矩，只需将抗滑力矩设为负值，进行力矩和计算即可。在计算切应力 ${\tau _i}$ 时也如此，不需要区分，只需根据实际发生的计算应力进行计算即可。不需要区分任意弧段的 ${\sigma _i}$ 和 ${\tau _i}$ 产生的力矩的正负，更符合“和”方法的本意，而且在力学上，滑动面上的 ${\sigma _i}$ 和 ${\tau _i}$ ，等效于滑动体的重力和作用于滑动体上的全部外力，求和时结果是一致的。

（2）抗转动力矩和计算

岩土体的材料强度服从Mohr-Coulomb准则，滑动体沿滑动面发生剪切，其抗转动力矩则由滑动面上的切向抗剪强度和法向反力提供。对于滑动面上的任一微段 ${\Delta }{l_i}$ ，岩土体的抗剪强度 ${\tau _{{\text{f}}i}}$ 由莫尔库仑强度理论进行计算，方向与位移趋势方向相反，法向反力 ${\sigma '_i}$ = ${\sigma _i}$ ，方向相反。

${\tau _{{\text{f}}i}} = {c_i} + {\sigma _i}\tan {\varphi _i}。$

(9)

滑动面上切向抗剪强度和法向反力对转动点的力矩积分（或求和）构成了抗转动力矩和，即为式（3）中的分子项。

${\tau _{{\text{f}}i}}$ 所在直线与 ${\tau _i}$ 一致，故转动点(x_p, y_p)到 ${\tau _{{\text{f}}i}}$ 所在直线的距离也为 ${d_{\tau i}}$ ，提供抗转动力矩之和为 $\sum _{i=1}^{n}{\tau }_{\text{f}i}\text{Δ}{l}_{i}{d}_{\tau i}$ 。

转动点(x_p, y_p)到 ${\sigma '_i}$ 所在直线的距离也为 ${d_{\sigma i}}$ ，提供抗转动力矩之和为 $\sum\limits_{i=1}^{n}{{\sigma }^{\prime }}_{i}\text{Δ}{l}_{i}{d}_{\sigma i}$ 。

抗转动力矩和即为

$\sum\limits_{i=1}^n M_{\mathrm{r}}\left(x_{\mathrm{p}}, y_{\mathrm{p}}\right)=\sum\limits_{i=1}^n \tau_{\mathrm{f}i} \Delta l_i d_{\tau i}+\sum\limits_{i=1}^n \sigma_i^{\prime} \Delta l_i d_{\sigma i}。$

(10)

（3）安全系数的计算

由式（8），（10），并代入式（3），得到力矩和法安全系数的表达式（1）或（3）的计算式。

$F({x_{\text{p}}},{y_{\text{p}}}) = \frac{{\sum ({c_i} + {\sigma _i}\tan {\varphi _i}){\Delta }{l_i}{d_{\tau i}} + \sum {{\sigma '}_i}{\Delta }{l_i}{d_{\sigma i}}}}{{\sum {\tau _i}{\Delta }{l_i}{d_{\tau i}} + \sum {\sigma _i}{\Delta }{l_i}{d_{\sigma i}}}}\text{，}$

(11)

写成积分形式为

$F({x_{\text{p}}},{y_{\text{p}}}) = \frac{{\int {\left[ {({c_i} + {\sigma _i}\tan {\varphi _i}){d_{\tau i}} + {{\sigma '}_i}{d_{\sigma i}}} \right]{\text{d}}l} }}{{\int {({\tau _i}{d_{\tau i}} + {\sigma _i}{d_{\sigma i}}){\text{d}}l} }}。$

(12)

1.5 力矩和方法求解步骤

力矩和方法的基本计算步骤如下：

（1）根据边坡地质勘查资料，确定边坡的潜在滑动面和对应的滑动体。

（2）建立考虑所分析边坡和滑动面的力学（应力）分析模型（如有限元、离散元、条分法等可以获得滑动面应力状态的方法均可），根据边坡的载荷、边界条件和岩土体物理力学参数，得到潜在滑动面的计算应力分布状态。

（3）将潜在滑动面分为若干段，从分析模型中获得每段编号、位置坐标和几何参数，从计算所得应力分布状态获得每滑面段的应力。

（4）根据滑动面位置，拟定试算的转动点位置区域，依次以区域内节点为转动点，根据式（8），（10）计算滑动面上的转动力矩和与抗转动力矩和，再按式（11）或（12）计算安全系数。

（5）对所试算区域的所有转动节点，通过计算找到最小安全系数，如果最小安全系数对应节点位于区域内部，可以认为此时的最小安全系数为所分析滑动面的基于力矩和法分析的安全系数，及其对应的转动节点。若最小安全系数对应节点位于区域边界，则回到第（4）步骤，调整试算的转动点位置区域。

2. 基于条分法模型的力矩和方法求解

力矩和法具有物理含义明确、概念清晰的特点，它可以直接采用显式计算，在已经确定了滑动面位置前提下，只要得到边坡在荷载和自重下的计算应力场，就可以通过力矩和法计算边坡的稳定性安全系数。为了获得滑动面上的应力，力矩和法是需要结合其它分析方法的。在边坡稳定性分析方法中，条分法是目前应用最广泛的分析方法之一。为了获得滑动面的应力，本文将采用条分法计算模型，通过条分法的计算获得滑动面上的应力，从而进行力矩和方法的计算。基于条分法的滑动面应力求解，其一便于方法之间的对比特别是与广泛应用的条分法的结果之间的对比分析；其二是更易于实现，避免了有限元计算时对模型和程序的较多依赖。

条分法中，Morgenstern-Price方法（M-P法）是一种同时满足力和力矩平衡条件的通用条分法^[6]，它是边坡极限平衡条分法^[7-9]中最具普遍性和假设条件的一种分析方法。下面将介绍M-P法与力矩和法相结合，求解安全系数。

如图 3所示，将滑体划分成n个竖直条块，取任意条块i，其中条块i的底面长度是 ${l_i}$ ，底边与X轴的夹角是 ${\alpha _i}$ ，转动点坐标是(x_p, y_p)，条块i底面左下脚坐标是(x_i1, y_i1)，条块i底面右下脚坐标是(x_i, y_i)。W_i是条块i的重力，E_i1，E_i分别是条块i左右两侧的水平条间力，X_i1，X_i分别是条块i左右两侧的条间剪力，h_i-1为作用在条块i左侧水平条间力E_i-1离基底的距离，h_i为作用在条块i右侧水平条间力E_i离基底的距离，H是坡高，β是坡角，P_i是底部的法向力, T_i是底部的切向剪力, u_i是底部的孔隙水压， ${\varphi '_i}$ 是滑动表面的有效内摩擦角，有效黏聚力为 ${c'_i}$ ，绕点(x_p, y_p)转动的安全系数F(x_p, y_p)。

图 3 条分法安全系数求解示意图

Figure 3. Schematic diagram of solving factor of safety by slices method

下载: 全尺寸图片幻灯片

M-P法假定条间力函数关系为 ${X_i} = \lambda f(x){E_i}$ 。式中 $\lambda$ 为任意常数， $f(x)$ 可取半正弦函数、常数等。 $f(x) =$ ${\text{sin}}\left[ {(x - a)/(b - a)} \right]{\text{π }}$ （正弦曲线），其中a，b分别为滑动面两端的x坐标。如果 $f(x)$ =1，则条间力( ${X_i}$ 和 ${E_i}$ 的合力)方向平行，即为斯宾塞法对条间力的假设。当力与力矩的平衡被严格满足时，在一个合理的范围内，对条间力的不同假设对其安全系数的影响不大^[10]。

根据莫尔-库仑准则及安全系数F的定义，建立条块的水平和垂直方向力平衡、滑体整体水平向力平衡以及将作用在条块上的力对土条条块底中点取矩建立力矩平衡。为了便于计算力矩和安全系数，得到 ${P_i}$ 和 ${T_i}$ 的表达式，有如下公式：

${E}_{i}={E}_{i-1}+ {P}_{i}\mathrm{sin}{\alpha }_{i}-{T}_{i}\mathrm{cos}{\alpha }_{i} \text{，}$

(13)

${P}_{i}=\frac{{W}_{i}-({X}_{i}-{X}_{i-1})-\frac{1}{F}({{c}^{\prime }}_{i}{l}_{i}-{u}_{i}{l}_{i}\mathrm{tan}{{\varphi }^{\prime }}_{i})}{\mathrm{cos}{\alpha }_{i}+\frac{1}{F}(\mathrm{sin}{\alpha }_{i}\mathrm{tan}{{\varphi }^{\prime }}_{i})} \text{，}$

(14)

${T}_{i}=\frac{{{c}^{\prime }}_{i}{l}_{i}+({P}_{i}-{u}_{i}{l}_{i})\mathrm{tan}{{\varphi}^{\prime }}_{i}}{F} \text{，}$

(15)

$\sum {({P_i}\sin {\alpha _i} - {T_i}\cos {\alpha _i})} = 0 \text{，}$

(16)

${E_{i - 1}}\left( {{h_{i - 1}} + {y_{i - 1}} - \frac{{{y_{i - 1}} + {y_i}}}{2}} \right) - {E_i}({h_i} - \frac{{{y_{i - 1}} + {y_i}}}{2} + {y_i}) -\\ ~~~~~~~~~~ {X}_{i-1}\left(\frac{{x}_{i-1}+{x}_{i}}{2}-{x}_{i-1}\right)-{X}_{i}({x}_{i}-\frac{{x}_{i-1}+{x}_{i}}{2})\text{=}0\text{ }。$

(17)

岩土体的抗剪强度 ${\tau _{{\text{f}}i}}$ 提供抗滑力，其值由莫尔库仑强度准则进行计算，滑动面上条块i的 ${\tau _{{\text{f}}i}}$ 的斜率为 ${k_{\tau {\text{f}}i}}$ ，所在直线方程为

$y={k}_{\tau \text{f}i}x+{b}_{i} \text{，}$

(18)

式中： ${k_{\tau {\text{f}}i}} = \frac{{{y_i} - {y_{i - 1}}}}{{{x_i} - {x_{i - 1}}}}$ ； ${b_i} = {y_i} - \frac{{{y_i} - {y_{i - 1}}}}{{{x_i} - {x_{i - 1}}}}{x_i}$ 。

转动点(x_p, y_p)到 ${\tau _{{\text{f}}i}}$ 所在直线的距离为

${d}_{\tau \text{f}i}=\frac{\mid {k}_{\tau \text{f}i}{x}_{\text{p}}+(-1){y}_{\text{p}}+{b}_{i}\mid }{\sqrt{{k}_{\tau \text{f}i}{}^{2}+{(-1)}^{2}}} 。$

(19)

抗剪强度 ${\tau _{{\text{f}}i}}$ 提供抗转动力矩，之和为

$\sum {\tau _{{\text{f}}i}}{l_i}{d_{\tau {\text{f}}i}} = \sum [{c'_i}{l_i} + ({P_i} - {u_i}{l_i})\tan {\varphi '_i}]{d_{\tau {\text{f}}i}} 。$

(20)

式中： ${l_i}$ 为条块i的基底长度。

滑动面上条块i的基底反力 ${P}_{i}$ 的斜率为 ${k_{Pi}}$ ，所在直线方程为

$y={k}_{Pi}x+{B}_{i} 。$

(21)

式中： ${k_{Pi}} = - \frac{1}{{{k_{\tau {\text{f}}i}}}} = - \frac{{{x_i} - {x_{i - 1}}}}{{{y_i} - {y_{i - 1}}}}$ ； ${B_i} = \frac{{{y_i} + {y_{{\text{i}} - 1}}}}{2} - {k_{Pi}} \cdot$ $\frac{{{x_i} + {x_{i - 1}}}}{2}$ 。

转动点(x_p, y_p)到 ${P}_{i}$ 所在直线的距离为

${d}_{Pi}=\frac{\mid {k}_{Pi}{x}_{\text{p}}+(-1){y}_{\text{p}}+{B}_{i}\mid }{\sqrt{{k}_{Pi}{}^{2}+{(-1)}^{2}}} \text{，}$

(22)

${P}_{i}$ 提供的抗转动力矩和为 $\sum {{P_i}{d_{Pi}}}$ 。

抗转动力矩和即为

$\sum {M}_{\text{r}}\text{(}{x}_{\text{p}},{y}_{\text{p}}\text{)}=\sum {\tau }_{\text{f}i}{l}_{i}{d}_{\tau \text{f}i}+\sum {P}_{i}{d}_{Pi}。$

(23)

${W_i}$ 提供转动力矩，之和为 $\sum {{W_i}{d_{Wi}}}$ 。其中 ${d_{Wi}} = {x_{\text{p}}} -$ ${x_i}$ 。

力矩和安全系数表达式如下：

$F{\text{(}}{x_{\text{p}}},{y_{\text{p}}}{\text{)}} = \frac{{\sum {M_{\text{r}}}{\text{(}}{x_{\text{p}}},{y_{\text{p}}}{\text{)}}}}{{\sum {M_{\text{s}}}{\text{(}}{x_{\text{p}}},{y_{\text{p}}}{\text{)}}}} = \frac{{\sum {\tau _{{\text{f}}i}}{l_i}{d_{\tau fi}} + \sum {P_i}{d_{Pi}}}}{{\sum {W_i}{d_{Wi}}}}。$

(24)

3. 算例

3.1 圆弧滑动面边坡算例

选用ACADS考题1中EX1(c)，此边坡为非均质边坡，由3层土组成，其计算模型见图 4，材料参数见表 1。

图 4 ACADS考题EX1(c)计算模型

Figure 4. Model of problem of EX1(c) of ACADS

下载: 全尺寸图片幻灯片

表 1 ACADS考题EX1(c)材料参数

Table 1. Material parameters of EX1(c) of ACADS

土层	c/kPa	φ/(°)	γ/(kN·m^-3)
1	0	38	19.5
2	5.3	23	19.5
3	7.2	20	19.5

下载: 导出CSV

| 显示表格

本文暂时不对边坡稳定力矩和法的临界滑动面的搜索问题进行讨论，图 4中的滑动面为边坡稳定性分析软件slide求得的临界滑动面，本算例用此滑动面计算边坡稳定的力矩和法安全系数，其中考题EX1(c)的滑动面同时穿过3种土层。

表 2为条分法计算的条块数据。表 3为多种方法计算结果。计算均采用同一个滑动面，从表中可以看出，本文方法计算得到的安全系数最小，其值为1.043，且对于圆弧滑动面，当转动点为圆心时，本文方法计算结果为1.378，与裁判答案基本一致，误差小于2%。

表 2 条方法计算过程

Table 2. Calculation process by slices method

条块i	基底长l_i/ m	条重W_i/ (kN·m^-3)	抗滑力 ${\tau _{{\text{f}}i}}{l_i}$ /(kN·m^-1)	基底反力P_i/ (kN·m^-1)
1	5.23	110.11	75.25	96.31
2	3.97	232.52	143.43	183.58
3	3.30	261.32	110.98	220.29
4	3.35	284.60	119.98	263.44
5	3.19	243.01	117.80	260.58
6	3.13	167.48	95.96	201.65
7	3.17	59.87	50.78	76.88

下载: 导出CSV

| 显示表格

表 3 安全系数计算结果

Table 3. Calculation results of safety factor

计算方法	圆心坐标/m	转动点坐标/m	安全系数
费伦纽斯法	(35.77, 23.54)		1.263
毕肖普法	(35.77, 23.54)		1.372
简布法	(35.77, 23.54)		1.296
M-P法	(35.77, 23.54)		1.379
矢量和法			1.375
力矩和法		(33.5, 4.6)	1.043

下载: 导出CSV

| 显示表格

3.2 非圆弧滑面边坡算例

非圆弧滑面边坡算例选用ACADS提供的考题3中EX3(b)，EX3(b)是一个含软弱夹层的非均质边坡考核题目，边坡计算模型见图 5，岩土参数见表 4。

图 5 ACADS考题EX3(b)计算模型

Figure 5. Model of problem of EX3(b) of ACADS

下载: 全尺寸图片幻灯片

表 4 ACADS考题EX3(b)材料参数

Table 4. Material parameters of EX3(b) of ACADS

土层	c/kPa	φ/(°)	γ/(kN·m^-3)
1	28.5	20	18.84
2	0	10	18.84

下载: 导出CSV

| 显示表格

对于具有软弱夹层的边坡，其滑面必定有一段会沿软弱夹层滑动，此时若依然采用圆弧滑动面，得到的结果往往是不正确的，故这里采用指定的折线滑动面进行计算边坡稳定的力矩和法安全系数。图 5给出剖面几何特征和指定滑动面位置，将滑动体划分为7个条块，滑动面同时穿过2种土层，底部一段沿软弱夹层滑动，滑动面控制点坐标见表 5。

表 5 滑动面控制点坐标

Table 5. Coordinates of control points for sliding surface

点号	X/m	Y/m
1	10.69	20.00
2	20.50	7.00
3	40.00	6.50
4	42.15	7.75

下载: 导出CSV

| 显示表格

表 6为条分法计算的条块数据。表 7为多种方法计算结果。

表 6 条分法计算过程

Table 6. Calculation process by slices method

条块i	基底长l_i/m	条重W_i/ (kN·m^-3)	抗滑力 ${\tau _{{\text{f}}i}}{l_i}$ / (kN·m^-1)	基底反力P_i/ (kN·m^-1)
1	8.14	300.33	299.76	185.96
2	8.14	825.64	491.09	711.65
3	6.50	1158.27	207.78	1178.36
4	6.50	780.68	144.07	817.07
5	6.50	403.10	75.24	426.69
6	0.99	20.82	5.24	29.71
7	1.49	9.21	56.75	39.08

下载: 导出CSV

| 显示表格

表 7 安全系数计算结果

Table 7. Calculated results of factor of safety

计算方法		转动点坐标	安全系数
简布法			1.195
陆军工程师团法			1.364
M-P法			1.286
矢量和法			1.221
力矩和法	条间力函数为常数	(40, 6.5)	1.051
力矩和法	条间力函数为半正弦函数	(40, 6.5)	1.049

下载: 导出CSV

| 显示表格

从表 7可以看出，在非圆弧滑面的情况下，本文方法得到的安全系数较小，而不同的条间力假设对力矩和安全系数的影响并不明显，而对于含软弱夹层的非均质边坡，滑面偏离圆弧较大，本文方法的优越性是能够处理各种非圆弧滑面。

4. 结论

（1）提出了分析边坡稳定性的新方法即力矩和方法。该方法对一确定的滑坡体和滑动面，假定滑坡体绕某点转动，求出滑坡体在自重和外荷载作用下的转动力矩和，同时求出滑动面上基岩法向反力和剪切强度对滑动体的抗转动力矩和，定义抗转动力矩和与转动力矩和的比值为对应的安全系数。本文中，转动点由最小安全系数决定，并通过试算搜索得到。

（2）力矩和方法是对矢量和方法的一个补充。滑动体视为刚体时，它的位移应该包括刚体的平动和转动两部分（只有外力合力通过重心时没有外力引起的转动），矢量和方法认为滑动体是沿某一直线方向平动，不考虑其转动，而力矩和方法考虑了转动，两个方法互为补充。当然，若能同时考虑滑体的平动和转动，将会更合理，也是本文后续将延续的工作。

（3）本文在求解滑动面应力时，采用了更简便的条分法计算（M-P法）。并结合M-P法模型，推导了力矩和方法求解公式。

（4）算例1结果显示，对于所分析的滑动面，由于所选择的转动点不同，力矩和方法所得到的安全系数较其它几种方法更小。这主要是因为最终所得最小安全系数对应的转动点不是所分析滑动面圆弧的圆心点。滑动面圆弧的圆心点只是力矩和方法所分析转动点的一个，其对应的安全系数不一定是滑动体绕某点发生转动时的最小安全系数。

这一结果也提示一点，具有圆弧滑动面的滑坡体，其初始失稳转动时所对应的转动点（或圆心）不一定是圆弧滑动面对应的圆心。进一步可认为，任意形状的滑动面，对于滑动体，滑动面只是其约束边界条件，与滑动体初始失稳或破坏的位移趋势无直接关系。这与人们习惯的认知是不同的，值得后续进一步的研究。

算例1中，当转动点取为圆心时，本文方法求得的安全系数与裁判答案基本一致，说明力矩和方法计算安全系数是可行的。

（5）力矩和方法可以分析任意滑动面。滑动体所对应的滑动面，只是滑动体的边界条件，其形状不会直接影响滑动体作为刚体自身的位移趋势，而是通过边界提供反力影响滑动体的后续位移轨迹。算例2给出了本文方法分析非圆弧滑动面的情况，并可以搜索到比原结果较小的安全系数，说明假定滑动体绕某点转动的位移趋势是值得对边坡进行稳定性分析时考虑。

图 1 广义非线性Baker准则在Mohr平面τ_n-σ_n

Figure 1. Generalized nonlinear Baker criterion in Mohr's plane

下载: 全尺寸图片幻灯片

图 2 地震作用下Baker岩土体硐室拱顶坍塌机制

Figure 2. Collapse mechanism of cavity crown under earthquake and Baker rock/soil medium

下载: 全尺寸图片幻灯片

图 3 深埋岩土体隧道有限元模型

Figure 3. Finite element model for deeply-buried tunnel

下载: 全尺寸图片幻灯片

图 4 无竖向地震荷载下数值解与解析解对比

Figure 4. Comparison of numerical and analytical solutions under no vertical earthquake

下载: 全尺寸图片幻灯片

图 5 竖向地震荷载下数值解与解析解对比

Figure 5. Comparison of numerical and analytical solutions under vertical earthquake effects

下载: 全尺寸图片幻灯片

图 6 参数A影响下的塌落形状及范围

Figure 6. Collapse shapes and dimensions under varying values of A

下载: 全尺寸图片幻灯片

图 7 参数n影响下的塌落形状及范围

Figure 7. Collapse shapes and dimensions under varying values of n

下载: 全尺寸图片幻灯片

图 8 参数T影响下的塌落形状及范围

Figure 8. Collapse shapes and dimensions under varying values of T

下载: 全尺寸图片幻灯片

图 9 参数γ影响下的塌落形状及范围

Figure 9. Collapse shapes and dimensions under varying values of γ

下载: 全尺寸图片幻灯片

图 10 参数k_v影响下的塌落形状及范围

Figure 10. Collapse shapes and dimensions under varying values of k_v

下载: 全尺寸图片幻灯片

图 11 工程地质横断图

Figure 11. Diagram of engineering geological cross section

下载: 全尺寸图片幻灯片

表 1 有限元模型几何参数

Table 1 Geometric parameters of finite element model

H/m B/m H₀/m L_c/m

40 50 5 10

下载: 导出CSV

表 2 不同破坏准则下的岩土体参数值

Table 2 Rock/soil parameters under different failure criteria

γ/
(kN·m^-3) Baker强度参数 HB强度参数 MC强度

A n T p_a/
kPa A B σ_t/
kPa σ_c/
kPa c/
kPa ϕ/
(°)

25 2.08 0.7 0.3 100 0.75 0.7 30 3000 115 51.6

下载: 导出CSV

表 3 隧道围岩力学参数

Table 3 Mechanical parameters of surrounding rock of tunnel

γ/(kN·m^-3) c/kPa ϕ/(°) A n T

18 30 24 0.45 0.7 0.67

下载: 导出CSV

表 4 围岩压力结果对比

Table 4 Comparative results of pressure of surrounding rock

监测点桩号埋深/
m 本文结果/kPa 规范计算结果/
kPa 误差/
%

ZK3+610 41.41 118.63 113.66 4.4

下载: 导出CSV

参考文献(22)

[1]	CHEN W F. Limit Analysis and Soil Plasticity[M]. Amsterdam; New York: Elsevier Scientific Pub. Co, 1975.
[2]	MICHALOWSKI R L. Slope stability analysis: a kinematical approach[J]. Géotechnique, 1995, 45(2): 283-293. doi: 10.1680/geot.1995.45.2.283
[3]	赵炼恒, 李亮, 杨峰, 等. 加筋土坡动态稳定性拟静力分析[J]. 岩石力学与工程学报, 2009, 28(9): 1904-1917. doi: 10.3321/j.issn:1000-6915.2009.09.023 ZHAO Lianheng, LI Liang, YANG Feng, et al. Dynamic stability pseudo-static analysis of reinforcement soil slopes[J]. Chinese Journal of Rock Mechanics and Engineering, 2009, 28(9): 1904-1917. (in Chinese) doi: 10.3321/j.issn:1000-6915.2009.09.023
[4]	SLOAN S W. Geotechnical stability analysis[J]. Géotechnique, 2013, 63(7): 531-571. doi: 10.1680/geot.12.RL.001
[5]	孙志彬, 潘秋景, 杨小礼, 等. 非均质边坡上限分析的离散机构及应用[J]. 岩石力学与工程学报, 2017, 36(7): 1680-1688. SUN Zhibin, PAN Qiujing, YANG Xiaoli, et al. Discrete mechanism for upper bound analysis of nonhomogeneous slopes[J]. Chinese Journal of Rock Mechanics and Engineering, 2017, 36(7): 1680-1688. (in Chinese)
[6]	QIN C B, CHIAN S C. Kinematic analysis of seismic slope stability with a discretisation technique and pseudo-dynamic approach: a new perspective[J]. Géotechnique, 2018, 68(6): 492-503. doi: 10.1680/jgeot.16.P.200
[7]	FRALDI M, GUARRACINO F. Limit analysis of collapse mechanisms in cavities and tunnels according to the Hoek–Brown failure criterion[J]. International Journal of Rock Mechanics and Mining Sciences, 2009, 46(4): 665-673. doi: 10.1016/j.ijrmms.2008.09.014
[8]	FRALDI M, GUARRACINO F. Analytical solutions for collapse mechanisms in tunnels with arbitrary cross sections[J]. International Journal of Solids and Structures, 2010, 47(2): 216-223. doi: 10.1016/j.ijsolstr.2009.09.028
[9]	FRALDI M, GUARRACINO F. Evaluation of impending collapse in circular tunnels by analytical and numerical approaches[J]. Tunnelling and Underground Space Technology, 2011, 26(4): 507-516. doi: 10.1016/j.tust.2011.03.003
[10]	YANG X L, QIN C B. Limit analysis of rectangular cavity subjected to seepage forces based on Hoek-Brown failure criterion[J]. Geomechanics and Engineering, 2014, 6(5): 503-515. doi: 10.12989/gae.2014.6.5.503
[11]	YANG X L, HUANG F. Three-dimensional failure mechanism of a rectangular cavity in a Hoek-Brown rock medium[J]. International Journal of Rock Mechanics and Mining Sciences, 2013, 61: 189-195. doi: 10.1016/j.ijrmms.2013.02.014
[12]	QIN C B, LI Y Y, YU J, et al. Closed-form solutions for collapse mechanisms of tunnel crown in saturated non-uniform rock surrounds[J]. Tunnelling and Underground Space Technology, 2022, 126: 104529. doi: 10.1016/j.tust.2022.104529
[13]	BAKER R. Nonlinear Mohr envelopes based on triaxial data[J]. Journal of Geotechnical and Geoenvironmental Engineering, 2004, 130(5): 498-506. doi: 10.1061/(ASCE)1090-0241(2004)130:5(498)
[14]	ZHANG D B, MA Z Y, YU B, et al. Upper bound solution of surrounding rock pressure of shallow tunnel under nonlinear failure criterion[J]. Journal of Central South University, 2019, 26(7): 1696-1705. doi: 10.1007/s11771-019-4126-3
[15]	刘智振. 非线性Baker破坏准则下地下硐室围岩压力上限解研究[D]. 湘潭: 湖南科技大学, 2017. LIU Zhizhen. Study on Upper Bound Solution of Surrounding Rock Pressure in Underground Cavity under Nonlinear Baker Failure Criterion[D]. Xiangtan: Hunan University of Science and Technology, 2017. (in Chinese)
[16]	HOEK E, BROWN E T. Empirical strength criterion for rock masses[J]. Journal of the Geotechnical Engineering Division, 1980, 106(9): 1013-1035. doi: 10.1061/AJGEB6.0001029
[17]	HOEK E, CARRANZA-TORRES C, CORKUM B. Hoek-Brown failure criterion: 2002 edition[C]// Proceedings of the North American Rock Mechanics Symposium. Toronto, 2002.
[18]	XU J S, YANG X L. Seismic stability analysis and charts of a 3D rock slope in Hoek-Brown media[J]. International Journal of Rock Mechanics and Mining Sciences, 2018, 112: 64-76. doi: 10.1016/j.ijrmms.2018.10.005
[19]	黄阜, 杨小礼, 赵炼恒, 等. 基于Hoek-Brown破坏准则的浅埋条形锚板抗拔力上限分析[J]. 岩土力学, 2012, 33(1): 179-184, 190. doi: 10.3969/j.issn.1000-7598.2012.01.028 HUANG Fu, YANG Xiaoli, ZHAO Lianheng, et al. Upper bound solution of ultimate pullout capacity of strip plate anchor based on Hoek-Brown failure criterion[J]. Rock and Soil Mechanics, 2012, 33(1): 179-184, 190. (in Chinese) doi: 10.3969/j.issn.1000-7598.2012.01.028
[20]	JIANG J C, BAKER R, YAMAGAMI T. The effect of strength envelope nonlinearity on slope stability computations[J]. Canadian Geotechnical Journal, 2003, 40(2): 308-325. doi: 10.1139/t02-111
[21]	LIU Z Z, CAO P, LIN H, et al. Three-dimensional upper bound limit analysis of underground cavities using nonlinear Baker failure criterion[J]. Transactions of Nonferrous Metals Society of China, 2020, 30(7): 1916-1927. doi: 10.1016/S1003-6326(20)65350-X
[22]	FRALDI M, CAVUOTO R, CUTOLO A, et al. Stability of tunnels according to depth and variability of rock mass parameters[J]. International Journal of Rock Mechanics and Mining Sciences, 2019, 119: 222-229. doi: 10.1016/j.ijrmms.2019.05.001