基于结构极限分析的锚杆加固隧道掌子面稳定性改进评估方法

靳子贤 ,  张成平 ,  王赛旭 ,  孟雪 ,  李巍

科学技术与工程 ›› 2026, Vol. 26 ›› Issue (21) :9260-9271

PDF (10409KB)
科学技术与工程 ›› 2026, Vol. 26 ›› Issue (21) :9260-9271 DOI: 10.12404/j.issn.1671-1815.2505437
交通运输

基于结构极限分析的锚杆加固隧道掌子面稳定性改进评估方法

作者信息 +

Improved Assessment Method for Stability of Tunnel Face Reinforced by Anchor Rods Based on Structural Limit Analysis

Author information +
文章历史 +

摘要

中国战略性基础设施建设正逐步向西部复杂艰险山区推进。在深埋长大山岭隧道建设过程中,当穿越复杂不良地质条件时,为保障施工安全与结构稳定,通常采用主动防控措施对隧道掌子面进行加固。其中,锚杆加固是一种常见且有效的主动防护手段。为了实现对锚杆拉伸状态下极限承载力的快速评估,改进既有基于常规极限分析的锚杆加固隧道掌子面稳定性分析方法,提出了一种将Von-Mises屈服准则与弹性补偿法(elastic compensation method, ECM)相结合的结构极限分析数值方法。详细阐述了在主应力空间中利用Haigh-Westergaard坐标系表达Von-Mises屈服面的坐标变换方法,并进一步说明了在该屈服面上实现弹性补偿法的关键迭代技术。通过Python编程语言对商用有限元软件ABAQUS进行二次开发,实现了基于结构极限分析的锚杆极限承载力评估,为锚杆极限承载力的快速评估提供了一种高效且可行的新途径。在此基础上,引入结构极限分析方法,改进了传统的基于极限分析的隧道掌子面锚杆加固稳定性分析方法。研究结果表明,该改进方法具有高度自动化的特点,迭代计算过程稳定可靠。计算结果进一步显示,锚杆在较弱地层中的加固效果更为显著,且其安全系数对地层强度的变化具有一定的敏感性。

Abstract

Strategic infrastructure construction in China is currently being intensified towards the complex mountainous regions in the western part of the country. During the construction of deep, long mountain tunnels, active reinforcement measures are routinely implemented to ensure construction safety and structural stability when unfavorable geological conditions are encountered. Among these measures, anchor bolt reinforcement is widely regarded as a prevalent and effective protective measure. To enable rapid assessment of the ultimate bearing capacity of anchor bolts under tensile conditions and to improve conventional limit-analysis-based stability assessment methods for tunnel faces reinforced with anchor bolts, a novel numerical approach integrating the Von-Mises yield criterion with the ECM (elastic compensation method) was developed. The coordinate transformation methodology for expressing the Von-Mises yield surface in the principal stress space using Haigh-Westergaard coordinates was systematically described. Critical iterative techniques for implementing ECM on this yield surface were thoroughly elucidated. Secondary development of the commercial finite element software ABAQUS was carried out using Python programming, thereby enabling the automated evaluation of the ultimate capacity of anchor bolts through structural limit analysis. Subsequently, the conventional limit analysis method for tunnel face stability assessment was refined through the incorporation of structural limit analysis. The research results demonstrate that the improved method is highly automated, with a stable and reliable iterative calculation process. It is revealed that the effectiveness of anchor bolt reinforcement is significantly greater in weaker strata. Moreover, safety factors are shown to be notably sensitive to variations in stratum strength parameters.

Graphical abstract

关键词

锚杆 / 隧道 / 极限分析 / 安全系数 / 数值模拟

Key words

anchor bolt / tunnel / limit analysis / safety factor / numerical simulation

引用本文

引用格式 ▾
靳子贤, 张成平, 王赛旭, 孟雪, 李巍. 基于结构极限分析的锚杆加固隧道掌子面稳定性改进评估方法[J]. 科学技术与工程, 2026, 26(21): 9260-9271 DOI:10.12404/j.issn.1671-1815.2505437
JIN Zi-xian, ZHANG Cheng-ping, WANG Sai-xu, MENG Xue, LI Wei. Improved Assessment Method for Stability of Tunnel Face Reinforced by Anchor Rods Based on Structural Limit Analysis[J]. Science Technology and Engineering, 2026, 26(21): 9260-9271 DOI:10.12404/j.issn.1671-1815.2505437

登录浏览全文

4963

注册一个新账户 忘记密码

随着我国经济飞速发展以及“西部大开发”和“一带一路”倡议的深入推进,我国隧道工程建设正逐步向地质条件复杂的西部地区延伸[1]。受复杂地质构造、不良地质体空间分布不确定性及施工经济性等因素的制约,深埋长大山岭隧道往往不可避免地穿越多种不良地质单元。这些特殊地质条件显著增加了隧道施工过程中发生塌方、岩爆等地质灾害的风险,给工程建设带来了严峻挑战。
岩土工程稳定性研究近年来受到广泛关注[2-5],而在隧道工程中,掌子面稳定性是关系到安全施工的重要问题之一。在常规隧道开挖过程中,隧道所处地层强度大,掌子面自稳能力较强,在隧道掌子面上不需要施加支护压力。但在软弱围岩隧道的施工过程中,为了实现对隧道掌子面的加固,保证围岩的稳定性,近年来在隧道施工过程中往往会对隧道掌子面进行超前预加固。超前预加固措施已经在很多工程中取得了良好的效果,通过预加固,掌子面稳定性会得到提高,从而降低了施工风险。
许多研究学者通过数值模拟、模型试验和理论分析等方法对加固隧道掌子面稳定性问题进行了研究。张霄等[6]研究了隧道掌子面预注浆参数的定量化设计方法,建立了考虑预注浆加固作用的开挖面稳定性计算模型,采用极限平衡法以及强度折减法相结合的方式评估了隧道掌子面的稳定性。Pan等[7]基于极限分析法,将隧道掌子面的破坏机制构建为三维牛角形,并进一步将超前锚杆简化为嵌入隧道掌子面中的土钉结构,考虑了锚杆-土体之间相互作用的破坏模式,建立了超前锚杆加固作用下的隧道掌子面稳定性分析理论模型。王明年等[8]基于极限平衡方法,创建了一种由棱柱体与楔形体结合的隧道掌子面破坏机制,并分析了隧道掌子面稳定性。Dias等[9]采用均匀化方法对锚杆加固隧道掌子面进行了三维数值分析,结果表明均匀化方法是水平锚杆加固隧道工作面设计的有效方法。Shin等[10]开发了大比例尺模型试验装置,并据此研究了颗粒土中超前管棚加固条件下的隧道掘进施工过程。王凤云等[11]通过室内相似模型试验,研究了长短锚杆对深埋隧道软弱围岩锚固效果的影响。王惠贤等[12]推导了非均匀应力场中深埋圆形隧洞在锚杆支护条件下的应力与位移弹性解析解。实践证明,在工程地质条件较差的隧道施工中,隧道掌子面加固是提高施工过程中隧道掌子面稳定性的一种经济有效的加固技术。
尽管已有大量学者针对锚杆加固条件下隧道掌子面的稳定性开展了广泛研究,采用了多种研究方法,并提出了诸多具有重要参考价值的掌子面稳定性分析理论模型,但关于锚杆这一关键加固构件的极限承载力研究仍较为缺乏。目前对锚杆极限承载力的取值主要依赖于常规检测或经验取值,通常直接采用规定的标准屈服强度值或其他理论计算值。然而,锚杆的截面形状、长度以及材料属性等均会对锚杆的极限承载力产生影响,进而加大了评估锚杆加固条件下隧道掌子面稳定性的难度。
现提出一种将Von-Mises屈服准则与弹性补偿法(elastic compensation method, ECM)相结合的结构极限分析数值方法。基于Python编程语言对有限元分析软件ABAQUS进行二次开发,提出一种用于锚杆极限承载力快速评估的新方法。并进一步在常规极限分析方法中引入ECM,完成对锚杆加固隧道掌子面稳定性评估方法的改进。该方法具有较高的适用性和稳定性,分析过程实现高度自动化,以期能够更准确地反映锚杆对隧道掌子面的加固效果,并对地层强度的变化表现出一定的敏感性。

1 问题阐述与理论模型

用于隧道掌子面稳定性分析的三维旋转破坏机制是由Mollon等[13]开发的,用于计算圆形断面隧道的安全系数和搜索其破坏面。这种旋转破坏机制已被证明能够为摩擦土中隧道掌子面稳定性分析提供精确上限解[13-14]。Pan等[7]对三维旋转破坏机制进行了改进,提出了适用于非圆形隧道掌子面(如马蹄形隧道)的三维旋转破坏机制,并假设有一个厚度有限的相互作用区来计算只考虑拉伸作用的锚杆失效时的能量耗散。

在实际工程中,隧道开挖时掌子面前方的弯曲效应很小,可以忽略不计[15]。基于这一现象,用于加固隧道的锚杆的失效模式可分为拉伸失效和拔出失效。因此,如何快速且有效地评估锚杆的拉伸极限承载力和拔出极限承载力变得尤为重要。本文研究将结构极限分析方法与马蹄形隧道的三维旋转破坏机制相结合,运用结构极限分析对锚杆的极限承载力进行精确评估,并基于三维旋转破坏机制对隧道掌子面的破坏区域进行搜索及安全系数的计算,旨在提出一种更为准确且高效的锚杆加固隧道掌子面稳定性分析方法。图1为高D、埋深C的隧道的旋转破坏机制纵剖面示意图,A点、B点分别代表隧道掌子面的顶部和仰拱,F点代表破坏机制的顶端。整个破坏机制以角速度ω绕通过O点的水平轴旋转,并假设该破坏机构采用圆柱形旋转速度场,对于破坏机制中的任何一点,速度等于角速度与旋转中心与考虑点之间的距离的乘积。锚杆在两种失效模式下的受力情况也在图1中给出。锚杆拉伸破坏主要关注锚杆自身的极限抗拉强度,而锚杆拔出破坏则侧重于地层与锚杆之间的黏结性能。

2 结构极限分析数值方法

2.1 结构极限分析基础理论

极限分析是一种基于塑性力学理论发展而来的研究材料在极限状态下的分析方法。该方法计算过程简便,物理意义清晰,便于在实际工程中应用,能够为工程师评估地下与隧道工程的稳定性提供有力的理论支持和实践参考。这一方法已经在岩土工程稳定性问题中得到了广泛的应用[7,16-18]。

对于结构极限分析,Hodge等[19]首次提出了针对结构的基于有限元程序的极限分析方法。这一方法已在工业设计上使用多年,涉及在有限元分析过程中调整材料单元的弹性模量,使应力在固定应变分布下处于屈服条件内,然后进一步利用材料单元的弹性模量的新空间分布来求解弹性问题。在每个阶段,可以通过缩放得到的弹性解来找到极限载荷的下界,从而使应力位于当前弹性解的屈服范围内。既有的研究结果表明,结构极限分析的计算过程通常可以得到一个单调递增的下界序列,直到计算得到的结果满足设计约束时可以终止该过程。该方法优于涉及编程技术的直接求解方法[20-22],因为它涉及传统的有限元方法,并且该方法适用于解决具有大量自由度的问题,规避了收敛缓慢的技术缺陷。结构极限分析方法也比简化的分析方法[23]有优势,因为它有更广泛的应用范围。

如图2所示,锚杆中的主体构件为拉杆,用于抵抗隧道掌子面失稳时产生的拉力,保持掌子面稳定性。拉杆的屈服破坏是锚杆失效的重要标志。实验结果表明,对于锚杆这类受拉构件,塑性屈服准则中最接近实际条件的是Von-Mises屈服准则[24]。Ponter等[25]证明了材料在Von-Mises屈服条件下的精确极限载荷计算结果可以用剪切模量具有空间变化的不可压缩线性弹性解精确模拟,这为基于Von-Mises屈服条件的结构极限分析提供了理论基础。

2.2 Von-Mises屈服准则

根据Von-Mises屈服准则给出的屈服条件,材料在主应力空间中的3个主应力(σ1,σ2,σ3)共同作用下达到屈服,这一准则可表示为

$\sqrt{\frac{1}{2}\left[\right({\sigma }_{1}-{\sigma }_{2}{)}^{2}+({\sigma }_{2}-{\sigma }_{3}{)}^{2}+({\sigma }_{3}-{\sigma }_{1}{)}^{2}]}$=σS

式(1)中:σS为材料的屈服强度。

式(1)表示的Von-Mises屈服面在主应力空间中为一个空心的非实体圆柱面,如图3(a)所示。静水应力轴穿过Von-Mises屈服面中心处,屈服面沿着静水应力轴方向在π-平面上的投影是标准圆形。在材料屈服问题的研究过程中,如果只考虑2个主应力,则可以假设σ3的取值为0,屈服面在二维应力空间中的投影可以表示为

${\sigma }_{1}^{2}$-σ1σ2+${\sigma }_{2}^{2}$=${\sigma }_{\mathrm{S}}^{2}$

如图3(b)所示,基于各向同性硬化准则,当材料到达屈服阶段时,塑性应变的增加会进一步导致屈服面增大。本文研究重点关注材料的极限承载力,即材料临近屈服的极限状态,这种塑性硬化现象不在考虑的范围内。

2.3 Haigh-Westergaard坐标系

数值极限分析方法涉及有限元分析技术,为了最大限度地减少有限元分析时所需的输入、输出参数,优化计算过程,采用圆柱坐标系对主应力空间中的Von-Mises屈服面进行应力状态转换,这一过程涉及Haigh-Westergaard坐标系。

在当前步骤中,主应力直接从柯西应力张量σij计算,以便于转换到主应力空间。主应力(σ1,σ2,σ3)均为σij的特征值,表达式为

det(σij-λδij)=0
λ3-I1λ2+J2λ-J3=0

式中:I1、J2和J3为应力不变量;λ为特征值变量。

I1=σ11+σ22+σ33
J2=$\frac{1}{2}$(${I}_{1}^{2}$-σijσij)
J3=det(σij)

为了描述该空间中的应力状态,采用Haigh-Westergaard坐标系(ξ,ρ,θ)。这些坐标直接由应力不变量导出,可以被定义为

ξ=$\frac{{I}_{1}}{\sqrt{3}}$=$\frac{{\sigma }_{1}+{\sigma }_{2}+{\sigma }_{3}}{\sqrt{3}}$
$\rho=\frac{1}{\sqrt{3}} \sqrt{\left(\sigma_{1}-\sigma_{2}\right)^{2}+\left(\sigma_{2}-\sigma_{3}\right)^{2}+\left(\sigma_{3}-\sigma_{1}\right)^{2}}$
θ=$\frac{1}{3}$arccos$\left(\frac{3\sqrt{3}}{2}\frac{{J}_{3}}{{J}_{2}^{3/2}}\right)$

式中:ξ为静水应力的大小;ρ为偏应力的大小;θ定义了主应力轴在偏平面内的方向。

Haigh-Westergaard坐标(ξ,ρ,θ)与主应力坐标之间的变换方程为

$\left[\begin{array}{l}{\sigma }_{1}\\ {\sigma }_{2}\\ {\sigma }_{3}\end{array}\right]$=$\frac{1}{\sqrt{3}}$ξ+$\sqrt{\frac{2}{3}}$ρ$\left[\begin{array}{l}\mathrm{c}\mathrm{o}\mathrm{s}\theta \\ \mathrm{c}\mathrm{o}\mathrm{s}\left(\theta -\frac{2\mathrm{\pi }}{3}\right)\\ \mathrm{c}\mathrm{o}\mathrm{s}\left(\theta +\frac{2\mathrm{\pi }}{3}\right)\end{array}\right]$

Haigh-Westergaard坐标系的优势在于它为主应力空间中的每一种应力状态提供了独特而全面的表示。具体来说,任何应力张量都可以映射到一个唯一的点(ξ,ρ,θ),这确保了整个应力空间被覆盖,没有冗余或间隙。

2.4 弹性补偿法

采用被称为弹性补偿法(ECM)的极限分析下限方法实现对锚杆极限承载力的预测,这是一种极限分析的静态方法,需要建立适合分析的容许应力场。ECM通过对结构进行一系列基于线性弹性有限元分析的迭代计算来逐步逼近结构的极限状态。在这一过程中,结构的高负载区域(应力大于材料屈服值的区域)的强度会被削弱(即弹性模量会被降低),结构所受荷载的值也会更新。这实现了对结构内部应力的重新分配,使整个结构在给定荷载条件下的最大允许应力值可以被搜索。

数值计算过程的基本原理比较简单,即在给定的荷载条件下,计算结构的线性弹性解,提取所有单元中各高斯点(GP)处的应力值,并在主应力空间中采用Haigh-Westergaard坐标系表示。为了简化计算,#ith单元的弹性解的应力值被定义为该单元上所有高斯点处主应力的平均值,在主应力空间中应力点被标记为${{P}^{e}}_{\#i}$(ξe,ρe,θe),如图4所示。在${{\overrightarrow{OP}}^{\mathrm{e}}}_{\#\mathrm{i}}$/$\left|{{\overrightarrow{OP}}^{\mathrm{e}}}_{\#i}\right|$方向上的射线与Von-Mises屈服面的交点即为屈服点,表示为${{P}^{\mathrm{Y}}}_{\#i}$。在进行分析时,首先需要设置材料参数[E(k),v]及荷载${P}_{\mathrm{D}}^{\left(s\right)}{\stackrel{-}{p}}_{\mathrm{i}}$,其中E(k)为材料的弹性模量,v为材料的泊松比,${\stackrel{-}{p}}_{i}$为指定的参考荷载。进行第一次线性弹性有限元分析时,k和s的取值为1,即E(k)=E(1),${P}_{\mathrm{D}}^{\left(s\right)}$=${P}_{\mathrm{D}}^{\left(1\right)}$。当#ith单元的线性弹性有限元分析得到的应力满足${\left|{{\overrightarrow{OP}}^{\mathrm{e}}}_{\#\mathrm{i}}\right|}^{\left(\mathrm{k}\right)}$>${\left|{{\overrightarrow{OP}}^{\mathrm{Y}}}_{\#i}\right|}^{\left(k\right)}$,单元的弹性模量将根据式(12)进行缩减。

${E}_{\#i}^{(k+1)}$=${E}_{\#i}^{\left(k\right)}\frac{{\left|{{\overrightarrow{OP}}^{\mathrm{Y}}}_{\#\mathrm{i}}\right|}^{\left(\mathrm{k}\right)}}{{\left|{{\overrightarrow{OP}}^{\mathrm{e}}}_{\#\mathrm{i}}\right|}^{\left(\mathrm{k}\right)}}$

值得注意的是,在迭代过程中,v的值是不会变化的。在进行第k次分析时,需要重点关注所有单元的线性弹性解中的最大值,即最大应力值,定义为${P}_{\mathrm{R}}^{\left(k\right)}$。如图5所示,在主应力空间中,为所有应力解中距离Von-Mises屈服面最远的点,与其相对应的屈服点为${P}_{\mathrm{R}}^{\mathrm{Y}\left(k\right)}$。当${\left|{\overrightarrow{OP}}_{\mathrm{R}}\right|}^{\left(\mathrm{k}\right)}$大于或者等于${\left|{\overrightarrow{OP}}_{\mathrm{R}}^{\mathrm{Y}}\right|}^{\left(\mathrm{k}\right)}$时,开始进入新一轮线性弹性有限元分析,保持荷载${P}_{\mathrm{D}}^{\left(s\right)}{\stackrel{-}{p}}_{i}$不变,采用式(12)得到所有单元的更新的弹性模量${E}_{\#i}^{(k+1)}$。直到最大主应力PR小于或者等于与其对应的屈服值${P}_{\mathrm{R}}^{\mathrm{Y}}$时,弹性模量停止更新。以上对弹性模量${E}_{\#i}^{\left(k\right)}$进行更新的过程被包含在对荷载${P}_{\mathrm{D}}^{\left(s\right)}{\stackrel{-}{p}}_{i}$进行更新的迭代过程中,即每次迭代过程中都增加${P}_{\mathrm{D}}^{\left(s\right)}{\stackrel{-}{p}}_{i}$的值,并重复对弹性模量的更新过程,直到进一步的载荷增加不允许最大应力PR的值低于与其对应的屈服应力${P}_{\mathrm{R}}^{\mathrm{Y}}$的值。最终,假设迭代过程结束时k=K和s=S,极限荷载乘数的下界PLB可以被计算,表达式为

PLB=${\left|{\overrightarrow{OP}}_{\mathrm{R}}^{\mathrm{Y}}\right|}^{\left(K\right)}\frac{{P}_{\mathrm{D}}^{\left(S\right)}}{{\left|{\overrightarrow{OP}}_{\mathrm{R}}\right|}^{\left(K\right)}}$

式(13)中:PLB对应于每次迭代中最高载荷时的荷载乘数,使得弹性模量在更新过程中的最终补偿应力${P}_{\mathrm{R}}^{\left(K\right)}$正好达到屈服。准确地说,是在极限荷载PLB${\stackrel{-}{p}}_{i}$下,结构网格单元的最大应力恰好达到屈服值${P}_{\mathrm{R}}^{\mathrm{Y}}$(至少在有限元网格的一个区域),所有其他应力值不超过屈服应力。迭代过程中,不考虑材料的塑性硬化,否则Von-Mises屈服面将会扩展,导致ECM方法难以接近极限分析下限解。

3 基于ABAQUS的结构极限分析

3.1 结构极限分析流程

极限分析的数值方法基于有限元分析,所有迭代过程中的线性弹性分析都使用商用有限元软件ABAQUS进行,弹性模量和荷载的更新以及极限荷载乘数的下界的计算使用Python实现。为了更清楚地解释ECM的基本步骤与迭代方案,绘制了ECM在ABAQUS中的流程图,如图6所示。

3.2 结构极限分析实施

Pan等[7]记录了一个采用水平锚杆加固隧道掌子面的实际案例,即Toulon隧道。Toulon隧道是一个实际存在的马蹄形隧道,属于非圆形隧道,其隧道截面几何形状与圆形隧道之间存在明显差异,采用传统的隧道开挖方法施工。

Toulon隧道截面如图7所示。为凸显马蹄形隧道截面几何形状与圆形隧道的差异,明确界定其边界轮廓特征,并体现锚杆加固范围,绘制了与马蹄形隧道截面面积相等的虚线圆作为参考线。为了防止施工过程中隧道掌子面出现坍塌等问题,掌子面的锚杆加固密度达到了0.7 根/m2。加固隧道掌子面的锚杆参数如表1所示。

在实际施工过程中,观察到隧道开挖过程中掌子面前方的弯曲效应很小,可以忽略不计。根据Pan等[7]的研究,在构建数值模型时,假设锚杆为一维结构单元,仅考虑其轴向抗拉性能。为了模拟锚杆真实受力情况,得到更可靠的锚杆极限承载力评估结果,采用实体单元构建锚杆数值模型。基于实体单元的数值模型可以精细化模拟不同类型锚杆,而不是将锚杆假设为一个仅受轴向力的简单的一维杆件。ECM是在每个单元的每个高斯点上进行的,这使得结构极限分析的计算量非常庞大,对计算机的运算能力提出了很高的要求。实际应用中,为提高计算速度,降低对计算机运算能力的要求,采用的ABAQUS内置的缩减积分网格单元对数值模型进行离散。相比于完全积分,缩减积分网格单元在每个方向上都减少了一个积分点,一阶缩减积分单元只有一个积分点在网格单元的质心处。这极大地减少了ECM的计算量,也使ECM在大型精细化数值模型上的使用成为可能。值得注意的是,由于缩减积分单元引入了“人工沙漏刚度”,该方法在网格细化时效果更佳,因此适当细化网格能解决计算过程中额外的沙漏问题(即发生变形但无应变的现象)。本文研究使用的单元类型均通过了形状、尺寸和分析校核标准。经过初步的网格敏感性研究后,选择有限元的数量,以确保精确的有限元弹性解。建立的锚杆数值模型及网格如图8所示,数值模型被离散为23 564个三维实体网格单元。

为了模拟锚杆承受轴向拉力时的受力情况,在计算拉伸极限承载力时,在杆件两端网格单元节点上施加轴向荷载,荷载为${P}_{\mathrm{D}}^{\left(s\right)}{\stackrel{-}{p}}_{i}$。值得注意的是,${P}_{\mathrm{D}}^{\left(1\right)}$被设置为1,参考荷载${\stackrel{-}{p}}_{i}$决定了荷载初始值的大小,如果荷载初始值过大,那么结构在第一次进行线性有限元计算时就会被破坏,导致ECM的迭代无法进行,本文研究中,参考荷载${\stackrel{-}{p}}_{i}$取10 kN。在计算锚杆拔出极限承载力时,在ABAQUS中为锚杆结构表面施加摩擦力,这一过程通过在接触属性中定义切向行为来实现。

3.3 结构极限分析计算结果与讨论

将锚杆数值模型通过ECM得到的拉伸极限承载力Pmax(PLB${\stackrel{-}{p}}_{i}$)与Pan等[7]记录的锚杆拉伸屈服强度Tm对比,如图9(a)所示。二者的计算结果比较接近,但是与ECM方法得到峰值荷载结果相比,Pan等[7]给出的锚杆拉伸屈服强度相对保守,这可能是因为该强度值作为结构设计参数,在实际工程中出于安全考虑通常取值偏保守,以确保施工过程的可靠性。锚杆的拔出力Tp取决于锚杆在土体中的有效长度等因素,其理论计算公式为

Tp=πτmΦmin(Le,L-Le)

式(14)中:Tp为锚杆发生拔出失效时的拔出力;τm为锚杆与土之间接触面极限黏结强度;L为锚杆的总长度;Le为锚杆的有效长度;Φ为锚杆的直径。选取了有效长度Le为6 m的锚杆,采用ECM计算其被拔出时的极限承载力,将计算结果与Tp的理论计算公式给出的结果进行对比,如图9(b)所示。与理论方法相比,ECM方法得到的锚杆拔出失效时的极限承载力大于由式(14)得到的理论值,这是由于ECM方法在计算拔出极限承载力时,以接触面的方式更好地模拟了土与锚杆之间的相互作用,并且在拔出过程中考虑了锚杆的受力及变形情况,而非将锚杆简单假设为刚性杆件。

ECM得到的荷载乘数${P}_{\mathrm{D}}^{\left(s\right)}$随着迭代次数s的增加而逐渐增大,并接近锚杆拉伸失效和拔出失效时的极限承载力,这符合极限分析下限方法的规律。事实上,在ABAQUS中采用ECM得到的结构极限承载力的值可以精确到1×10-12,但是在实际应用中,过高的精确度是没有必要的,因此本研究中计算结果精确度取1×10-3。在计算迭代初期,荷载乘数PLB增长速度很快,此后开始减缓;在接近结构极限状态前,荷载乘数增长速度变得非常缓慢。这是由于荷载初始值设置较小,当施加在结构上的荷载与峰值荷载之间的差距较大时,弹性模量的变化程度较大。

ECM在对结构极限承载力预测的同时也给出了锚杆结构中可能发生破坏的薄弱区域。以计算Toulon隧道中锚杆发生拉伸破坏时的极限承载力Pmax为例,本文研究的破坏机制通过ECM最后一次迭代产生的应力云图表达,如图10所示。ECM方法的最终结果显示了较大应力出现的位置,反映出了结构的薄弱区域。值得注意的是,前文提到,ECM方法是通过削弱结构的高负载区域,对结构的应力进行重分配来实现对结构的极限承载力真实值的逼近。由于在迭代过程中大多数网格单元的弹性模量都被更新,这导致被分析结构的材料属性被改变,成为一种现实中并不存在的“虚拟结构”。这一过程使计算得到的应力云图可以对结构的破坏机制进行定性分析,而云图中计算结果的数值没有指导意义。不过,ECM方法已经为结构的极限承载力提供了非常可靠的数值解,这使云图中计算结果的数值变得并不关键。

4 改进的稳定性评估方法

4.1 极限分析方法改进

基于Pan等[7]提出的水平锚杆加固隧道工作面安全系数评价方法,将马蹄形隧道掌子面的破坏区域假设为“三维牛角”,并假设锚杆穿过破坏区域,如图11所示。在采用极限分析方法评估隧道掌子面安全性时,将锚杆的拉伸破坏和拉出破坏分别考虑,并假设了一个相互作用区来表达锚杆与土体之间的相互作用[7]。锚杆发生拉伸破坏及拔出破坏时产生的总抵抗力的做功速率可以分别表示为

${W}_{\mathrm{b}}^{\mathrm{t}\mathrm{f}}$=ω$\stackrel{m}{\sum _{1}}$TmRkcosβk$\left(\frac{1}{2}+\frac{2\alpha \delta }{{\mathrm{e}}^{2\alpha \delta }-{\mathrm{e}}^{-2\alpha \delta }}\right)$
${W}_{\mathrm{b}}^{\mathrm{p}\mathrm{f}}$=ω$\stackrel{m}{\sum _{1}}$TpRkcosβk$\left[1+\frac{\delta }{2\mathrm{m}\mathrm{i}\mathrm{n}({L}_{\mathrm{e}},L-{L}_{\mathrm{e}})}\right]$

式中:ω为破坏区域的旋转角速度;Tm为锚杆拉伸屈服强度;Rk和βk为锚杆与破坏区域交点的极坐标;α为锚杆与土体相互作用区相关的参数;δ为相互作用区的宽度;m为锚杆的数量。更详细的做功速率计算方法可参见Pan等[7]的研究。

基于本文研究提供的锚杆极限承载力评估方法,式(15)与式(16)中的Tm与Tp将由采用ECM计算得到的Pmax替代。本文研究采用Python编程语言实现三维旋转破坏机制对隧道掌子面稳定性的评估,确保了结构极限分析与岩土极限分析之间的有效衔接,实现了计算过程的自动化,并保障了锚杆加固条件下隧道掌子面稳定性评估方法的运算稳定性。

4.2 改进计算结果对比与验证

本文研究以Toulon隧道为依托,将分别采用数值模拟[7]、常规极限分析[7]和ECM改进极限分析得到的隧道掌子面安全系数(FS)进行对比。对比结果如图12所示。ECM改进极限分析得到的安全系数与常规极限分析得到的安全系数(FS)非常接近。在不同围岩摩擦角(ϕ)条件下,采用ECM改进极限分析得到的安全系数均大于常规极限分析得到的安全系数,这是由于结构极限分析给出的锚杆拉伸失效与拔出失效时的极限承载力均大于常规极限分析[7]中Tm与Tp的取值。值得注意的是,采用ECM改进的极限分析在较小的摩擦角(ϕ)条件下,计算出的安全系数与常规极限分析得到的安全系数之间的差异较大,随着摩擦角(ϕ)增大,二者之间的差异在逐渐减小。这说明在地层强度较弱的情况下,锚杆对隧道掌子面的加固作用更加明显,此时锚杆承载力的提高对隧道掌子面安全系数有不可忽视的影响。虽然常规极限分析得到的安全系数更加保守,可以给实际工程设计提供更具安全性的理论依据,但快速得到更精确的锚杆极限承载力评估结果,进而预测锚杆加固条件下隧道掌子面稳定性同样具有研究意义。该技术的发展有助于研究人员更深入地理解不同类型锚杆的力学特性及其对隧道掌子面的加固效果,具有重要的实际意义。

Pan等[7]对采用水平锚杆加固的圆形截面隧道进行了稳定性分析,该圆形截面具有与Toulon隧道马蹄形截面相当的等效面积,具体形式可参考图7。采用常规极限分析[7]和ECM改进极限分析得到的圆形隧道掌子面安全系数(FS)进行对比,结果如图13所示。与马蹄形隧道类似,ECM改进极限分析所得的安全系数略高于常规极限分析结果。在摩擦角(ϕ)较小的条件下,两种方法计算得到的安全系数差异较大;随着摩擦角(ϕ)增大,该差异逐渐减小。

为进一步验证改进方法的科学性和有效性,本研究参考并对比了Zhang等[26]的研究成果。该研究系统分析了锚杆加固对隧道掌子面稳定性的影响,并提出了适用于实际工程的简化理论模型(以下简称简化理论模型)。此外,该研究对比了理论模型与数值模拟计算所得的安全系数(FS),相关隧道参数及锚杆参数如表2所示。

依据表2中的锚杆参数,在ABAQUS中建立了锚杆数值模型,用于评估锚杆的拉伸极限承载力与拔出极限承载力,并进一步采用引入ECM改进的极限分析方法对隧道掌子面稳定性进行评估。图14给出了采用数值模型[27]、简化理论模型[26]和ECM改进极限分析所得锚杆加固条件下隧道掌子面安全系数(FS)的对比结果。

图14(a)显示,随着黏聚力c的增大,不同方法所得的安全系数均呈线性增长趋势。采用ECM改进极限分析得到的安全系数高于简化理论模型[26]和数值模型[27]的结果,且与数值模型结果更为接近。此外,随着c的增大,ECM改进极限分析与数值模型及简化理论模型之间的安全系数差异逐渐减小。图14(b)显示,随着摩擦角ϕ的增大,不同方法所得的安全系数均呈现非线性增长趋势。ECM改进极限分析与数值模型所得的安全系数十分接近,且二者均高于简化理论模型的结果。对比中所涉及计算方法在结果上产生差异的原因在于其基本原理的不同。简化理论模型方法的运算基于静力平衡条件,且对地层滑移破坏模式的假设较为理想化,同时在确定锚杆极限承载力时,未充分考虑锚杆截面形状、长度及材料属性等关键因素的影响。相比之下,数值模型方法能够更真实地模拟锚杆加固条件下隧道掌子面临界破坏状态,从而提供更为合理的安全系数评估结果。然而,该方法未能计入锚杆失效过程中的能量耗散机制,这正是ECM改进极限分析法与数值模型所得安全系数存在差异的主要原因。总体而言,ECM改进极限分析方法能够提供合理的安全系数评估结果,在地层抗剪强度参数变化条件下,其安全系数的变化趋势与既有研究成果高度一致。该方法更全面地考虑了锚杆对隧道掌子面稳定性的影响,所得评估结果具有较高的实际工程应用价值。

4.3 改进计算结果分析与讨论

图15展示了锚杆长度比(L/D)为1.2,锚杆密度为0.7 根/m2时,不同黏聚力(c)(15、20、25、30 kPa)条件下,摩擦角(ϕ)与安全系数(FS)的关系曲线。不同黏聚力曲线基本呈现出线性增长趋势,即随着摩擦角的增大,安全系数也相应增加。同时,随着黏聚力的增加,安全系数的基线水平也相应提升。围岩自身强度增加,隧道开挖面安全系数随之增加,即隧道的稳定性逐渐增强,这与以往的研究结论一致。同时,在隧道开挖面锚杆密度为0.7 根/m2,锚杆长度比为1.2时,图中所示不同围岩参数下的隧道开挖面安全系数均大于1,即隧道面处于稳定状态。

图16展示了锚杆长度比(L/D)为1.2,围岩黏聚力为20 kPa时,不同围岩摩擦角(15°、20°、25°、30°、35°)条件下,锚杆密度与安全系数(FS)的关系曲线。不同围岩摩擦角曲线基本呈现出线性增长趋势,即随着锚杆密度的增大,安全系数也相应增加。同时,随着围岩摩擦角的增加,安全系数的基线水平也相应提升。围岩摩擦角增加,隧道开挖面安全系数随之增加,且增加趋势逐渐变缓。同时,在隧道开挖面锚杆密度为0根/m2时,图中存在隧道开挖面安全系数小于1,即隧道面处于不稳定状态的情况,随着锚杆密度的增加,隧道开挖面安全系数均变为大于1,这体现了在实际隧道工程中开挖面锚杆加固措施的重要性。

图17展示了围岩黏聚力为20 kPa,锚杆密度为0.7 根/m2时,不同围岩摩擦角(15°、20°、25°、30°、35°)条件下,锚杆长度比(L/D)与安全系数(FS)的关系曲线。当锚杆长度比小于1时,不同围岩摩擦角曲线基本呈现出线性增长趋势,即随着锚杆长度比的增大,安全系数也相应增加。当锚杆长度比大于1时,安全系数基本不再随锚杆长度比变化,这是因为锚杆长度超过隧道前方破坏区土体厚度,锚杆超出长度部分无法对破坏区土体产生有效作用。同时,在隧道开挖面锚杆密度为0.7 根/m2,c=20 kPa时,图中所示不同围岩参数下的隧道开挖面安全系数均大于1,即隧道面均处于稳定状态。

5 结论与建议

首次将基于弹性补偿法的结构极限分析方法应用于锚杆极限承载力的评估,为解决传统依赖常规检测或经验取值来确定锚杆拉伸极限承载力与拔出极限承载力的问题提供了新途径,并在此基础上改进了用于锚杆加固隧道掌子面稳定性评估的常规极限分析方法,获得了对地层强度更具敏感性且具有重要研究意义的计算结果。

(1)基于既有研究,明确了Von-Mises屈服准则模拟锚杆进行数值结构极限分析的可行性,并在主应力空间中采用Haigh-Westergaard坐标系对Von-Mises屈服面进行了坐标转换。

(2)提出了基于Von-Mises屈服准则的弹性补偿法(ECM),给出了该方法所涉及迭代过程的详细说明,并进一步阐述了该方法的适用性与局限性等注意事项。

(3)采用编程技术常用的Python语言对商用有限元分析软件ABAQUS进行了二次开发,实现了ECM对锚杆极限承载力的预测,与既有锚杆极限承载力进行了对比,分析了ECM的迭代过程与计算结果。

(4)引入结构极限分析以改进基于常规极限分析的锚杆加固隧道掌子面稳定性分析方法,采用改进方法计算发现锚杆加固效果在较弱地层中更明显,且安全系数对地层强度变化具有敏感性。

在后续研究中,可进一步探讨结构极限分析技术,并研究基于Von-Mises屈服准则与弹性补偿法(ECM)的结构极限分析静力学方法是否适用于大型地下工程结构,从而不断完善隧道掌子面稳定性分析体系,为实际工程提供更为坚实的理论支持,提升施工过程的安全性。

参考文献

原文顺序 | 出版日期 | 本文引用
[1]

李术才, 王康, 李利平, 等. 岩溶隧道突水灾害形成机理及发展趋势[J]. 力学学报, 2017, 49(1): 22-30.
[2]

Li Shucai, Wang Kang, Li Liping, et al. Mechanical mechanism and development trend of water-inrush disasters in karst tunnels[J]. Chinese Journal of Theoretical and Applied Mechanics, 2017, 49 (1): 22-30.
[3]

Jin Z X, Zhang C P, Li W, et al. Stability analysis for excavation in frictional soils based on upper bound method[J]. Computers and Geotechnics, 2024, 165: 105916.
[4]

Tu S Q, Li W, Zhang C P, et al. Seepage effect on progressive failure of shield tunnel face in granular soils by coupled continuum-discrete method[J]. Computers and Geotechnics, 2024, 166: 106009.
[5]

Jin Z X, Zhang C P, Li W, et al. Probabilistic analysis for excavation stability in spatially variable soil based on upper bound method[J]. Computers and Geotechnics, 2024, 176: 106769.
[6]

Li W, Zhang C P. Face stability analysis for a shield tunnel in anisotropic sands[J]. International Journal of Geomechanics, 2020, 20(5): 04020043.
[7]

张霄, 王明年, 于丽, 等. 隧道开挖面预注浆加固稳定性理论计算模型[J]. 现代隧道技术, 2024, 61(5): 42-51.
[8]

Zhang Xiao, Wang Mingnian, Yu Li, et al. Theoretical calculation model for stability of tunnel excavation face with pre-grouting reinforcement[J]. Modern Tunnelling Technology, 2024, 61(5): 42-51.
[9]

Pan Q J, Dias D. Safety factor assessment of a tunnel face reinforced by horizontal dowels[J]. Engineering Structures, 2017, 142: 56-66.
[10]

王明年, 张霄, 赵思光, 等. 软弱围岩隧道机械化全断面施工超前支护体系设计方法研究[J]. 铁道学报, 2020, 42(8): 146-154.
[11]

Wang Mingnian, Zhang Xiao, Zhao Siguang, et al. Research on design method of advance support system for mechanized full-face excavation of tunnels in soft and weak surrounding Rock[J]. Journal of the China Railway Society, 2020, 42(8): 146-154.
[12]

Dias D, Kastner R, Jassionnesse C. Sols renforcés par boulonnage-Etude numérique et application au front de taille d’un tunnel profond[J]. Geotechnique, 2002, GE52(1): 15-27.
[13]

Shin J H, Choi Y K, Kwon O Y, et al. Model testing for pipe-reinforced tunnel heading in a granular soil[J]. Tunnelling and Underground Space Technology, 2008, 23(3): 241-250.
[14]

王凤云, 潘锐, 蔡毅, 等. 深埋隧道长短锚杆锚固围岩破坏模型试验[J]. 科学技术与工程, 2025, 25(6): 2480-2489.
[15]

Wang Fengyun, Pan Rui, Cai Yi, et al. Experimental of failure process of long and short rock bolts for a week deep tunnel[J]. Science Technology and Engineering, 2025, 25(6): 2480-2489.
[16]

王惠贤, 吴庆良, 鲍安红, 等. 锚杆支护下两向不等压深埋隧洞受力与变形弹性分析[J]. 科学技术与工程, 2024, 24(2): 757-765.
[17]

Wang Huixian, Wu Qingliang, Bao Anhong, et al. Elastic analysis of stress and deformation of deep buried tunnel supported by bolts under two-direction unequal pressure[J]. Science Technology and Engineering, 2024, 24(2): 757-765.
[18]

Mollon G, Dias D, Soubra A H. Rotational failure mechanisms for the face stability analysis of tunnels driven by a pressurized shield[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 2011, 35(12): 1363-1388.
[19]

Ibrahim E, Soubra A H, Mollon G, et al. Three-dimensional face stability analysis of pressurized tunnels driven in a multilayered purely frictional medium[J]. Tunnelling and Underground Space Technology, 2015, 49: 18-34.
[20]

Kamata H, Mashimo H. Centrifuge model test of tunnel face reinforcement by bolting[J]. Tunnelling and Underground Space Technology, 2003, 18(2): 205-212.
[21]

Li W, Zhang C P, Zhang D L, et al. Face stability of shield tunnels considering a kinematically admissible velocity field of soil arching[J]. Journal of Rock Mechanics and Geotechnical Engineering, 2022, 14: 505-526.
[22]

Li W, Zhang C P, Tu S Q, et al. Face stability analysis of a shield tunnel excavated along inclined strata[J]. Underground Space, 2023, 13: 183-204.
[23]

Wang S X, Zhang C P, Li W, et al. Partial collapse mechanism of a horseshoe-shaped tunnel face in layered soils[J]. Computers and Geotechnics, 2024, 167: 106114.
[24]

Hodge P G, Belytschko T. Numerical methods for the limit analysis of plates[J]. Journal of Applied Mechanics, 1968, 35: 796-802.
[25]

Karadeniz S, Ponter A R S. A linear programming upper bound approach to the shakedown limit of thin shells subjected to variable thermal loading[J]. The Journal of Strain Analysis for Engineering Design, 1984, 19(4): 221-230.
[26]

Nguyen D H. Direct limit analysis via rigid-plastic finite elements[J]. Computer Methods in Applied Mechanics and Engineering, 1976, 8(1): 81-116.
[27]

Ponter A R S, Karadenix S. A linear programming upper bound approach to the shakedown limit of thin shells subjected to variable thermal loading[J]. Journal of Strain Analysis for Engineering Design, 1984, 19: 221-230.
[28]

Zarka J, Engel J J, Ingebert G. On a simplified inelastic analysis of structures[J]. Nuclear Science and Engineering, 1980, 57(2): 333-368.
[29]

Khan A S, Huang S. Continuum theory of plasticity[M]. New York: John Wiley & Sons, 1995.
[30]

Ponter A R S, Carter K F. Limit state solutions, based upon linear elastic solutions with spatially varying elastic modulus[J]. Computer Methods in Applied Mechanics and Engineering, 1997, 140: 237-258.
[31]

Zhang X, Yu L, Wang M N, et al. Investigation of the failure mechanism and theoretical model of bolt-reinforced shallow tunnel faces with different bolt lengths[J]. Underground Space, 2024, 16: 126-142.
[32]

Zhang X, Wang M N, Wang Z L, et al. A limit equilibrium model for the reinforced face stability analysis of a shallow tunnel in cohesive-frictional soils[J]. Tunnelling and Underground Space Technology, 2020, 105: 103562.

基金资助

国家自然科学基金面上项目(52278388)

国家自然科学基金青年科学基金(52208383)

中央高校基本科研业务费专项资金(2024QYBS002)

北京交通大学人才基金(2025JBRC006)

PDF (10409KB)

54

访问

0

被引

详细
导航
相关文章

/

〈 〉

京ICP备05035734号-4
版权所有 © 《科学技术与工程》编辑部
地址： 北京市海淀区学院南路86号 《科学技术与工程》编辑部　邮编：100081
电话：(010)62118920 　E-mail：stae@vip.163.com
技术支持：北京玛格泰克科技发展有限公司 　
网站管理员登录