强度折减有限元法中锚杆计算模型研究

文章编号：1000－7598 (2010) 12－4021－06

强度折减有限元法中锚杆计算模型研究

刘文连1, 2，韦立德3

（1.昆明理工大学 建筑工程学院，昆明 650224；2.中国有色金属工业昆明勘察设计研究院，昆明 650051

3.中科院武汉岩土力学研究所 岩土力学与工程国家重点实验室，武汉 430071）

摘 要：针对目前岩土工程中锚杆数值模型存在的局限性，受商用软件中锚杆强度模型的启发而提出了一种适用于强度折减有限元法的锚杆强度模型，并建议了考虑锚杆强度折减的强度折减方法，研制了相应程序。该模型应用于岩坡锚固计算，计算结果说明了模型的有效性；应用于土质边坡稳定性评估计算，结果表明模型方便可行。 关 键 词：边坡；锚杆；计算模型；有限元模型；强度折减 中图分类号：TU 470 文献标识码：A

Study of calculation model of anchors in strength reduction FEM

LIU Wen-lian, WEI Li-de3

(1. College of Architectural Engineering, Kunming University of Science and Technology, Kunming 650224, China;

2. Kunming Prospecting Design Institute of China Nonferrous Metals Industry, Kunming 650051, China; 3. State Key Laboratory of Geomechanics and

Geotechnical Engineering, Institute of Rock and Soil Mechanics, Chinese Academy of Sciences, Wuhan 430071, China)

1, 2

Abstract: Due to the limitation of the numerical models of anchors by finite element method (FEM) in geotechnical engineering at present, an anchor model which is fit for the strength reduction FEM is put forward. Strength reduction FEM considering the anchor model is proposed. The corresponding program is also developed; and the proposed model and method are tested by a case study calculation results of a simple rock slope reinforced with anchors. At last the model is used in evaluation of a soil slope stability; and the results of this example indicate that the model is convenient and feasible.

Key words: slope; anchor; calculation model; finite elements model; strength reduction

1 引 言

由于计算结果和二维极限平衡法计算结果有很大的一致性，强度折减有限元法得到了承认和采用。在三维边坡稳定分析方面，三维极限平衡法还没有成熟的理论成果，而三维强度折减有限元法是目前三维边坡稳定分析比较理想的方法，因此发展和完善强度折减有限元法具有重要的理论价值和工程应用意义。锚杆加固是边坡主要加固方法之一，但当前强度折减有限元法中锚杆计算模型不能体现边坡失稳过程锚杆的作用特征，因此当前强度折减有限元算法不适用于需要考虑锚杆锚固作用的边坡稳定性评估问题。

本文讨论了当前强度折减有限元法中锚杆计算模型存在的问题，探索适用的强度折减有限元法中

锚杆计算模型，研制了相应程序，并应用于岩体边坡和土质边坡稳定性评估计算中。

2 锚杆计算模型现状

2.1 不考虑锚杆塑性屈服

用于边坡支护的锚杆失效形式有卸剥、拔出、拉断等形式，其中卸剥和拔出为主要形式，而卸剥（或拔出）对应的锚杆屈服条件为穿过滑动面的锚杆单元轴力大于等于卸剥（或拔出）段沿壁摩阻力之和。在采用强度折减有限元法计算边坡稳定性中，由于不能预知边坡滑动面的准确位置，计算前卸剥（或拔出）段不能确定，从而卸剥（或拔出）段沿壁摩阻力之和不能确定，因此在采用强度折减有限元法计算边坡稳定性中不考虑锚杆失效，而按照弹性模型计算锚杆的作用，这样计算相当于假设锚杆

收稿日期：2009-12-25

第一作者简介：刘文连，男，19年生，博士研究生，教授级高级工程师，主要从事岩土工程、工程勘察等方面勘察、设计和科研工作。 E-mail: LWENL702@sina.com

4022 岩 土 力 学 2010年

不失效。文献[2－4]采用了这样的处理过程。 2.2 考虑锚杆塑性屈服但不折减锚杆强度参数

这种计算模型，对于全程黏结锚杆考虑锚杆塑性屈服，强度参数采用的是锚杆材料抗拉强度参数，

对应锚杆失效形式为拉断形式；对于预应力锚杆考 虑锚杆塑性屈服，强度参数采用的是锚杆锚固段承

载力，对应锚杆失效形式为拔出形式，且滑动面通过自由段。其计算特点为在采用强度折减有限元法分析边坡稳定性过程中，只对岩土体强度参数进行折减，不对锚杆强度参数进行折减。文献[5]采用了这样的处理过程。

3 锚杆计算模型探索

在这里主要回顾成熟的和得到承认的二维极限平衡法中的锚杆计算模型，通过类比建立适合于强度折减有限元法的锚杆计算模型。 3.1 二维极限平衡法中的锚杆计算模型

用于边坡稳定分析的二维极限平衡法商用软件Slide、Geo-slope等软件是普遍得到承认的商用软件，滑坡时锚杆失效形式有卸剥、拔出、拉断等形式。对于均匀岩土体中全程黏结锚杆（如图1所示），在商用软件Slide中锚杆承载力F计算过程如下：

图1 边坡中的锚杆 Fig.1 Anchor in slope

F=min(F1,F2,F3)

⎫F⎪

1=BL0/S,拔出⎪

F=T/S,拉断⎬ （1）

2⎪F3=(P+BLi)/S,卸剥 ⎪⎭

式中：Li为滑坡体内锚杆段长度；L0为滑坡体以外锚杆段长度；B为1 m长的锚杆所具有的沿壁摩阻力；S为锚杆间距；T为锚杆材料抗拉承载力；P为锚杆外端托盘抗拉承载力。按照此式计算，包含三段的典型的锚杆承载力曲线如图2所示。在实际情况中，一般锚杆材料强度比较高而不出现拉断，锚杆承载力曲线对应的锚杆材料抗拉强度部分长度为0；如果锚杆外段不安装托盘等固定装置，锚杆承载力曲线托盘承载力为

0；对不均匀岩土体中全程黏

结锚杆承载力计算过程按此类推。

图2 典型锚杆承载力

Fig.2 A typical anchor force diagram

圆弧滑面是一般土体边坡常见的破坏型式，而瑞典条分法又是圆弧滑动模式下进行边坡稳定性分析的最简单方法。根据Bishop安全系数定义，可将瑞典法中的边坡安全系数写为

fs=

∑cijlij+∑⎡⎣(wijcosαij+Fksinβk)tanϕij+Fkcosβk⎤⎦

∑wijsinαij

（2）式中：fs为滑体安全系数；i、j分别为滑体的大分块数循环变量和大分块中的小土条数循环变量；k为锚杆层数循环变量；cij、ϕij分别为土的黏结力和摩擦系数；lij为滑弧上的小土条底长；wij为小土条浮重度，无水时为湿重度；αij为滑面切线倾角；βk为锚杆轴线与滑面交点切线间夹角；Fk为锚杆在滑面处的轴向拉力。而对于岩体边坡平面滑动可认为

是式（2）的一个特例，即只有1个大分块，且大分块中只有1个小分块。

3.2 适用于边坡稳定性分析的锚杆强度模型

在这里，滑坡时锚杆失效形式考虑有卸剥、拔出、拉断等形式。首先定义卸剥、拔出情况锚杆强度模型如下。为了考虑锚杆的杆壁和岩土体之间剪切屈服性质在每个锚杆单元的两个节点都分别设有一个屈服应力，它由该节点到锚杆两端能够承受沿壁摩阻力的最小值来确定。锚杆单元i节点屈服应

力σi

bp可表示为[6]

min⎛⎜∑2πrkkmm⎞bcFL,σi

bp

=⎝k∑2πrmbcFL⎟

⎠πr2

（3）b

式中：ck

F为i节点到锚杆两端经历的第k个锚杆单

元所在灌浆和岩土介质之间接触面的凝聚力；rb为锚杆半径；Lk为i节点到锚杆两端经历的第k个锚杆单元的长度。

而锚杆单元高斯积分点的屈服应力由两节点屈

第12期 刘文连等：强度折减有限元法中锚杆计算模型研究 4023

服应力线性插值而得。在锚杆两端对应的节点屈服应力为0。为了提高锚杆的承载能力和边坡整体稳定性，经常在锚杆自由端部加有托盘（或板条）结构，这时锚杆单元i节点到锚杆外端的承载力由

∑2πrkkkk

bcFL变为P+∑2πrk

bcFL（P是板条结构能够

k

给锚杆端部提供的最大拉力）。而同时考虑卸剥、拔

出、拉断形式的锚杆单元i节点屈服应力σi

bp

可表示为

i

i

i

i

σi

min(F1

,F2,F3)bp=F⎫πr2=2

⎪bπrb⎪Fi=∑2πrmLm

, 拔出⎪1bcF⎪

m⎬ （4）

Fi, 拉断⎪

2=T⎪

Fikk

⎪3=P+∑2πrkbcFL,卸剥⎪

⎭

3.3 考虑锚杆强度折减的强度折减法

参照文献[7]，式（2）可做如下变换：

∑⎡⎢c⎛wijcosα+Fksinβ⎤ij ⎞k⎟ tan ϕij⎥lij+Gf⎢ij+⎜⎣⎜⎝lijlij⎟⎠⎥⎦s=

∑⎛⎜w≈ij⎜sinα⎞

⎝lij⎟lijij⎟⎠

∑(cij+σijtanϕij)lij+G

∑τ≈

ijlij ∫l

0(c+σtanϕ)dl+G

∫l

0τij

dl （5）

式中：G=∑Fkcosβk。

将式（5）两边同除以fs，则变为

∫ l⎛ 0

⎜ctanFk1=⎝f+σϕ⎞

s

f⎟dl+∑cosβks⎠fs

∫l

l

=0τijd

l

∫ 0(c′+σtanϕ′)dl+∑Fk′cosβk

（6）

∫l

0τijdl

c′=

c,tanϕ⎫

f tanϕ′=sfs

⎪⎪σiFi′min(Fi1′,Fi2′,Fi3′)⎪bp′=πr2=

⎪bπr2b⎪

mm⎪

⎬′∑2πr （7） Fi=m

bLcF

⎪1fs

⎪2πrkk

⎪FibLcF2′=Ti

f, F3′=P∑⎪

k

f+

⎪sfs⎪⎭

式（6）左边等于1，表明当强度折减fs以后，坡体达到极限状态。式（7）表明，考虑锚杆锚固效果的强度折减法只要按照式（7）进行折减，有限元强度折减法在本质上与传统方法是一致的。式（7）即是考虑锚杆强度折减的强度折减方法。

4 计算比较

按照以上原理研制得到考虑锚杆强度折减的三维强度折减有限元法程序，程序中采用了锚杆单元节点可以和实体单元节点不重合技术，岩土体本构模型采用非关联准则的弹塑性本构模型，使用Mohr-

Coulomb准则作为屈服函数，用Drucker-Prager方程作为塑性势函数，具体见文献[6－8]。下面介绍本文算法在一个简单的岩坡锚杆支护方案计算中的应用实例，以证实本文算法的有效性。

某岩石边坡工程，边坡比例为1:1，边坡高 20 m，由岩块A、岩块C和软弱夹层B组成，软弱夹层从边坡正中部延伸到坡顶，夹层和水平面夹角为33.51°，水平宽度为0.1 m，几何模型见图3。计划采用一排锚杆支护，锚杆长度为6 m，直径为

0.1 m，锚芯间距为4 m，锚杆外端点坐标为（15.1，2，15），锚固角为22.51°（不考虑施工因素影响）。

(20,0,20) (25,0,20) CB (10,0,10)(0,0,0)A

图3 几何模型

Fig.3 Geometry model

边坡的有限元网格模型见图4。边坡宽4 m，约束条件如下：坡体侧面约束相应的水平位移，坡体底面约束3个方向的位移。岩体是各向同性体，参数见表1，锚芯弹性模量为29 GPa，半径为0.05 m。

图4 三维有限元模型

Fig.4 3D finite element model

4024 岩 土 力 学 2010年 表1 岩体力学参数表

Table 1 Parameters of rock mass for mechanical analysis

区域

c ϕ ψ E γ / kPa

/ (°)

/ (°)

/ GPa

/ (kN/m3)

ν

A 800 52.4 0 36.000 27 0.26B 10 20.0 0 0.196 20 0.30C 300 26.0 0 1.500 26 0.30

计算过程中锚杆杆壁和岩土介质之间接触面的凝聚力和所穿过岩石凝聚力相等。本文采用笔者研制的三维强度折减有限元法程序进行了3个方案计算。

现在采用极限平衡法进行分析。假设边坡滑动面沿着软弱岩层B层中面，滑动面以上岩体重量经计算为258 822.40 kg，滑动面凝聚力提供的抗滑力为702 886.38 kN，滑动角为33.51°。对于不锚固情况，采用本文研制的强度折减有限元法程序计算结果稳定性安全系数为1.034，而极限平衡法计算为

F=F

sQ

=(702 886.38+258.822 4×1 000×

9.8cos33.51°tan20°)/(258.822 4×1 000×

9.8sin33.51°)=1.061 6

式中：F为抗滑力；Q为下滑力。可见对不锚固情况本文程序计算结果和极限平衡法计算结果近似相等。

对锚固角为22.51°情况，本文程序计算结果安全系数为1.256。按照图5成果通过插值得到B区锚杆单元中点处锚杆承载力为311 738.00 kN，极限平衡法计算为

FF

s=Q

=[702 886.38+258.822 4×1 000×

9.8cos33.51°tan20°+311 738.00×

(cos22.51°+sin22.51°tan20°)]/(258.822 4×1 000×9.8sin33.51°)=1.061 6

本文程序计算结果误差为3.36%，可见本文建议的方法和解析解结果一致。在有限元法计算中，采用高斯点积分代替精确积分，把位移复杂分布曲线在单元内简单化为与坐标的一次或者二次函数、本文提出的锚杆沿壁摩阻力近似考虑方法等近似处理过程，这些近似处理过程是产生系统数值计算误差的主要原因，而考虑强度折减以前按照本文算法锚杆节点抗拉屈服应力分布如图5所示。这个锚固情况下的水平位移等值线图见图6，可见滑移面在

B区，与实际情况符合。

计算结果表明，不对锚杆屈服应力进行折减时，本文程序计算结果安全系数为1.324，在假设锚杆为弹性材料（不考虑锚杆屈服）时本文程序计算结果

安全系数大于5.0（安全系数等于5.0时还没有不收

敛迹象），因此可以下定论，本文第二部分提到的当前强度折减有限元法中锚杆计算模型都高估了锚杆的作用，以此计算模型计算稳定性安全系数偏高，以此计算结果指导边坡加固容易导致边坡滑坡事故。因此，本文第二部分提到的当前强度折减有限元法中锚杆计算模型是有缺陷的，不适用于边坡稳定性评价。

6050 aPM40 /力应30服屈201000

1

2

3

4

5

6

离锚杆头距离/ m

图5 屈服应力分布图

Fig.5 Distribution of yield stresses along anchor

810.0 -510. 00- 71-00 .150×00-.10--

(a) 等势线分布

(b) 对应的有限元网格图

图6 水平位移等值线图(单位: m)

Fig.6 Contours of horizontal displacements(unit: m)

5 计算模型应用

某路堑土质边坡工程，几何模型见图7，从上往下依次是黏土层、碎石土层和基岩层，从下往上三段人工开挖边坡比例分别为1:0.75、1:1和1:1.25。土体力学参数如表2所示，因为要求评估边坡在天然状态和饱和（降雨）状态下边坡的稳定性，因此对可能在滑坡体内的岩土材料提供天然状态和饱和

第12期 刘文连等：强度折减有限元法中锚杆计算模型研究 4025

状态2组参数。要求评估边坡在天然状态和饱和（降雨）状态下边坡的稳定性，如不稳定采用锚杆支护。经计算边坡不稳定，采用锚杆加固。加固方案是：锚杆直径为0.11 m，排拒为1 m，倾角为15°，最下面一段人工开挖边坡锚杆长度分别为10、15、

20 m，第2段人工开挖边坡最下面两根锚杆长度为20 m，其余锚杆长度为25 m，如图7所示。锚杆壁灌浆和黏土、碎石土间接触面凝聚力分别取为30 kPa和60 kPa。

图7 几何模型 Fig.7 Geometry model

表2 力学参数表

Table 2 Mechanical parameters

土层名称

c ϕ ψ E γ

/ kPa

/ (°)

/ (°)

/ MPa / (kN/m3)

ν

黏土层 30 10 0 25 18 0.35碎石土层 23 25 0 50 19 0.30饱和黏土层 20 9 0 25 19 0.35饱和碎石土层 17 21 0

50 20 0.30

基岩 36 000 27 0.26

边坡的有限元网格模型见图8，边坡宽1 m。约束条件：坡体侧面约束相应的水平位移，坡体底面约束3个方向的位移。锚芯弹性模量为0.725 GPa，半径为0.055 m（模量已适当降低）。在三维强度折减有限元计算过程中，先对包含锚杆的边坡体进行重力作用下的弹性变形计算，以求取地应力分布；再进行折减强度情况下的塑性变形计算，以确定稳定性安全系数及对应的极限平衡状态下的塑性变形、应力分布等特征。如果采用锚杆真实弹性模量，将由于锚杆弹性模量和土体弹性模量差异很大而使得第一阶段计算土坡得到的地应力分布和实际地应力分布有比较大的差距而失真，导致稳定性安全系数计算结果产生比较大的误差。为了减少这种误差，这里采用了降低锚杆弹性模量的处理方法，即降低锚杆弹性模量进行强度折减计算。当再次降低锚杆弹性模量进行强度折减，计算结果变化不大，即认为前一次采用的锚杆弹性模量为计算应当使用的锚杆弹性模量。采用传统极限平衡法和笔者研制的三

维强度折减有限元程序进行了3个方案计算，Fs的增量取0.002，最大迭代次数是800，所得结果列于表3，锚固滑弧面如图9、10所示。经比较可见，本文有限元法和极限平衡法计算结果误差比较小。计算结果表明，边坡稳定，建议锚固方案可行。

图8 三维有限元模型

Fig.8 3D finite element model

表3 各种分析方法计算的边坡稳定性安全系数

Table 3 Safety factors with different analytical methods

状态 极限平衡法 (Slide软件) 毕肖普法

简布法

强度折减有限元法

不锚固天然

0.988 0.956

0.970 锚固天然 1.331

1.270 1.310 锚固饱和 1.043

1.008

1.048

2.0991.331 1.766 2.178(a) 用Bishop简化确定滑弧面

(b) 强度折减有限元法根据屈服区域确定滑弧面

图9 天然状态下滑弧面

Fig.9 The failure faces of slope without draining

6 结 语

本文总结当前强度折减有限元法中锚杆计算模型特点并指出不足之处，介绍了得到公认的边坡稳定性分析极限平衡法商用软件中锚杆强度模型，由

4026 岩 土 力 学 2010年

1.4061.043 1.320 1.638

(a) Bishop简化

(b) 强度折减有限元法根据屈服区域确定滑弧面

图10 饱和状态下边坡滑坡滑动面

Fig.10 The failure faces of slope under draining

商用软件中锚杆强度模型和瑞典法启发建立了适用于强度折减有限元法的锚杆强度模型和考虑锚杆强度折减的强度折减方法，并通过和解析解计算结果、商用软件Slide计算结果进行对比验证了本文提出方法的适用性。

参 考 文 献

[1] 蔡永昌, 朱合华, 李晓军. 一种用于锚杆支护数值模拟

的单元处理方法[J]. 岩石力学与工程学报, 2003, 22(7): 1137－1140.

CAI Yong-chang, ZHU He-hua, LI Xiao-jun. New element form in the numerical simulation of bolt- reinforced structure[J]. Chinese Journal of Rock Mechanics and Engineering, 2003, 22(7): 1137－1140. [2] 崔颖波, 潘懿. 某露天边坡加固方案的稳定性数值分

析[J]. 化工矿物与加工, 2008, 11: 33－39.

CUI Ying-bo, PAN Yi. Numerical analysis on stability of slope consolidated in a open-pit mine[J]. Industrial Minerals & Processing, 2008, 11: 33－35.

[3] 宋开灼, 郭志华. 水电站尾水出口边坡平面有限元分

析[J]. 土工基础, 2008, 22(1): 59－63.

SONG Kai zhuo, GUO Zhi hua. Plane finite element

analysis to the slope at tail water outlet of hydro-power station[J]. Soil Engineering and Foundation, 2008, 22(1): 59－63.

[4] 王剑波, 何朋立, 陈新胜, 等. 隧道锚喷参数优选正交

有限元强度折减分析[J]. 青岛理工大学学报, 2007, 28(2): 1－5.

WANG Jian-bo, HE Peng-li, CHEN Xin-sheng, et al. Research on optimization of supporting parameters of road tunnel engineering based on the orthogonal finite element strength reduction experimental method[J]. Journal of Qingdao Technological University, 2007, 28(2): 1－5.

[5] 荣冠, 王思敬, 王恩志, 等. 强降雨下元磨公路典型工

程边坡稳定性研究[J]. 岩石力学与工程学报, 2008,

27(4): 704－711.

RONG Guan, WANG Si-jing, WANG En-zhi, et al. Stability analysis of Yuanmo highway slope under intensive rainfall[J]. Chinese Journal of Rock Mechanics and Engineering, 2008, 27(4): 704－711. [6] 韦立德, 叶志华, 陈从新, 等. 一种锚杆计算模型及其

在岩坡锚杆支护中的应用[J]. 岩土工程学报, 2008, 30(5): 732－738.

WEI Li-de, YE Zhi-hua, CHEN Cong-xin, et al. Effective

calculation method of bolts and its application to rock slopes anchored by bolts[J]. Chinese Journal of Geotechnical Engineering, 2008, 30(5): 732－738. [7] 赵尚毅, 郑颖人, 时卫民, 等. 用有限元强度折减法求

边坡稳定安全系数[J]. 岩土工程学报, 2002, 24(3): 343－346.

ZHAO Shang-yi, ZHENG Ying-ren, SHI Wei-min, et al. Analysis of safety factor of slope by strength reduction FEM[J]. Chinese Journal of Geotechnical Engineering, 2002, 24(3): 343－346.

[8] 韦立德, 杨春和, 高长胜. 基于三维强度折减有限元的

抗滑桩优化探讨[J]. 岩土工程学报, 2005, 27(11): 1366－1368

WEI Li-de, YANG Chun-he, GAO Chang-sheng. Optimization of slide-resistant piles based on strength reduction method with 3D FEM[J]. Chinese Journal of Geotechnical Engineering, 2005, 27(11): 1366－1368.