基于3DEC节理岩体隧道开挖与支护的数值模拟

周 斌，何明卫，祝方才，邬明亮，刘思远，过 江

（1.湖南工业大学 土木工程学院，湖南 株洲 412007；2.中南大学 资源与安全工程学院，湖南 长沙 410012）

0 引言

近十几年来，随着国民经济的稳步增长和工程技术的不断进步，大型山岭隧道、水工隧道和矿山巷道建设均取得了突破性成果[1-3]。但是，隧道工程往往会穿越发育良好的节理裂隙岩体地区，由于节理裂隙岩体等复杂结构面的存在，使得在此类岩体中进行隧道开挖的难度极大。在开挖失去完整性、整体强度和自稳能力的岩体隧道时，特别容易发生结构失稳、坍塌等安全事故。

已有研究表明，对节理裂隙岩体隧道围岩的变形及稳定性影响较大的因素，主要有隧道上覆岩层的厚度[4]、节理面的特征[5]、施工技术[6]和支护手段[7]等。目前，研究此类问题的学者和工程人员，常用有限单元法、有限差分法和离散单元法等方法对隧道、巷道等地下工程的稳定性进行分析[8-9]，且发现对于含有复杂结构面的地下工程，采用离散单元法更能模拟岩体的真实情况[10-11]。因此，本研究拟以国内某山岭高速公路隧道为研究对象，基于离散单元3DEC软件，考虑裂隙水压力的作用，建立节理裂隙岩体隧道的开挖与支护模型，并在有无裂隙水两种工况下，分别分析了隧道在开挖与支护过程中岩体及节理的响应情况，以期为今后研究节理岩体等相关问题提供有益的参考。

1 基本概况

1.1 工程概况与模型建立

本研究以国内某山岭高速公路隧道为研究对象，该隧道长1 538.9 m，隧道垂直埋深为50.0～200.0 m，局部埋深超过200.0 m，隧道围岩以Ⅱ、Ⅲ、Ⅳ、Ⅴ类围岩为主，并且分别占隧洞长度的1%,50%,38%,11%。隧道断面呈马蹄形，隧道断面的具体尺寸如图1所示，图中尺寸单位为cm。

图1 隧道断面尺寸图Fig.1 Sectional dimensions of the tunnel

模拟段隧道存在3 组节理穿过，计算模型岩石节理的倾向分别取5°,30°,255°，倾角分别取83°,5°,79°，节理间距分别取13.0,9.0,11.0 m，3 组节理的空间位置形态如图2所示。隧道计算模型的尺寸为40 m×40 m×32 m。

图2 节理模型图Fig.2 Joint model illustration

在所建立的计算模型中，隧道岩体采用弹塑性模型，屈服条件符合Mohr-Coulomb 强度准则，其屈服函数如下：

式（1）（2）中：fs为剪切应力；

ft为张拉应力；

σ1、σ3分别为最大、最小主应力；

c为黏聚力；

Rm为抗拉强度；

Nφ为特征参数，且

其中φ为内摩擦角。

当岩体内部某点的应力满足fs＜0 时，其发生剪切屈服；当岩体内部某点的应力满足ft＞0 时，其发生张拉屈服。结构面（节理）采用面接触-库伦滑移模型来描述其滑移与破坏，其计算原理如下。

在模拟过程中的每一计算步内，节理接触面的法向力矢量Fn和切向力矢量Fs的计算方法如下：

kn、ks分别为节理面的法向刚度和切向刚度；

ΔXn、ΔXs分别为接触面的法向位移增量和切向位移增量；

Ac为接触面面积。

当节理接触面没有发生滑移（弹性阶段）时，节理面承受的最大拉应力和最大剪应力分别如下：

式中Rq为结构面抗剪强度。

当节理面被拉断或者发生剪切破坏时，节理的最大拉应力、黏聚力和抗拉强度均为零，最大剪应力为Rqtanφ，而节理面的法向力和切向力将发生相应的变化。当节理受剪切破坏时，节理面的切向力为残余强度切向力，法向力不变。计算模型的岩体与节理裂隙参数取值见表1。

表1 模型岩体与节理参数与取值Table 1 Mechanical paremeters and values of rock mass and joints

岩体中的渗流本质为水力耦合问题，即水的渗透压力与岩体的力学变形反复作用，进而达到一种平衡状态。裂隙渗流遵循立方定律，即认为岩体不透水，渗流只在裂隙中发生，立方定律认为单裂隙流量Qi与裂隙开度e的三次方成正比，即：

式中：μ为裂隙水的动黏滞系数；

Ji为水力梯度；

g为重力加速度；

ρ为裂隙水密度。

对于实际岩体而言，节理面的特性、填充物质以及岩体所处应力状态等因素往往会对裂隙开度产生影响，并且会进一步影响裂隙水流的运动，导致立方定律不满足实际工程要求。在此基础上，C.Louis 等[12]提出利用节理粗糙度对立方定律进行修正，即

式中：Ra为节理的相对粗糙度，且

其中，N为裂隙总数，b为裂隙等效隙宽。

1.2 3DEC 软件介绍及其特点

3DEC 软件是ITASCA 公司开发的，以离散单元法为基础模拟非连续介质力学行为的三维计算分析程序，也称非连续性方法。其以拉格朗日为算法基础，适合模拟静态或动态作用下离散介质（如节理）的受力以及破坏状态，将非连续介质视为离散块体集合，以间断面充当块体间的边界条件，允许块体沿间断面产生大位移或者旋转。3DEC 软件以其能更好地模拟块体间的分离、旋转和块体塌落过程，并且在块体塌落后能与新的块体实现接触和进行力学解算等优势而被国内外学者广泛使用，但是在其实际工程应用中还存在以下不足：

1）由于3DEC 软件缺少曲线建模命令，在针对砌体单元开挖内侧壁模拟时，以折线代替曲线，因而难以确保其合理性，如图3所示；

2）利用较多的节理单元（jset）对模型进行切割后，往往会产生较多的不规则异形体，以至于混淆块体，也会对模拟隧道开挖造成较大的影响，甚至会出现计算错误，如图4所示。

图3 以折线代替曲线示意图Fig.3 Illustration of polylines instead of curves

图4 不规则异形体示意图Fig.4 Illustration of the irregularly-shaped body

针对上述3DEC 软件的不足，科研工作者们通常借助ANSYS、Rhino 等软件的前处理功能、自身的单元划分、忽略微小节理的影响等手段，以减小因模拟而造成的误差。本研究主要根据现场的地质情况，对开挖临界面进行细分，以减小折线对砌体单元开挖时的影响。同时，选取对隧道开挖影响较大的节理，忽略其微小节理的影响，以消除不规则异形体的产生。

2 数值模拟及结果分析

2.1 数值模拟

在常规静压力作用下开挖隧道后，应力往往会进行重新分布，由于不连续节理的存在和地下水压力的影响，随着裂隙水压力的增大，节理面的剪切强度将会降低，内摩擦角减小[13]，从而会进一步影响隧道开挖时围岩的稳定性能[14]。因此，根据工程的实际情况，针对其在常规静压力荷载作用下和在地下水压力作用下，分别进行隧道开挖与支护模拟。利用节理单元（jset）建立三组节理裂隙隧道模型，通过整体式中空锚杆支护单元（cable）和混凝土衬砌单元（liner）对开挖后的隧道进行支护模拟，隧道锚杆与衬砌模型如图5所示。

图5 隧道锚杆与衬砌模型图Fig.5 Model of tunnel anchor with its lining

在建立的计算模型中，对模型四周以及下部采用位移、应力边界条件，上部为自由边界，并且选取拱顶（点A）、拱肩（点B）、拱腰（点C）、拱脚（点D）和拱底（点E）5 个特征点作为岩体的变量研究监测点。以实际施工顺序，即采用先加载再开挖的施工顺序，模拟隧道开挖与支护的整个施工过程。得到初始地应力后，再按照“台阶法”分为4 组，上下双台阶先后开挖隧道围岩，并且同时对开挖临空面向外以梅花形设置长为3.5 m、直径为22 mm 的砂浆锚杆（HRB400），且进行10 cm 厚度的C20 混凝土衬砌支护模拟，直至隧道围岩变形稳定或者失稳，锚杆参数及其取值如表2所示。

表2 锚杆参数及取值Table 2 Bolt parameters and values

2.2 结果分析

2.2.1 常规静压力荷载下隧道开挖后的响应分析

图6所示为常规静压力荷载下，隧道开挖后未进行任何支护时的围岩竖向位移分布云图。

图6 静压力荷载下未支护时围岩竖向位移分布云图Fig.6 Nephogram of vertical displacement distribution of the surrounding rock without support under the static pressure load

由图6可以得知，常规静压力荷载下，隧道开挖后未进行支护时，由开挖引起的竖向位移主要集中在隧洞的上下两侧，并且由于节理的存在，会使得岩体的竖向位移发生突变。

图7所示为常规静压力荷载下，隧道开挖后未进行支护时的围岩水平位移分布云图。

图7 静压力荷载下未支护时围岩水平位移分布云图Fig.7 Nephogram of horizontal displacement distribution of the surrounding rock without support under the static pressure load

由图7可知，常规静压力荷载下隧道开挖后未进行支护时，由开挖引起的围岩水平位移主要集中在隧道左右两侧，水平位移未因存在节理而发生突变。

表3给出了模拟所得常规静压力荷载下，隧道开挖后，未进行支护与进行锚杆衬砌支护时，围岩的整体位移量。

表3 静压力荷载下围岩整体的位移量Table 3 Total displacement of the surrounding rock under the static pressure load cm

分析表3中的模拟结果数值可以得知，在开挖后对隧道围岩施加砂浆锚杆与混凝土衬砌支护时，隧道围岩的竖向位移、水平位移较未支护时约分别减小了27.74%,24.49%。这一结果说明，本文所采用的支护手段能有效减小围岩位移，保证了开挖时岩体隧道的稳定性。

图8所示为静压力作用下，隧道开挖后未进行任何支护时围岩的最大主应力分布云图。

图8 静压力作用下未支护围岩最大主应力分布云图Fig.8 Nephogram of maximum principal stress distributionof the surrounding rock without support under the static pressure load

由图8所示围岩最大主应力分布云图可以得知，隧道开挖后未进行支护时，围岩的最大主应力主要分布于隧洞周围，大小为0.38 MPa，且最大主应力在节理处发生了突变。

图9所示为静压力作用下，隧道开挖后未进行支护时围岩的受力状态分布云图。由图9所示围岩受力状态分布云图可以得知，隧道开挖后未进行支护时，隧道围岩主要以受剪、受压状态为主，部分围岩存在受拉的情况，图中的shear-p 和shear-n 分别表示岩体单元已经发生剪切破坏和正在发生剪切破坏，而tension-p 和tension-n 分别表示岩体单元已经发生拉伸破坏和正在发生拉伸破坏。

图9 静压力作用下未支护围岩受力状态分布云图Fig.9 Nephogram of stress state distribution of the surrounding rock without support under the static pressure load

图10和图11分别为隧道开挖后未进行支护和开挖后采取锚杆衬砌进行支护时，对隧道围岩开挖特征点拱顶（A点）的竖向位移曲线和拱腰（C点）的水平位移数据的监测曲线。

图10 静压力作用下拱顶（点A）竖向位移曲线Pig.10 Vertical displacement curves of vault (point A)under the static pressure load

图11 静压力作用下拱腰（点C）水平位移曲线Pig.11 Horizontal displacement curves of spandrel (point C)under the static pressure load

由图10和图11可以得知，在开挖计算前期，两个监测点的未支护与支护两组位移数值以相同的速率增加，随着迭代计算的进行，未支护组的围岩位移进一步增大，而采取锚杆衬砌支护组的围岩位移增加率小于未支护组的。

表4给出了隧道开挖后未进行支护和开挖后采取锚杆衬砌进行支护时，在静压力作用下A、C两点的位移量。

表4 静压力作用下A、C 点的位移量Table 4 Displacement of points A and C under the static pressure load cm

由表4中的数据可以得知，在采用锚杆衬砌支护后，拱顶（A点）的竖向位移和拱腰（C点）的水平位移较未支护时分别约减小了20.10%,27.11%。这一结果与围岩整体位移量的变化率基本一致，进一步印证了围岩整体位移与局部位移量的变化关系。

2.2.2 裂隙水压力作用下隧道开挖后的响应分析

由于隧道开挖地区中裂隙水压力的存在，当进行隧道模拟开挖时，需要对围岩应力与裂隙水渗流进行耦合，这会使得隧道开挖后围岩的位移量显著增加，模拟所得裂隙水压力作用下的隧道围岩位移量如表5所示。

表5 裂隙水压力作用下隧道围岩的位移量Table 5 Displacement of tunnel surrounding rock under the fracture water pressure cm

由表5中的数据可以得知，由于受地下水压力的影响，隧道围岩的整体竖向位移与水平位移较未考虑裂隙水压力时增加了2～3 倍。在对开挖隧道围岩施加砂浆锚杆与衬砌支护后，围岩整体的竖向位移与水平位移分别约减小了27.77%,25.67%。这一结果说明，本文模拟所采用的支护手段在考虑裂隙水的围岩中同样适用。

图12所示为裂隙水压力作用下，未支护时围岩的最大主应力分布云图。由图可知，裂隙水压力作用下未支护时围岩的最大主应力主要分布于隧洞周围，大小为12.0 MPa，可见，最大主应力值较无裂隙水作用时大幅度增加。从图12中还可以看出，围岩的最大主应力在节理两侧也发生了突变。

图12 裂隙水压力作用下未支护围岩最大主应力分布云图Fig.12 Nephogram of maximum principal stress distribution of the surrounding rock without support under the fracture water pressure

图13所示为裂隙水压力作用下，未支护时围岩的受力状态分布云图。

图13 裂隙水压力作用下未支护围岩受力状态分布云图Fig.13 Nephogram of stress state distribution of the surrounding rock without support under the fracture water pressure

由图13可以得知，隧道围岩主要以受压、受剪状态为主，部分围岩也存在受拉状态，整体呈现为“蝴蝶”型的受力状态。

对比图9和图13可以得知，其与静压力作用下的受力状态基本一致。

图14和图15分别为在裂隙水压力作用下，对隧道围岩开挖特征点拱顶（A点）竖向位移和拱腰（C点）水平位移数据监测曲线。由图14和图15可以得知，在开挖计算前期，两特征点处的未支护和支护后的两组位移均以相同的速率增加；随着迭代计算的进行，未支护组的位移进一步增加，而采用锚杆衬砌支护组的位移增加速率相对较小，但也在持续增大。

图14 裂隙水作用下拱顶（点A）竖向位移曲线Fig.14 Vertical displacement curves of vault (point A)under fracture water pressure

图15 裂隙水作用下拱肩（点C）水平位移曲线Fig.15 Horizontal displacement curves of spandrel (point C)under fracture water pressure

表6给出了模拟所得裂隙水作用下A、C两测点的位移数值。

表6 裂隙水压力作用下A、C 点的位移量Table 6 Displacement of point A and C under fracture water pressure cm

由表6可以得知，采用锚杆衬砌支护手段后，拱顶（A点）的竖向位移和拱腰（C点）的水平位移较未支护时分别约减小了32.85%,24.19%，相对整体位移量变化率稍小。由于裂隙水压力的影响，A、C两点的位移量较无裂隙水作用时分别约增加了1.52倍、4.40 倍，说明围岩整体位移与局部位移的一致性，也可进一步说明裂隙水对围岩整体位移和局部位移有较大的影响。

3 结论

本研究以国内某山岭高速公路隧道为背景，通过3DEC 离散元软件，对山岭裂隙岩体隧道进行开挖与支护模拟，主要得到以下结论：

1）对于含有复杂节理裂隙的岩体地下工程，采用离散单元法能更有效地模拟工程的实际情况，但3DEC 软件本身也存在不足，本研究主要通过临空面单元细分、忽略微小节理等手段予以减小软件对模拟结果的影响。

2）利用静压力荷载作用、静压力荷载与裂隙水压力耦合作用进行模拟对比，发现后者围岩的整体位移和局部位移量都较前者增加了2～3 倍，进一步说明围岩裂隙水对隧道开挖具有较大的影响；本文所模拟的砂浆锚杆单元（cable）和混凝土衬砌单元（liner）等支护手段能减小20%～32%的围岩位移量。

本文的研究成果可为今后研究裂隙水压力对隧道围岩开挖支护提供一定的参考。