Research on the Mechanism and Safe Thickness of Karst Tunnel-Induced Water Inrush under the Coupling Action of Blasting Load and Water Pressure

Abstract

When the drilling and blasting method is used to construct a tunnel through the karst stratum, the coupling effect of the blasting load and the karst water pressure in front of the tunnel face exposes the tunnel face to the risk of water inrush, which threatens the safety of personnel and property. It is very important for the design and construction of related tunnels to study the evolution mechanism of water inrush in karst tunnels and determine the minimum thickness of outburst prevention under blasting. Relying on the Dejiang tunnel Project in Tongren City, this paper adopts the Smoothed Particle Hydrodynamic–Finite Element Method (SPH-FEM) coupling calculation method to study the evolution process of water inrush in karst tunnels under blasting, analyzing the results of water inrush in tunnels under different rock wall thicknesses under blasting. Then, according to the regression of rock wall stress peak data, the analysis determines the minimum outburst prevention thickness of the karst tunnel. The research results show that there is a superposition effect between the blasting stress wave and the gravitational interaction of the karst water itself, and that the Smoothed Particle Hydrodynamic (SPH) particles in the aquifer cause damage and cracks to the rock wall under the coupling action of the blasting load and the karst water pressure, further leading to the expansion of the cracks and the formation of inrush channels. the stress, vibration velocity, and displacement of the unit at the junction of the aquifer and the rock layer show a trend that first decreases, then increases, and then decreases with an increase in the thickness of the rock wall. Based on the actual geological conditions of the Dejiang tunnel project parameters, when the thickness of the rock wall is 3.08 m, the peak stress of the rock formation unit at the junction with the aquifer reaches the maximum value. In order to avoid water inrush during blasting, the minimum outburst prevention thickness should be greater than 3.08 m. Based on the analysis results, a corresponding water inrush prevention plan was formulated on site which effectively guaranteed construction safety and, at the same time, verified the reliability of the analysis results. The relevant research results can provide useful references for similar projects.

1. Introduction

With the continuous development of traffic construction, many traffic tunnels inevitably pass through the karst strata. At present, drilling and blasting is still the main method of tunnel excavation. The coupling effect of blasting and karst water pressure leads to more cracks in the rock wall of the karst tunnel and, at the same time, will cause the cracks to expand and connect with each other, eventually forming a water inrush channel and resulting in a sudden water inrush accident [1,2,3]. Carrying out research on water inrush in karst tunnels under blasting and obtaining the minimum outburst prevention thickness will help to actively prevent water inrush disasters in karst tunnels.

Regarding previous research on water inrush in karst tunnels, most scholars have mainly carried out research from the main aspects of tunnel water inrush risk assessment, evolution mechanism, stability, and the safety thickness of water-resistant rock walls. In order to assess the risk of water inrush in tunnels, some scholars have proposed a mechanical model depending on geology, a mathematical model featuring the quantification of margin and uncertainty (QMU) [4], and two types of fuzzy comprehensive evaluation models with multiple factors affecting water inrush [5], based on the material element theory and ideal point method, to improve the risk assessment model of water inrush risk assessments in karst tunnels [6]. Other proposals have included the use of a geographic information system (GIS) for the dynamic prediction of a water inrush risk assessment model [7], a three-level, three-staged water inrush risk assessment model and water inrush interval risk assessment model [8], integrated standardization processes and an Analytic Hierarchy Process (AHP) inrush water risk assessment model [9], a karst water inrush (KWBF) conceptual model based on statistical analysis and phenomenological induction [10], a fuzzy comprehensive evaluation model of karst water inrush based on grey theory [11], and an evaluation model of the water inrush mechanism based on the hydromechanical coupling response behavior of geology and surrounding rocks [12]. In order to study the evolution mechanism of water inrush, some scholars have established model test systems for water inrush in high seepage and high in situ stress tunnels [13], a large-scale, three-dimensional model including two stages of tunnel excavation and hydraulic loading [14], a large buried depth and high water pressure tunnel water inrush model test system [15], a tunnel water inrush evaluation model suitable for fluid–structure coupling to test new materials [16], a multi-type water inrush model test system that can be serialized and visualized [17], and research on model tests for filling cracks such as the water inrush process test system [18] to explore the mechanism of water inrush in karst tunnels. At the same time, some scholars have established numerical models to analyze the inrush water caused by the gradual failure of the rock mass [19], the effect of cracks on the water flow [20], the gas–liquid two-phase flow mechanism after the water inrush in the tunnel [21], and the effect of seepage force on the expansion of the water inrush channel [22]. For the stability study of water-resistant rock walls, Li [23] proposed the upper bound theorem based on limit analysis and a Hoek–Brown failure criterion considering the influence of seepage force; Liu [24] conducted a real triaxial mechanical geological mechanical model test to evaluate rock wall stability; Wu [25] proposed a karst cave water pressure model to study the effect of cave water pressure on tunnel lining; Li, L. [26] established a simplified mechanical model to study the stability of rock wall caused by water-rich karst caves at the top of the tunnel sexual influence. Regarding the thickness of the water-resistant rock wall, some scholars have used the catastrophe theory and strength theory [27], the orthogonal numerical simulation test [28], combined the numerical simulation results and multiple linear regression [29], and the semi-analytical solution theory of the minimum safe thickness of the karst cave [30], minimum safe thickness calculation theory based on two-band theory and critical water pressure [31], minimum thickness calculation theory combining pseudo-static and dynamic theory [32], minimum safe thickness calculation theory considering energy dissipation [33], the critical safety thickness calculation theory based on the upper theorem of limit analysis [34], among other methods, to study rock wall safety thickness. The above analyses have mainly studied the mechanism of water inrush in karst tunnels under water pressure, while less research has involved the mechanism of water inrush in karst tunnels under the coupling action of blasting load and water pressure. At present, some scholars have established dynamic models to conduct theoretical analysis and numerical analysis on the mechanism of water inrush in karst tunnels under blasting loads [35,36,37,38]. However, there is a lack of scholars to conduct a new numerical analysis on the mechanism of water inrush in karst tunnels under the coupling action of blasting load and karst water pressure. Numerical analysis of the mechanism of water inrush in karst tunnels under blasting load requires the realization of groundwater inrush deformation under blasting disturbance, which cannot be solved by conventional finite element technology. As a new discrete element and finite element coupling calculation method, the Smoothed Particle Hydrodynamic–Finite Element Method (SPH-FEM) technology can solve the large deformation analysis that finite element technology cannot handle [39,40,41,42]. J. Wang [43] adopted the SPH-FEM coupling calculation method, in which the evolution process of water inrush under the blasting of karst tunnel is simulated by dividing Smoothed Particle Hydrodynamic (SPH) particles into groundwater. However, the above article did not study the minimum outburst prevention thickness of karst tunnels under blasting loads.

Relying on the karst section of the Dejiang tunnel in Tongren City, this paper adopts the SPH-FEM coupling calculation method to study the water inflow mechanism of the karst tunnel under the coupling action of blasting load and karst water pressure. Moreover, this study compares and analyzes the damage of rock wall under the independent action of blasting load and karst water pressure, analyzes the coupling effect of the blasting load and water pressure, analyzes the evolution results of water inflow in karst tunnels under different rock wall thicknesses, and performs regression analysis on the stress data of the numerical model to determine the minimum outburst prevention thickness before comparing the results with the minimum outburst prevention thickness obtained from the theoretical calculation formula of the minimum outburst prevention thickness of a karst tunnel under the blasting load. This has been carried out in order to provide a new suggestion for prevention and control in karst tunnel water inrush defense which has great practical significance.

2. SPH-FEM Coupling Analysis Method

2.1. Introduction to SPH Method

As a discrete element calculation method, the SPH method has gradually been applied to large analysis areas such as blasting and penetration due to its ability to handle large deformation properties. The system state in SPH is represented by a set of particles that can be assigned material properties and interact within a range of smoothing or weighting functions. These particle motions are controlled by the SPH formulation of the physical partial differential equation characterizing fluid flow; moreover, the SPH-controlled equation of motion can be transformed into the NS equation by replacing the SPH approximation with the function and its derivative with respect to the NS equation.

$$\frac{D{\rho }_i}{Dt}={\displaystyle \sum _{j=1}^N{m}_j{v}_{ij}^{\beta }\frac{\partial W_{ij}}{\partial x_i^{\beta }}}$$

(1)

$$\frac{D{v}_i^{\alpha }}{Dt}=-{\displaystyle \sum _{j=1}^N{m}_j(\frac{{\sigma }_i^{\alpha \beta }}{{\rho }_i^2}+\frac{{\sigma }_j^{\alpha \beta }}{{\rho }_j^2}})\frac{\partial W_{ij}}{\partial x_i^{\beta }}+F_i$$

(2)

$$\frac{D{e}_i}{Dt}=\frac{1}{2}{\displaystyle \sum _{j=1}^N{m}_j\left(\frac{P_i}{{\rho }_i^2}+\frac{P_j}{{\rho }_j^2}\right)v_{ij}^{\beta }\frac{\partial W_{ij}}{\partial x_i^{\beta }}}+\frac{{\mu }_i}{2{\rho }_i}{\epsilon }_i^{\alpha \beta }{\epsilon }_j^{\alpha \beta }$$

(3)

In the formula: ρ_i, $v_i^{\alpha }$, ε_i are the density, velocity component, and internal energy of particle i, respectively; N denotes the number of particles; α and β denote the coordinate directions; $x_i^{\beta }$ is the coordinate component of x along β; P_j denotes the isotropic pressure; W_ij denotes the smooth function between particles i and j; m_j denotes the particle mass.

From Equations (1)–(3), it can be found that SPH can calculate the velocity and energy changes of the particles inside the system and that the calculated data between the particles can be connected and act on each other; moreover, the whole action process varies with time. Further, the mesh-free nature of the SPH particles allows for the direct treatment of large deformations, avoiding the problem of large mesh distortions in numerical calculations.

Specifically, when the SPH calculation method is applied to karst tunnel surge water under the coupled action of blasting load and water pressure, by dividing the groundwater in front of the tunnel palm face into SPH particles and applying gravity to the SPH particles, the groundwater in front of the palm face of the karst tunnel can be realized as a result of the blasting dynamic action and its own gravity action to produce a water inrush deformation along the tunnel interior, as shown schematically in Figure 1. Compared to the finite element, it is more intuitive to represent the evolution of the water inrushing while avoiding large deformation distortions in the mesh.

Although SPH has the advantage of simulating large deformation, there is a serious boundary defect in the far field of the computational model due to particle inconsistency. In order to overcome this defect, a coupled analysis method using the SPH method in the near field and the Finite Element Method (FEM) in the far field is developed and applied, which combines the advantages of SPH (namely, being able to simulate large deformation) and FEM boundary calculation (being more stable and able to effectively overcome the problem of large deformation distortion and boundary instability).

2.2. SPH and FEM Contact Settings

The key to coupling SPH with FEM is that the intersection can effectively transfer SPH particle information to the finite element mesh. The common method is to set up interfacial contacts to bond the SPH particles to the finite element mesh to achieve displacement coordination between the two methods in a point–surface glued contact. The contact parameters are adjusted by setting the keyword *CONTROL_CONTACT. The contact schematic is shown in Figure 2.

3. Project Background

The Dejiang tunnel referred to in this study is located in Tongren City, Guizhou Province. It is a separate super-long tunnel with a total length of 5500 m. The surrounding rocks are moderately weathered argillaceous siltstone, limestone, mudstone, and mudstone intercalated with argillaceous limestone. The rock mass is developed with joints and fissures, and the rock mass is relatively broken, showing a mosaic-fragmented structure. The Dejiang tunnel approximately traverses the Wenxin fault, the Changfeng anticline, and the Shixiang line in an approximately orthogonal manner, as shown in Figure 3. Among them, an inclined karst catchment spanning 1050 m across the rock-facing area has been constructed for the water-rich karst section. According to the hydrogeological survey, the catchment area is 29.6 km². At the same time, because the rock faces the oblique catchment and the bottom limestone blocks the bottom water, a huge underground river system is formed above the tunnel.

The YK9 + 896~YK968 karst section is located in the core of the rock-facing slope, and the surrounding rock has developed fissures and is relatively fragmented, mostly with argillaceous deposits. Affected by the syncline catchment, there is high-pressure karst water. During the construction process, phenomena such as water inrush, mud inrush, and landslides occurred, as shown in Figure 4.

4. SPH-FEM Model

The section near the water-rich section of YK942~968 was used as the study object to establish the SPH-FEM numerical analysis model, and the cross and longitudinal section of the water-rich section is shown in Figure 5.

The model is a single layer grid model with a length of 55 m and a height of 60 m. The explosive is located in the center of the upper step of the tunnel and is 4 cm in diameter and 1.5 m in length. The aquifer in front of the palm face is divided by SPH particles, with a thickness of 4 m and a head height of 34 m. The rest of the model is meshed with a finite element mesh. Normal constraints and reflection-free boundaries are set at the bottom and left and right of the model, and the top of the model is free. Gravitational interaction is defined by the *LOAD_BODY keyword. SPH-FEM contact is achieved by setting *CONTROL_CONTACT to enable SPH-FEM interaction. The parameters contained in *CONTROL_CONTACT are shown in

Table 1. Contact Control Parameters.

Slsfac	Rwpnal	Islchk	Shlthk	Penopt	Thkchg	Orien
0.1	0.0	1.0	0.0	1.0	1.0	1.0
Enmass	Xpene	Rwgaps	Rwgdth	Rwksf	Icov	Swradf
0.0	4.0	1.0	0.0	1.0	0.0	0.0

.

The explosive material is simulated using *MAT_HIGH_EXPLOSIVE_BURN and the pressure–volume relationship during the explosive blast is described by the Jones–Wilkins–Lee (JWL) equation of state. The combustion fraction is multiplied by the high explosive’s equation of state to control the release of chemical energy to simulate the explosive detonation process. The air material is simulated using *MAT_NULL, while the pressure change process of the air is described using the EOS_LINEAR_POLYNOMIAL equation of state. The aqueous media material is also simulated using *MAT_NULL, and the pressure change process of the aqueous media material is represented by *EOS_GRUNEISEN. The initial support was simulated using the *MAT_HJC principal structure, which is a damage material model with applicability to concrete and rock at large deformations and high strain rates. The limestone is simulated using *MAT_HJC; it is similar to the widely used Johnson–Cook material model in metal materials, and its equivalent yield strength is a function of pressure, strain rate, and damage; moreover, the damage amount is a function of plastic body strain, equivalent plastic strain, and pressure. The constitutive relation of the JHC model is:

$${\sigma }^{\ast }=[A(1-D)+B{p}^{\ast N}][1+C\ln ({\overline{\epsilon }}^{\ast })]$$

(4)

In the formula, ${\sigma }^{*}$ is the normalized equivalent stress, ${\sigma }^{*}=\sigma /f_c^{\prime }$, where $\sigma$ is the actual equivalent stress, $f_c^{\prime }$ is the quasi-static compressive strength; A is the normalized viscosity coefficient; D is the damage parameter; B is the normalized pressure hardening coefficient; $P^{*}=P/f_c^{\prime }$ is the normalized pressure; ${\overline{\epsilon }}^{*}=\overline{\epsilon }/{\overline{\epsilon }}_0$ is the dimensionless strain rate; c is the strain rate coefficient.

Model cumulative damage consists of equivalent plastic strain and equivalent plastic volumetric strain, where the damage parameter D is:

$$D={\displaystyle \sum \frac{\Delta {\epsilon }_p+\Delta {\mu }_p}{D_1{(P^{*}+T^{*})}^{D_2}}}$$

(5)

In the formula, $\Delta {\epsilon }_p$ and $\Delta {\mu }_p$ are the equivalent plastic strain and equivalent plastic volumetric strain, respectively; $D_1$ and $D_2$ are the material constants; $T^{*}=T/f_c^{\prime }$ are the normalized maximum tensile hydrostatic pressures. When the stress on the element exceeds the defined tensile strength, it can be determined that the element is fractured. The parameters are shown in

Table 2. Explosives parameters.

ρ (kg/m3)	D (m/s)	A (GPa)	B (GPa)	R1	R2	w	E0 (GPa)	Pcj (GPa)
1.0	3600	229	0.55	6.5	1.0	0.35	3.2 × 109	14.4

,

Table 3. Air parameters.

ρ (kg/m3)	C0	C1	C2	C3	C4	C5	C6	E0	V0
1.29	0	0	0	0	0.4	0.4	0	2.50 × 105	1.0

,

Table 4. Water medium parameters.

ρ (kg/m3)	C (m/s)	S1	S2	S3	${\mathsf{\gamma }}_0$	a	e0	v0
1020	1650	2.56	−0.096	0	0.35	0	0	0

,

Table 5. Initial support parameters.

ρ (kg/m3)	E (GPa)	G (GPa)	A	B	C
2.37	28	19.8	0.79	1.6	0.007
N	T (GPa)	K1 (GPa)	K2 (GPa)	K3 (GPa)	v
0.63	0.004	86	−173	204	0.18

and

Table 6. Limestone parameters.

A	N	n	A1 (GPa)	A2 (GPa)	A3 (GPa)	B*	BQ	B1
2.50	0.85	3	4.25	4.12	4.12	1.6	0.0105	1.22
B2	T1 (GPa)	T2	pc (MPa)	pl (MPa)	D1	D2	M	Q2.0
1.22	4.387	0	13.3	60	1	0.001	0.85	0.685

.

In order to study the effect of rock wall thickness d on tunnel water inrush under blasting conditions, d was set to 0.6 m, 1.2 m, 1.8 m, 2.4 m, 3.0 m, 3.6 m, 4.2 m, and 4.8 m for the modeling analysis. The calculation model is shown in Figure 6. The calculated working conditions are shown in

Table 7. Calculation of working conditions.

Calculation of Working Conditions	1	2	3	4	5	6	7	8
d (m)	0.6	1.2	1.8	2.4	3.0	3.6	4.2	4.8

d—The width of the rock wall between the tunnel face and the aquifer.

.

5. Analysis of Numerical Model Results

5.1. Analysis of the Separate Effects of Blast Load and Water Pressure

Based on the above model parameters for a rock wall thickness of d = 1.2 m, separate models for blast load only (no aquifer) and water pressure only (no blast load) were developed to investigate the separate effects of the blast load and aquifer water pressure on the rock wall. As can be seen from Figure 7, the blast cracks produced by the blast load did not extend to the aquifer location when the original aquifer material was classified as a perimeter rock condition, indicating that the blast load alone did not directly cause a sudden surge of water in the tunnel, although the blast load alone caused cracking damage to the rock wall, the depth of which was d = 1.18 m. In addition, the stress cloud under separate water pressure indicates that, although the rock wall is not cracked under water pressure, there is a stress concentration at the shell hole with a pressure of P_S = 0.32 MPa, which will produce a stress release under blasting disturbance and is more likely to produce sudden water surges. During karst tunnel drilling and blasting construction, the blasting load is coupled with karst water pressure; therefore, it is necessary that the blasting load and karst water pressure are coupled and included in karst tunnel sudden surge water mechanism research.

Based on the above model parameters, for a rock wall thickness of d = 1.2 m, separate models for blast load only (no aquifer) and water pressure only (no blast load) were developed to investigate the separate effects of blast load and aquifer water pressure on the rock wall. As can be seen from Figure 7, in the original, aquifer material is dispersed into the surrounding rock conditions, and the blast cracking generated by the blast load on the rock wall caused blast damage with a depth of d = 1.18 m. While the width of the rock wall is 1.2 m, the damaged area does not extend to the aquifer location, which indicates that, although the blast load alone caused damage cracking on the rock wall, albeit in the absence of water pressure, it did not directly cause a water inrush in the tunnel. In addition, the stress cloud diagram under separate water pressure shows that the rock wall did not crack under water pressure, although the stress concentration at the bottom of the blast is obvious and there is a high stress area with a peak stress of P_S = 0.32 MPa. When the rock wall cracks due to the damage caused by blasting loads, the end of the resulting crack is close to reaching the stress concentration area at the bottom of the blast hole, in which case a stress release is likely to occur. This will likely result in a stress release, causing the damage fissure to penetrate the aquifer and cause water inrush in the tunnel. In the karst tunnel drilling and blasting construction process, the blasting load will be coupled with karst water pressure, and it is necessary to study the blasting load and karst water pressure coupling under the karst tunnel surge water mechanism.

Note the stress curve at measurement point #1 under the blast load and water pressure alone (see Figure 6 for an illustration of the measurement point selection). As can be seen from Figure 8a, the peak stress at measurement point #1 reached 51.69 MPa under the blasting load; however, the stress at the measurement point decayed rapidly with time. As can be seen from Figure 8b, the stress at measurement point #1 gradually increases under water pressure; from t = 5000 μs onwards, the stress value accelerates to P_S = 0.32 MPa, and then the stress value starts to decrease sharply to P_S = −0.45 MPa. After that, the peak stress shows an increasing trend again. The stress curve analysis shows that the blast load and water pressure have a sequential effect on the damage to the wall, with the blast load first causing damage fissures in the wall and water pressure acting as the dominant effect on the wall after the blast stress wave has decayed, leading to the expansion of the damage fissures and the formation of a water inrush channel.

5.2. Analysis of Sudden Water Surges in Tunnels under Coupled Blast Load and Water Pressure

Taking the calculation model of rock wall thickness d = 1.2 m as an example to analyze the evolution process of tunnel water inrush under blasting conditions, it can be seen from Figure 9 that, when t = 9 ms, the blasting stress wave is transmitted to the aquifer, and the SPH particles on the near-explosion side in the aquifer are squeezed against each other due to the impact of the stress wave. The SPH particles are distributed in a waveform. When t = 360 ms, under the action of gravity and the blasting stress wave, the SPH particles in the aquifer begin to move to the tunnel along the blasting crack and squeeze the surrounding rock of the tunnel face, accelerating the development of the cracks in the surrounding rock of the tunnel face, and consequently expanding the tunnel’s water-gushing channel. From the calculation diagram of t = 810 ms~1800 ms, it can be seen that, as the SPH water particles surge along the interior of the tunnel under the action of gravity, the cracks on the tunnel face evolve from the initial burst cracks to the thick water-gushing channels, and the SPH particles in the aquifer are accelerated. The surge into the tunnel is embodied in the accelerated increase in the cumulative particle number inside the tunnel and the emergence of larger and larger voids at the top of the aquifer.

The changing trend of the vibration velocity curve of the second measuring point under the coupling action of blasting load and hydraulic pressure is analyzed. The selection of measurement point #2 is shown in Figure 6. It can be seen from Figure 10 that, during the initial stage, the vibration velocity reaches a peak value of 1.15 cm/s in the first 10 ms, and the vibration velocity during the initial stage is mainly affected by the explosion load, and then decreases rapidly. However, the vibration velocity curve starts to increase rapidly from t = 500 ms, and the peak value increases continuously. When t = 2 s, the vibration velocity rapidly increases from 0.28 cm/s to 1.72 cm/s, which is mainly because the sudden change in the tunnel occurs at t = 500 ms. In the process of water inrush deformation, SPH water particles enter the tunnel along the crack channel, and the movement speed increases greatly. This process is reflected in the second half of the vibration velocity curve, and the vibration velocity of the second measurement point increases rapidly.

5.3. Analysis of Sudden Water Surges in Tunnels with Different Rock Wall Thicknesses

Figure 11 lists the cloud map of tunnel water gushing at t = 1800 ms with rock wall thickness d = 1.8 m, 2.4 m, 3.0 m, 3.6 m, and 4.2 m. It can be seen from Figure 11 that, when the thickness of the rock wall is d = 1.8 m, the SPH particles of the aquifer are still continuously pouring into the tunnel, and the influx quantity is reduced compared with that when the thickness of the rock wall is d = 0.6 m. When the thickness of the rock wall is 2.4 m, the influx of SPH water particles in the tunnel is still obvious; however, the number of SPH water particles is obviously reduced compared with that when the thickness of the rock wall is d = 1.8 m. At the same time, it can be observed that the rock wall is affected by the water pressure of the aquifer, resulting in extrusion fold deformation. The cloud map of when the thickness of the rock wall is d = 3.0 m shows that there is no influx of SPH particles into the tunnel at t = 1800 ms; however, the SPH particles in the aquifer are obviously deformed into the tunnel by the rock wall, creating conditions for a subsequent water influx. However, when the thickness of the rock wall is d = 3.6 m, the cloud map shows that, although the SPH particles of the aquifer squeeze the rock wall, the water inflow channel has not been formed, and only the aquifer particles are leaked along the bottom of the aquifer.

The SPH particles of different rock wall thickness models were extracted from adjacent SPH particles and finite element meshes on both sides of the rock wall boundary and at the position of the horizontal line of the blasthole for data analysis. The schematic diagram of measurement point extraction is shown in Figure 12a. It can be seen from Figure 12b that, when the thickness of the rock wall increases from 0.6 m to 4.2 m, the stress of measuring point #1 first decreases, then increases, and then decreases. The stress peaks of measuring point #1 corresponding to each rock wall width are 39.20 MPa, 18.32 MPa, 28.52 MPa, 16.31 MPa, 57.20 MPa, 8.67 MPa, and 9.82 MPa. When the thickness of the rock wall is 3.0 m, it reaches a peak value of 57.20 MPa and then begins to decrease. It can be seen from Figure 12c that the combined vibration velocity of the particles at measurement point #2 also shows a trend of first decreasing, then rising, and then decreasing by 1.26 cm/s, 1.00 cm/s, 2.01 cm/s, 0.48 cm/s, 0.47 cm/s. When the thickness reaches 3.0 m, it reaches a peak value of 1.73 cm/s, and then the peak vibration velocity begins to decrease with an increase in rock wall width. It can be seen from Figure 12d that the variation trend of the displacement map of the SPH particles at measuring point #2 is similar to the combined vibration velocity. For 0.065 m, 0.080 m, 0.029 m, 0.033 m, the peak displacement of 0.08 cm was also reached when the rock wall thickness was 3.0 m, and then the peak displacement began to decrease. From the perspective of stress peak, vibration velocity peak, and displacement peak, the minimum outburst prevention thickness of the rock wall of the karst tunnel under blasting should be greater than the thickness corresponding to the three peaks.

6. Data Analysis of Minimum Outburst Prevention Thickness for Karst Tunnels

After extracting the stress peak value of the #1 measuring point node of the rock wall thickness d = 0.6~4.8 m model for analysis and performing data fitting on the stress peak value, as shown in Figure 13, it can be seen from Figure 13 that the element stress peak shows a trend of first decreasing, then increasing, and then decreasing. The reason for this trend is that, when the thickness of the rock wall is small, the blasting stress wave and the water pressure of the aquifer cause a certain degree of offset, resulting in a decrease in the peak stress of the unit. Then, when the thickness of the rock wall continues to increase, the impact of the blasting stress wave on measuring point #1 is relatively weakened; moreover, since the position of the aquifer relative to measuring point #1 remains unchanged, its influence on the stress peak value of the unit can be regarded as unchanged; therefore, the stress peak value of measuring point #1 increases when the thickness of the rock wall increases and when the rock wall increases. When the distance reaches a certain distance, the disturbance to the aquifer by blasting becomes smaller, and no obvious stress redistribution is caused. At this time, the stress peak value of measuring point #1 decreases.

Data fitting on the stress peak value was performed to obtain a fitting curve. The fitting formula used was:

$$y=(-19.38x+101.74)\cos {(x+2.97)}^2$$

(6)

After using the derivative of the fitted curve to obtain the derivative change trend curve of the stress peak curve and finding the zero point of the stress peak derivative curve to obtain x = 3.08 m, that is, when the thickness of the rock wall is 3.08 m, the unit stress peak reached the maximum. At this time, the risk of water inrush caused by blasting is the greatest, and the thickness of the rock wall should be at least greater than this thickness during construction.

7. Theoretical Verification of Minimum Overburden Thickness

The minimum thickness L of a karst tunnel under the blasting load consists of two parts: the blasting damage area thickness L_c and the hydraulic fracture area thickness L_w. During blasting construction, the stress waves generally disturb the front of the palm face to a depth of no more than 1.5 m, as shown in Figure 14. In this paper, Because L_c represents the initial crack before the action of water pressure, it is considered as the crack length generated by blasting when L_c is taken, and this data has been analyzed and obtained in Section 5.1; therefore, L_c = 1.18 m.

The minimum outburst prevention thickness L of a karst tunnel under the action of the blasting load consists of two parts, namely, the thickness of the blasting damage area L_c and the thickness of the hydraulic fracturing area L_w. During blasting construction, the disturbance depth of the stress wave to the front of the tunnel face generally does not exceed 1.5 m. The calculation diagram is shown in Figure 14. In this paper, L_c is set as 1.18 m according to the actual blasting damage thickness in the calculation.

When P_w < P_lj, the thickness L_w of the hydraulically fractured zone is calculated as:

$$\left.\left\{\ln \lambda -\ln \left[\lambda -\frac{2f{P}_{\mathrm{w}}\sqrt{\pi a}-\sqrt{3}{K}_{\mathrm{II}C}-2c\sqrt{\pi a}+2{\tau }_1\left|K_{\mathrm{II}}^{(2)}\right|\sqrt{\pi a}}{\gamma H\sqrt{\pi a}(f-f\cos 2\beta +\sin 2\beta )}\right]\right.-\frac{f+f\cos 2\beta -\sin 2\beta }{f-f\cos 2\beta +\sin 2\beta }\right\}$$

(7)

In the formula, R is the radius of the tunnel equivalence circle, R = (B + h)/4, B is the span of the tunnel, h is the height of the upper step of the tunnel excavation, λ is the lateral pressure coefficient, f is the friction coefficient, a is one-half of the crack length, K_ⅡC is the type II compression–shear extension stress intensity factor, c is the cohesion of the crack surface, ${\tau }_1$ for the maximum shear stress on the wavefront surface, $K_{\mathrm{II}}^{(2)}$ for the dynamic stress intensity factor, γ is the surrounding rock weight, H is the tunnel depth, and β is the vertical direction to the center of the fissure to establish the angle of the coordinate system X-axis.

The buried depth of Dejiang tunnel at the YK9 + 915 section is 312 m, the tunnel section span B is 10.8 m, and the height h is 5.6 m. There is an aquifer in this area. The head height of the tunnel roof is 34 m, the water pressure P_w is 0.32 MPa, and the limestone type II pressure. The shear expansion stress intensity factor is selected according to the literature, which is taken as 7.94 $\times {10}^{6}N/{\mathrm{m}}^{3/2}$, the lateral pressure coefficient λ is taken as 1.5, the crack length a is taken as the initial crack thickness of 1.18 m, the friction coefficient of the crack surface is taken as 0.7, and the cohesion coefficient is taken as 0.7. The force c is taken as 0.3 MPa, the dynamic stress intensity factor $K_{\mathrm{II}}^{(2)}$ is taken as 0.4, and β is taken as 12°.

According to the parameters of the Dejiang tunnel project, the thickness L_w of the hydraulic fracturing zone under the blasting load can be obtained to be 1.60 m; then, the minimum outburst prevention thickness of the karst tunnel under the blasting load calculated by the formula is L = 2.78 m, and the numerical model L determined by the mathematical analysis of the data is 3.08 m, and the error between the two is 9.7%. This proves that the minimum outburst prevention thickness of a karst tunnel under blasting action obtained by numerical model analysis is reliable.

8. On-Site Remediation Measures

Manuscripts reporting large datasets that are deposited in a public YK9 + 928-YK9 + 968 include the section with developed lava fissures and water-rich sections of the Dejiang tunnel, which has been proven to have groundwater through advanced drilling. During the construction period, advanced grouting reinforcement measures were carried out for this section. The grouting range was 5 m outside the excavation outline, the grouting length of each cycle was 30 m, and the excavation was 26 m. According to the previous analysis, the minimum anti-outburst thickness was 3.08 m. For the sake of safety, a 4 m stop magma wall is reserved. The grouting holes are arranged according to a slurry diffusion radius of 2 m and a hole bottom spacing of 3 m. There are four rings of grouting holes in each cycle. The diameter of the grouting hole is not less than 90 mm. On-site diagnosis and treatment are shown in Figure 15. After grouting reinforcement, excavation was carried out by means of short footage and weak blasting. Drill-and-blast excavation was carried out under the condition that the rock wall thickness was strictly controlled to be 4 m, and grouting prereinforcement was carried out before the drill-and-blast excavation. In the end, the construction of the water-rich section did not have any accident of sudden water gushing caused by blasting, which ensured the tunnel’s safety during drilling and blasting construction.

9. Conclusions

In this paper, based on the Dejiang tunnel in Tongren City, Guizhou Province, the SPH-FEM coupling calculation method is used to analyze the water inrush conditions of the tunnel under different rock wall thicknesses, and the minimum outburst prevention thickness of the karst tunnel under blasting is determined according to the analysis of the unit stress data of the calculation model. The main conclusions are as follows:

(1)

Compared with the independent effects of the blasting load and water pressure, the coupling effect of the blasting stress wave and water pressure will increase the depth of damage to the tunnel rock wall and increase the risk of water inrush inside the tunnel;

(2)

Under the coupling action of the blasting load and water pressure, as the thickness of the rock wall increases within a certain range, the stress, vibration velocity, and displacement of the unit at the junction of the rock wall and the aquifer show a variation law that first decreases, then increases, and then decreases, and there is peak stress, peak vibration velocity, and peak displacement. In order to reduce the risk of water inrush during tunnel construction, the minimum thickness of the rock wall for outburst prevention should be at least greater than the width of the rock wall corresponding to the peak stress, peak vibration velocity, and peak displacement;

(3)

Comparing the minimum outburst prevention thickness analyzed based on the actual geological parameter modeling with the existing theoretical calculation formula for the minimum outburst prevention thickness of karst tunnels under drilling and blasting conditions, the error between the two is only 9.7%. The minimum anti-outburst thickness guides the construction and ensures the construction safety of the supporting project.