基于Hoek-Brown准则的隧道围岩塑性变形研究

付大喜，王小伟，翁效林

(1.河南省交通规划设计研究院股份有限公司，河南 郑州 450000；2.长安大学 公路学院，陕西 西安 710000)

地下隧道开挖打破了岩体中的初始应力平衡状态，围岩应力出现应力重分布现象。由于围岩应力重分布的影响会出现2种状况：一是围岩处于弹性状态，弹性理论依然适用，此时隧道可以不进行支护；二是围岩的重分布应力超过岩体的弹性极限状态进入塑性状态，围岩周围出现弹塑性区域，处在塑性区的围岩由于软化作用出现塑性流动区，围岩进步被破坏。因此对围岩弹塑性区的应力和位移分析是硐室支护的重要参考依据。

许多学者对深埋圆形隧道的弹塑性解进行了深入研究。张常光等[1]归纳了圆形隧道平面应变条件下8种强度准则新的围岩应力和位移解，并对新解和进行了比较分析。王宏伟等[2]为探究围岩破碎区分布范围和应力分布特点，建立非静水压力下围岩应力状态模型，通过数值模拟技术发现侧压力系数对围岩塑性区半径和影响范围随隧道埋深呈反比。张道兵等[3]为获得更精准的深埋隧道围岩压力，运用极限分析法所得结果与工程实际对比发现围岩压力随侧压力系数和黏聚力的增大而减小。关晓迪等[4]基于统一强度理论推导了非对称荷载作用下围岩弹塑性区应力表达式，结果表明屈服准则系数对围岩弹塑性区交界面处的径向和切向应力均有较大影响。常帅斌等[5]针对隧道开挖后产生的塑性区围岩稳定性问题，基于统一强度理论分析了围岩力学参数和支护压力对塑性区的影响。赵翔等[6]分析了由于开挖导致岩体物理力学参数改变对围岩弹塑性区分布的影响，发现隧洞开挖后硐室应力由单向受力状态转向双向受力状态，并以此为依据划分弹塑性应力区。

但上述成果大多是基于非线性的M-C准则，大量工程实验和数值模拟结果表明非线性的H-B准则能更好的描述围岩的应力和变形状态[7]。朱艳峰等[8]基于应变软化模型，推导了非圆形隧道塑性区半径的表达式，并通过现场声波试验验证其准确性。崔建锋等[9]采用室"楔形塌落体+圆弧转动体"破坏模型探究动荷载对硐室稳定性的影响，发现较大地质强度指标和围岩单轴抗压强度有利于硐室稳定，而围岩重度和隧道尺寸的增大会增大硐室垮塌的风险。王睿等[10]将现场实测的竖向应力和水平应力带入推导的理想化松动圈公式，得到了较为准确的松动圈分布规律。贺耕夫等[11]综合考虑了冻胀深度和衬砌的弹性模量对围岩塑性区的影响，推导了寒区隧道围岩的塑性区半径公式，结果表明在一定范围内岩体的单轴抗压强度对塑性区半径影响巨大。邢凯等[12]基于H-B准则推导了渗流条件下水下隧道开挖围岩的应力分布特征公式，并将应力分布区划分为强影响区、弱影响区、无影响区；徐强等[13]基于隧道衬砌支护与围岩相互作用的思想，将塑性区内的自重应力看作竖向松动压力，提出了新的围岩压力计算公式；朱明俊等[14]为同时考虑开挖面的空间效应和渗流作用，推导了综合考虑二者作用的隧道开挖弹塑性解，并根据计算结果建议提前施作支护可有效减小围岩变形。陆晓清等[15]将中间主应力引入到H-B弹塑性解析式中，并将新的弹塑性解析式与M-C准则和传统的H-B准则的弹塑性解进行对比，发现3种准则在围岩应力和围岩的差异主要集中在1～4倍的硐室半径范围内。

综上所述，目前基于H-B屈服准则推导的圆形硐室弹塑性解取得了丰富的研究成果，但对缺少对于综合考虑H-B准则参数地质强度指标GSI、岩石软硬强度指标 和扰动参数D隧洞围岩弹塑性区和塑形半径的影响。因此本文基于H-B屈服准则,考虑应力重分布作用(不考虑地下水作用)的影响，推导了深埋圆形隧洞弹塑性应力计算方法。在此基础上,对不同围岩支护力和不同初始地应力对围岩塑性区的位移和塑性半径的影响进行深入的分析。并将理论结果与现场试验结果进行对比分析，为后续类似施工提供借鉴参考。

1 基于Hoek-Brown强度准则求解圆形硐室围岩的弹塑性区应力及位移

1.1 Hoek-Brown强度准则

Hoek-Brown准则是霍克等通过大量的岩石三轴试验提出的非线性岩体破坏准则，其综合考虑了岩体结构面和应力状态对岩体强度的影响，对于各向异性的岩体或深埋破碎的岩体均适用。其表达式为：

(1)

式中：σ1、σ2分别为岩体破坏时的最大、最小有效主应力,MPa；σci为完整岩块的单轴抗压强度,MPa；mb、s、α为Hoek-Brown准则的相关参数，其中mb为岩石材料常数，取值范围在0.000 000 1～25之间，s为岩体质量有关常数，其取值范围为0～1，对破碎岩体取0，完整岩体取1.0。mi和s可依据室内三轴试验成果计算得到或直接按Hoek建议表格直接取经验值。mb、s、α参数由下式表示为：

(2)

(3)

(4)

式中：mi反映岩石的软硬程度；GSI为地质强度指标；D为扰动参数，值为0～1。施工无扰动时D取0，扰动严重时D取1。

1.2 硐室围岩的弹塑性区应力

基本假设：(1)不考虑围岩的蠕变和剪胀，假设岩体为均质的(2)假设侧压系数λ=1，原岩处于静水压力状态；(3)隧道断面为圆形,且隧道长度为无限长岩体质量保持均一，于是可以采用平面应变问题选取隧道内任意截面作为其代表研究；其简化的力学模型如图1所示。

图1 隧道围岩弹塑性区域

1.2.1 弹性应力解答

假定隧洞周围岩体为均质体，忽略自重，建立如下平衡微分方程：

(5)

式中：σr和σθ分别为径向应力和环向应力，规定拉应力为正，压应力为负。文中所有应力均为有效应力，后面不再重复说明。应力边界条件为：

σr=-pa,r=ra

(6)

σr=-p0,r=βra

(7)

式中：β为一很大的数。根据上述平衡微分方程和边界条件，由文献[16]当ha=h0时，即不考虑渗流作用，求得应力弹性解答为：

(8)

(9)

式中：K1、K2计算表达式如下：

(10)

(11)

(12)

(13)

式中：E为弹性模量，μ为泊松比。

1.2.2 弹塑性应力解答

考虑到洞室周岩体塑性变形性能发挥引起应力重分布，引入应力调整系数λ，重分布后弹性区(r≥Rp，Rp为塑性区半径)岩体应力为：

(14)

(15)

对于塑性区(r≤Rp)岩体，其大小主应力分别为σ1=-σθ，σ3=-σr，假定岩体服从Hoke-Brown准则，将大小主应力代入式(1)得：

(16)

由于塑性区岩体仍满足平衡微分方程，故将上式代入式(5)，可得：

(17)

由式(17)可转化为：

(18)

对上式积分后，将洞壁处应力边界条件式(6)代入上式得：

(19)

其中K5为待定系数。

(20)

(21)

(22)

(23)

(24)

(25)

(26)

2 对比分析

为了验证本文推导结果的准确性，将本文结论与既有研究成果进行比较。潘阳等[18]基于Hoek-Brown准则和平衡微分方程推导了圆形巷道弹塑性区的应力解。

参考文献中算例，在围岩质量不同的石灰岩中开挖圆形硐室，已知r0=4 m，隧洞覆盖深度h=200 m。因此p0= 5. 4 MPa。岩石材料常数mi=8,开挖扰动系数D=0. 6，完整岩石的单轴抗压强度σci=50 MPa，GSI=20。

通过图2可以看出在塑性区的切向应力本文采取的方法与文献[18]的方法绘制图像的变化趋势趋于一致，而塑性区的径向应力本文计算结果明显大于文献[18]的原因是由于参数a的取值影响，与已有的参考文献[18]相比，本文计算的塑性区径向应力较真实值更为接近。在弹性区无论是径向应力还是切向应力都比文献[18]更快趋近初始应力。而且本方法计算出的塑性区径向应力最大值(绝对值)明显小于文献[18]，切向应力在弹塑性交界处无突变发生(文献[18]在交界处发生突变是由于塑性区采用Hoek-Brown算法，弹性区采用厚壁圆筒法进行处理)。因此本方法更适合计算深埋隧洞的应力变形解。

图2 GSI=20时硐室弹塑性区应力变化

3 参数分析及讨论

隧道开挖半径R取为2.5 m，埋深H取600 m，初始地应力p0=0.025H=15 MPa,支护反力pi取为0。分析质量好、质量一般和质量差的岩体在岩体质量指标GSI分别取值30、50、70时围岩弹塑性区应力和孔壁处位移的变化规律。岩石材料常数mi分别取值为8、12、16，分别代表软岩、硬岩和中硬岩3种不同坚硬程度的岩石。计算时扰动参数D取0。

岩体变形模量由下式给出：

对于GSI分别取30、50、70的不同岩体，泊松比 分别取值为0.35、0.30和0.25。经计算，孔壁处围岩位移和塑性区半径分析所需的参数如表1所示。

表1 计算参数

3.1 初始地应力的影响

图3和图4分别为初始地应力分别取5 MPa、10 MPa、15 MPa、20 MPa、25 MPa、30 MPa条件(未支护)下塑性区半径和孔壁处围岩位移的变化图(计算参数从表1中获取)。从图3中可以看出在GSI=30相同时，硬岩和中硬岩的塑性区半径与初始应力呈一次线性函数变化，而软岩的塑性区半径与初始应力呈二次线性变化，且硬岩的塑性区半径变化速率最慢，中硬岩次之，软岩变化速率最快。对于硬岩与中硬岩在初始地应力小于15 MPa时其塑性半径变化速率缓慢，当初始应力大于15 MPa时，中硬岩的塑性半径增速开始变化而硬岩的增长速率基本不发生变化。出现这种现象的原因主要是岩石的材料指标与岩石的单轴抗压强度即岩体本身的性质不同引起的。对于同种围岩(σci=20 MPa)当GSI≥50时初始地应力对塑性区半径明显变小。

图3 不同初始地应力下塑性区半径的变化

图4 不同初始地应力下孔壁处塑性区位移变化

对于软岩的塑性半径变化速率过快的原因主要是相比硬岩，软岩隧道开挖后应力重分布调整幅度明显更大。而且在高地应力深埋的条件下，隧道未开挖时，围岩处于静水压力状态，在隧道开挖后，由于卸荷作用和应力重分布，岩体结构面发生剪切滑移破坏，软岩相比其他两种岩体更容易破碎(主要单轴抗压强度较低)，因此其塑性区半径变化速率随地应力的改变变化速率最大。

从图4中可以看出，随着初始应力的增大，隧道孔壁处的围岩也随着增大，软岩增长速度最快，中硬岩次之，硬岩增长速度最慢。在未支护条件下，地应力的变化对于孔壁处围岩都有较大影响，尤其是对于软岩(GSI=30时)，在p0=30 MPa时孔壁处围岩达到18 m多，且硬岩孔壁处的围岩也达到了0.2 m。由此表明对于GSI=30的岩体开挖后立即进行支护的必要。

3.2 支护力与初始地应力比例的影响

为探究支护比(围岩支护力与初始地应力比值)对塑性区半径和应力的影响(扰动参数D=0)，选取表1中计算参数行计算。

图5和图6分别为在支护比分别取0.00、0.05、0.10、0.15、0.20、0.25、0.30、0.35、0.40、0.45、0.50的情况下,孔壁处塑性区半径变化和位移变化图。不同支护力与初始应力比(以下简称支护比)下，中硬岩孔壁处的塑性半径随GSI的增大而减小。当支护比小于0.30时，其塑性区半径随支护比的增大而减小的速率明显加快，且GSI越小，加快速率越明显。当支护比大于0.30时，其塑性区变化速率明显减慢，而且此时GSI的变化对于塑性区半径影响较小。因此，在后期进行支护时当支护比取0.30时进行支护即可满足大部分设计要求。

图5 不同支护比下孔壁处塑性区半径变化图

图6 不同支护比下孔壁处塑性区位移变化图

当GSI=30时，支护力对软岩孔壁处位移影响最大，未支护时位移为1.52 m，支护比为0.30和0.50时，孔壁处位移分别为0.17 m和0.08 m。对于中硬岩和硬岩的影响较小。且GSI=50的软岩孔壁处位移与GSI=30中硬岩的位移大小基本趋于一致，由此岩体质量对隧道变形有较大影响。GSI=70支护比取0.00和0.50时，孔壁处位移为0.04 m和0.02 m,此时支护力对于围岩变形影响不大。

4 工程实例

栾卢高速某隧道位于三门峡市横涧乡，隧道左起ZK52+627—ZK58+124，全长5 508 m，右线里程K52+627—K58+124，全长5 497 m，隧道净宽10.25 m，净高5.0 m，最大开挖宽度12.3 m，最大埋深395 m，两线隧道间距为9 m～47 m，为分离式隧道。

以隧道ZK54+060为典型断面，围岩等级为Ⅳ级，该断面为中—微风化片岩，灰白色，变晶结构，片状构造，岩石裂隙发育，含水丰富，为地质断层影响带。围岩基本质量修正值[BQ]=274.6。隧道净空断面如下图7所示，隧道净宽10.25 m，净高5.0 m，断面处埋深为280 m，内轮廓半径R=5.45 m。初期支护参数如表2所示。

图7 隧道净空断面图(单位:mm)

表2 初期支护设计参数

4.1 理论计算

围岩等级S取4，[BQ]上和[BQ]下根据《公路隧道设计规范》(JTG 3370—2018)分别取330和260。竖向应力取初始应力，围岩压力按式(27)取值，围岩修正级别系数[S]取值按式(28)计算，Ⅳ级围岩承担初支荷载占比在60%～80%，取初支承担围岩载荷为70%将计算的到的q=0.202 MPa代入，假定支护压力为Pi=0.7q=0.141 MPa。

q=γh=26.5×0.45×2s-1[1+i(B-5)]

(27)

(28)

表3 H-B准则的围岩力学参数表

(1) 拱顶。在计算拱顶和拱腰的松动圈范围时，将隧道假定为圆形。同时考虑到侧压力系数对松动圈的影响，为保证公式适用性，在计算拱脚处塑性区半径时，将按实际侧压力系数λ=1.2进行计算。

(29)

围岩初始地应力P0=γ·H=7.42 MPa，地质强度标根据文献[21]将未修正的岩体基本质量BQ=314.6代入到式(38)中得出RMR=34(此处取整)，利用RMR求得修正后GSI=29，再将表3中参数带入式(2)中求得s和mb取值分别为7.74×10-5和0.680。将以上参数代入式(29)和式(35)求得塑性区半径RP和塑性区位移分别为8.03 m和71.1 mm。塑性区厚度为：

D=RP-R=2.58 m

(2) 拱腰。考虑到现场实测最大水平主应力范围为9.09 MPa～12.53 MPa，最小水平主应力值范围为7.04 MPa～8.41 MPa，拱脚处初始应力p0=λ·γH=8.90 MPa,同上可求得塑性区半径RP和塑性区位移分别为8.65 m和100.1 mm。

D=RP-R=3.20 m

4.2 现场监测及分析

为验证推导结果准确性，在栾卢高速某隧道ZK54+060断面拱顶及拱脚处布设多点位移计监测围岩深部位移。拱顶测点H1孔深6 m，共设置8个测点，分别为0.8 m、1.3 m、1.8 m、2.3 m、3.0 m、4.0 m、5.0 m、6.0 m；拱脚处布设测点H2、H3，孔深和测点布置和拱顶相同。

现场监测的围岩深部位移曲线如图8所示。由图8可得处围岩深部位移发展主要分为三个阶段，即剧烈发展阶段、快速发展阶段以及缓慢发展阶段。由图8(a)中可知，拱顶处围岩深部位移累计沉降最大处为1号测点位移为9.72 mm，最小累计沉降为8号测点位移为1.39 mm。按位移分布区域可划分为严重扰动区(测点累计沉降量h>7 mm)、中度扰动

图8 ZK54+060断面围岩深部累计位移时程曲线

区(测点累计沉降量7 mm≥h≥3 mm)和微扰动区(h< 3 mm)3个层次。如图9所示，根据围岩位移沉降曲线斜率可将松动圈划分为3个区域，靠近围岩洞壁的松动区，位移较大的升高区以及围岩变形量最小的弹性区。测点H1塑性区范围是3.0 m～4.0 m，左拱脚塑性区范围是4.0 m～5.0 m，右拱脚塑性区范围是4.0 m～5.0 m。

图9 围岩内部位移曲线

如图8(a)所示，1、2号测点变形速率和累计位移远大于其他测点，厚度为1.8 m。3、4、5号测点处于中度扰动区，其范围为2.3～4.0 m，塑性区半径就处于此范围。7号测点为微扰动区，也是隧道开挖后影响范围最大处，8号测点外为弹性区，基本不受开挖爆破的影响，弹性区处在5.0 m以外。

按上述分区方法，左右拱脚处围岩深部位移规律也可划分为严重扰动区、中度扰动区和微扰动区如图8(b)和图8(c)所示。与拱顶处相比，左右拱脚塑性区影响范围更大，且影响区域厚度并不一致，影响区厚度分别为5.0 m和4.0 m。微扰动区为6.0 m和5.0 m以外区域。左右拱腰塑性区范围不一致的主要原因是此断面处于断层破碎带影响区域，围岩节理裂隙发育不均衡，且受较大水平构造作用明显。

理论计算与现场实测结果对比如表4所示，拱顶和左右拱脚计算的塑性区范围均在现场实测范围内，表明推导的不考虑渗流作用的H-B准则在计算围岩塑性区范围具有一定的准确性。

表4 结果分析对比

此外，通过对围岩深部位移累计观测值与理论推导结果相比，发现拱顶和拱脚围岩内部0.8 m处变形量分别占理论计算总变形量的13.6%和10.9%，这表明围岩变形在距离孔壁0.8 m范围内变化剧烈，因此早期围岩加固时可重点对此范围内岩体进行加固。

5 结 论

(1) 在支护压力与初始比值较小时，隧道洞壁附近的围岩应力相对较小，但孔壁处位移和塑性区半径都很大。当支护压力较小时，尤其是支护比小于0.1时，GSI的改变对孔壁处位移和塑性区半径影响很大。支护比大于0.3时，塑性区半径和孔壁处位移变化速率明显减缓，因此支护比等于0.3可作经验参数指导施工。

(2) 初始地应力对于软岩隧道(对于本文σci<20 MPa)的变形尤为显著，尤其是初始地应力大于15 MPa时，对于硬岩(σci=80 MPa)与中硬岩(σci=50 MPa )在初始地应力小于15 MPa时其塑性半径变化速率缓慢，当初始应力大于15 MPa时，中硬岩的塑性半径和位移变化明显而硬岩基本不发生变化。

(3) 通过现场对围岩深部累计位移监测，按位移分布区域可划分为严重扰动区、中度扰动区和微扰动区3个层次。且左右拱脚处的塑性区厚度与拱顶塑性区厚度各不相同，这主要是由于围岩节理裂隙发育分布不均，水平构造应力对围岩变形影响明显。在后续初期支护时要加强对拱脚处支护。

(4) 围岩发生剧烈变形主要集中在距围岩表面0.8 m处。x=0.8 m处拱顶和左右拱腰的累计沉降值分别为9.72 mm、11.25 mm和10.58 mm，约占理论总变形量的11.5%，在初期支护时可重点对距围岩表面0.8 m范围内岩体进行加固。

本试验可以对围岩深部位移测点在靠近孔壁3.0 m～5.0 m范围内可适当增加测点，以判定具体围岩塑性区范围。