基于UDEC的岩体完整性系数对岩爆烈度影响的离散元分析

Discrete Element Analysis on the Influence of Rock Mass Integrity Coefficient on Rockburst Intensity Based on UDEC

摘要: 针对多数学者在岩爆预测指标选取中并未考虑岩体完整性系数K_V对岩爆烈度预测准确率的作用，为进一步明确岩爆发生时岩体完整性系数K_V对岩爆烈度的影响，本文以江边水电站引水隧洞为工程案例，基于弹性应变能理论，利用UDEC离散元软件建立隧洞开挖数值模型，通过在岩体中设置结构面间距分别为0.5 m、0.75 m、1 m、2 m的两组平行相交结构面，模拟隧洞开挖过程中岩体的力学行为特征，并分析围岩岩爆破坏时的弹性应变能云图、应力云图、位移云图以及弹射岩块速度云图。结果表明，随着结构面间距的增大，即岩体完整性的提高，岩体发生破坏时所需的能量越高，产生的应力、位移、弹射速度也越大，即岩爆烈度越大。

Abstract: Addressing the fact that most scholars have not considered the role of rock mass integrity coefficient (K_V) in improving the accuracy of rockburst intensity prediction in selecting predictive indicators, this study aims to further clarify the influence of rock mass integrity coefficient (K_V) on rockburst intensity during rockburst events. Taking the water diversion tunnel of Jiangbian Hydropower Station as a case study, based on the theory of elastic strain energy, a numerical model of tunnel excavation is established using the UDEC discrete element software. Two sets of parallel intersecting structural planes with spacings of 0.5 m, 1 m, 1.5 m, and 2 m are set within the rock mass to simulate the mechanical behavior of the rock mass during tunnel excavation. Subsequently, the elastic strain energy contours, stress contours, displacement contours, and ejection velocity contours of rock fragments during surrounding rockburst failure are analyzed. The results indicate that with the increase in the spacing between structural planes, which signifies an improvement in rock mass integrity, the energy required for rock failure increases, leading to higher stresses, displacements, and ejection velocities. In other words, the rockburst intensity intensifies.

1. 引言

岩爆是指地下岩体在经过人为开挖扰动或其他外加扰动后发生爆裂、松脱、抛射的一种动力性破坏现象，其现象的发生与岩体中应变能的集聚、突然释放密切相关。由于浅部资源日渐枯竭以及地下工程建设需要，人类开采矿物、开挖地下硐室逐渐从地表浅层转向地壳深层[1]。伴随着埋深与地应力的增加，岩体地质环境及扰动后的岩体力学行为更加复杂，地下工程建设中面临着越来越多的岩爆灾害问题，已成为限制大型水电、交通、矿山及国防地下工程进一步发展的难点问题[2]。

岩爆预测对岩爆防控至关重要，准确预测岩爆烈度等级可以帮助选择安全的工程选址和采取合适的施工方法，这能有效降低岩爆对工程造成的风险，并保障施工人员和设施的安全。近些年来，岩爆预测研究进一步深入，基于机器学习的多因素岩爆预测研究越来越受学者的青睐[3]-[5]。其中，岩石单轴抗压强度(σ_c)、岩石抗拉强度(σ_t)、岩石脆性系数(σ_c/σ_t)、最大切向应力(σ_θ)、围岩应力系数(σ_θ/σ_c)、岩石弹性变形能指数(W_et)等指标是学者在进行岩爆预测研究中常选取的评价指标[6]-[8]。吴川顺[9]等选取六指标(σ_θ, σ_c, σ_t, σ_θ/σ_c, σ_c/σ_t, W_et)并提出主成分分析法、概率神经网络对岩爆烈度进行分级预测；谢学斌[10]等基于四指标(σ_θ, σ_θ/σ_c, σ_c/σ_t, W_et)并采用改进CRITIC法预测岩爆等级；胡建华[11]等提出由三指标(σ_θ/σ_c, σ_c/σ_t, W_et)组成的评估体系，建立了岩爆倾向性有限元评价模型；齐石[12]等选择了四指标(σ_c/σ_t, σ_θ, σ_θ/σ_c, W_et)并基于L-XGB算法预测岩爆等级。上述学者在选择岩爆预测指标时，都并未考虑岩体完整性系数K_V对岩爆烈度的影响。

岩体完整性系数K_V是体现岩体稳定性的重要指标之一，一项[13]利用贝叶斯网络并基于429组岩爆案例建立的岩爆概率预测模型发现，岩体完整性系数对岩爆影响很大。张亭等[14]选取三指标(σ_θ /σ_c, σ_c/σ_t, W_et)建立BPM神经网络模型预测岩爆，准确率达84.78%；杨腾杰等[15]提出了一种基于改进TOPSIS-GRA的岩爆倾向性预测模型，该模型考虑了岩体完整性系数K_V，通过选取四指标(σ_θ/σ_c, σ_c/σ_t, W_et、K_V)构建评价体系，并在17组岩爆实例中验证，准确率达88.9%；赵国彦等[16]为高效预测岩爆烈度，结合主成分分析(PCA)与最优路径森林(OPF)算法，选取5项关键指标(σ_c, σ_θ/σ_c, σ_c/σ_t, W_et, K_V)构建PCA-OPF模型，该模型预测准确率达91.25%，优于其他模型，表明包含岩体完整性系数Kv的评价体系对岩爆预测具有可靠性和适用性。

为此，为了表明岩体完整性系数K_V对岩爆烈度的影响，本研究基于离散元数值模拟和弹性应变能理论，以江边水电站引水隧洞工程为工程案例，利用UDEC离散元软件建立隧洞开挖数值模型，设置两组间距不同的平行相交结构面，模拟岩体完整性系数K_V对岩爆的影响。模拟结果表明，随着岩体完整性的提高，岩爆烈度也越大，表明岩体完整性系数对岩爆有重要影响。

2. 工程背景与地质资料

江边水电站[17]位于踏卡背斜的西南翼，距离山层断裂处约918公里。该水电站地下结构穿越的地层为二叠系地层，地层产状结构相对稳定，结构面方向一般为N55˚~70˚W，N55˚~65˚E范围内。根据地质测绘，该水电站地质构造比较简单，没有大区域分布范围的断裂或者大规模的断层分布。其中该水电站场地的地质构造主要以小规模断层和结构面为主，主要地质岩性包括燕山期黑云母花岗岩(yB52)与二叠系甲黄沟群地层接触带的蚀变带。燕山期黑云母花岗岩地层区结构面较为发育，分为2组：N5˚~10˚W，N80˚~85˚E，结构面间距中充填片状岩层和石英石且结构面弯曲粗糙；N15˚~40˚E，S60˚~85˚E结构面间距充填岩块和石英石且结构面大多表面粗糙，产状弯曲延伸。根据现场施工的勘察记录，该水电站在引水隧洞施工中岩爆一般发生在两组平行相间结构面的围岩中，且在N80˚~85˚E和S60˚~85˚E范围组合的结构面围岩中发生最频繁，且发生岩爆烈度的大小不一样。因此，本文隧洞模型围岩设置结构面倾角为80˚~85˚和60˚~85˚，结构间距分别为2 m、1 m、0.75 m、0.5 m。

3. 弹性应变能理论

从能量的观点出发，岩爆发生时，岩体中的能量会释放出来。其原理是开挖造成的临空面会加剧围岩损伤，导致岩体内应变能累积，一旦此能量突破岩体承受极限，即发生岩爆。因此弹性应变能也是岩爆分级预测的重要指标，探讨不同岩体完整性下岩爆的起裂、裂纹扩展及岩块抛掷机制，从能量角度分析具有实际意义。在复杂应力状态下，岩体单元的能量满足下列关系：

$U=U^d+U^e$ (1)

岩体单元各部分可表示为：

$U={\displaystyle {\int }_{10}^{\epsilon }{\sigma }_1}\text{d}{\epsilon }_1+{\displaystyle {\int }_{10}^{\epsilon }{\sigma }_2}\text{d}{\epsilon }_2+{\displaystyle {\int }_{10}^{\epsilon }{\sigma }_3}\text{d}{\epsilon }_3$ (2)

$U^e={\sigma }_1{\epsilon }_1^e/2+\frac{1}{2}{\sigma }_2{\epsilon }_2^e\frac{1}{2}{\sigma }_3{\epsilon }_3^e$ (3)

${\epsilon }_i^e=\frac{1}{E_i}\left[{\sigma }_i-v_i\left({\sigma }_j+{\sigma }_k\right)\right]$ (4)

式中：U是外力对围岩所产生的总输入能量，$U^d$ 岩体单元的耗散能，$U^e$ 是单元能释放出的弹性应变能。${\epsilon }_i$ 是主应力${\sigma }_i$ 方向上的总应变，${\epsilon }_i^e$ 是弹性应变，$E_i$ 是卸载弹性模量，$V_i$ 是卸荷泊松比。

为适合开挖工况的岩体单元可释放应变能的统计，取初始的弹性模量$E_0$ 与泊松比v进行简化计算：

$U^e=\frac{1}{2E_o}\left[{\sigma }_1^2+{\sigma }_2^2+{\sigma }_3^2-2v\left({\sigma }_1{\sigma }_2+{\sigma }_1{\sigma }_3+{\sigma }_2{\sigma }_3\right)\right]$ (5)

式中：$U^e$ 是隧道围岩中原始应力场中积聚的弹性应变能，${\sigma }_1$ 、${\sigma }_2$ 、${\sigma }_3$ 分别是原始应力场中的3个方向的主应力。

4. 数值模拟方案

4.1. UDEC离散元软件介绍

UDEC程序被许多学者和工程技术人员广泛地应用于岩土体力学行为相关问题的研究[18]。岩石在自然状态下经历了各种构造运动和风化作用，形成了复杂的非均质结构，而离散元法将岩体视为由多个刚性块体组成的集合，每个块体之间通过接触面相互作用。与其他模拟软件不同，UDEC将岩石材料分割为块体，并将结构面视为块体的边界，模拟过程中，非连续的结构面必须贯穿整个岩体，否则将被自动删除。因此，虚拟结构面经常被引入以保持连通性，并根据实际情况进行调整。

4.2. 参数设置

在用离散元软件(UDEC)做隧道开挖模拟时，存在超级大岩土体难以进行模拟和处理问题，因此在设置开挖边界必须要考虑人工边界对于数值开挖模型计算结果的影响[19]。

1) 岩体破坏采用摩尔库伦破坏准则，并采用弹塑性本构模型。结构面本构采用摩尔库伦接触模型，模型上边界和左右边界为应力约束，下边界为固定位移边界。

2) 数值模拟的初始地应力参数根据江边水电站勘察资料[20]，初始垂直应力和水平应力设置为15 MPa和20 MPa，根据工况设置隧道开挖洞径半圆半径7 m，直墙为8 m，取模型计算尺寸为50 m × 50 m，详见图1。

Figure 1. Schematic diagram of numerical simulation scheme

图1. 数值模拟方案示意图

3) 设置两组平行相交的结构面，并设置四种不同间距平行相交结构面的分布工况，分别为2 m，1 m，0.75 m，0.5 m。四种工况下两组平行结构面的倾向与相交角度均保持一致，各倾向的角度与江边水电站统计的中结构面分布形式大致一致。

4) 先对模型施加初始地应力场，当计算达到应力平衡状态后，再进行开挖。

5) 根据统计江边水电站的工程地质勘察报告以及张德永对江边水电站地质参数的统计[17]，得到江水电站工程区岩石参数建议值，见表1和表2。

Table 1. Table of physical and mechanical parameters of surrounding rock

表1. 围岩的物理力学参数表

岩性	容重/t/m^3	弹性模型/GPa	泊松比	粘聚力/MPa	摩擦角/˚	抗拉强度/MPa
黑云母花岩岗	26.4	24.5	0.25	1.7	46	3.84

Table 2. Table of physical and mechanical parameters of structural plane

表2. 结构面物理力学参数表

结构面类型	倾角/˚	法向刚度/GPa	切向刚度/GPa	内摩擦角/˚	粘聚力/MPa	抗拉强度/MPa
结构面	0~90	20	6.6	25	0.025	0

4.3. 后处理

基于离散元软件UDEC建立不同结构间距分布的隧洞开挖模型，计算开挖平衡后的竖向应力、破坏位移、弹射岩块速度云图，分析不同岩体完整性条件下高地应力隧洞开挖的岩爆范围、力学行为特征。具体数值模拟方案示意图如图1所示。

5. 数值模拟结果分析

5.1. 弹性应变能分析

弹性应变能累积能量是从每个单元提取对应的变量，根据提取的变量，由公式(5)计算每个单元的弹性应变能密度，最后用fish函数的zone extra参数把累积结果计算出来。四种工况在隧道开挖后弹性应变能的释放与累积情况见表3。由表3可得，释放的能量数值随着节理间距的增大大致呈现增大的趋势，累积的能量数值呈现增大的趋势；不同间距节理岩体围岩中能量释放的单元面积与累计能量的单元面积皆大致相同。

Table 3. Release value and area of elastic strain energy

表3. 弹性应变能释放数值与面积

节理间距(m)	释放能量值	释放能量面积	累积能量值	累积能量面积(mm^2)
0.5	3.8e^6	3.6e^2	1.89e^7	1.95e^3
0.75	3.9e^6	3.45e^2	2.0e^7	1.96e^3
1	4.2e^6	3.43e^2	2.47e^7	1.96e^3
2	4.4e^6	3.42e^2	2.5e^7	1.96e^3

根据上述分析，隧洞开挖时，岩体的完整性程度基本不影响围岩的扰动范围，但影响围岩的储能性能、破坏区域大小与能量释放大小。结构面间距越大，即围岩完整性程度越好，岩体也能够承受更多的能量，同时也释放更多能量。

根据利用fish语言编写的弹性能计算公式，可以得到岩爆前后的各工况弹性能释放分布云图如图2所示。洞室的拱底和边墙有弹性应变能的释放，故可推断出该处发生了岩爆，随着结构面间距的减少，洞室围岩中释放出应变能的区域在增多，洞壁更多的地方发生岩爆破坏，可说明当间距越小时，岩体被切割成更小的岩块，小的岩块所能承受的临界能量更低，岩块更易达到它破坏临界能量条件。

(a) 0.5 m间距 (b) 0.75 m间距

(c) 1.0 m间距 (d) 2.0 m间距

Figure 2. Distribution cloud of elastic strain energy release after tunnel excavation with different joint spacing

图2. 不同间距节理岩体隧道开挖后的弹性应变能释放分布云图

图中绿色区域是弹性应变能积聚的区域，显然绿色区域多集中与两侧拱底和拱顶，因此当能量进一步积累或者受到扰动有发生岩爆的潜在可能。由于结构面的控制效应，沿着结构面的倾向积聚能量，并沿着结构面方向发生滑移型破坏，图中黑色线条区域为潜在的裂隙破坏区域。从图2可以看出，当间距较大时，潜在破坏区域集中在拱顶和拱底，当间距变小时，左右的拱肩和左边拱底也出现少部分潜在裂隙破坏区域，由此可得，结构面间距越大，岩爆发生的烈度越大。

5.2. 应力分析

隧道开挖后，会造成围岩应力的重分布并引起应力集中，当岩体承载的应力过大并超过围岩体的极限时，开挖后的围岩体便会发生岩爆。由图3可知，当结构面间距越大(0.5 m, 0.75 m, 1 m, 2 m)，应力集中程度越高，高应力的分布区域面积越大，同时岩爆发生时围岩体的最大竖向应力也越大，分别为29.6 MPa，30.3 MPa，31.1 MPa，32.4 MPa。这表明：结构面间距越大，岩体完整性越优，岩体所能承受破坏的极限强度也越大，符合表4关于岩体完整性系数的岩爆烈度分级标准。

(a) 0.5 m间距 (b) 0.75 m间距

(c) 1.0 m间距 (d) 2.0 m间距

Figure 3. Vertical stress cloud map of rock mass tunnel excavation with different spacing joints

图3. 不同间距节理岩体隧洞开挖后的竖向应力云图

Table 4. Rock mass integrity coefficient rock burst classification standard [21]

表4. 岩体完整性系数岩爆分级标准[21]

岩爆等级	无岩爆	弱岩爆	中等岩爆	强岩爆
Kv	<0.55	0.55~0.65	0.65~0.75	>0.75

5.3. 位移分析

根据隧道开挖过程中，显示出洞室围岩的位移破坏形式如图4所示，洞室开挖后围岩位移云图与图2的弹性应变能释放云图相对应，由于弹性应变能的释放，洞室的边墙发生明显的破坏位移，右边的位移比左边的位移更明显，围岩变形破坏主要在两侧的边墙，沿着结构面发生滑移破坏，因此产生的水平位移更大，大多数发生位移的岩块都是被结构面切割出来的，拱顶和拱底也有少量的碎屑弹出。随着节理间距的

(a) 0.5 m间距 (b) 0.75 m间距

(c) 1.0 m间距 (d) 2.0 m间距

Figure 4. Displacement cloud map of tunnel excavation of rock mass with different spacing joints

图4. 不同间距节理岩体隧洞开挖的位移云图

增大，变形破坏时发生的位移更大，被结构面切割的岩块更大，其破坏时能承受的临界能量越大，释放出来的弹性应变能越多，导致发生岩爆时岩块弹射出的位移更大。

5.4. 岩块弹射速度分析

图5为不同间距节理岩体隧洞开挖后发生岩爆时的岩块弹射速度云图，云图中显示被结构面切割的岩块在开挖过程中发生岩爆弹射现象，且随着节理间距的增大，发生岩爆弹射出来岩块的速度更大，同时质量也越大，其弹射的动能也会越大，预示着岩爆发生的烈度更大。

(a) 0.5 m间距 (b) 0.75 m间距

(c) 1.0 m间距 (d) 2.0 m间距

Figure 5. The cloud image of ejection velocity of damaged rock mass in surrounding rock after tunnel excavation with different spacing joints

图5. 不同间距节理岩体隧洞开挖后围岩破坏岩块的弹射速度云图

6. 结论

基于四种不同结构面间距岩体隧洞开挖过程的模拟分析，研究了岩体完整性程度对岩爆的影响规律。得到结论如下：

1) 随着结构面间距的增大，岩体完整性越强，岩体发生破坏时竖向应力越大，洞壁围岩应力集中程度越大。

2) 随着结构面间距的增大，岩体完整性越强，岩爆破坏产生的位移越大，岩爆越剧烈。由于结构面存在控制效应，多沿着结构面分割岩体的倾向发生位移。

3) 随着结构面间距的增大，岩体完整性越强，岩体发生破坏时，弹射出来的岩块速度更快，质量也越大，其弹射的动能也更大，岩爆烈度也越大。

4) 从弹性应变能的角度来说，结构面间距越大，围岩完整性程度越好，岩爆发生时岩体也能释放更多能量，岩爆烈度也越高。

综上所述，岩体完整性对岩爆范围与烈度有较强的控制性作用，表明了在进行岩爆预测时需要考虑岩体完整性系数对岩爆烈度的影响。