A Study on the Surface Vibration Effect of CO₂ Phase Transition Cracking Based on the Time-Domain Recursive Analysis Method

Abstract

Supercritical CO₂ phase change fracturing technology has been widely used in rock engineering. However, the surface vibration characteristics induced by phase change rock breaking have not been sufficiently studied. In this paper, a model for calculating the surface vibration velocity induced by supercritical CO₂ phase change was established based on the time-domain recursive method, and the reliability of the model was verified using LS-DYNA nonlinear finite element software; based on the model for calculating the blasting energy of a compressed gas and water vapor container, the surface vibration velocity and decay law induced by CO₂ phase change fracturing and equivalent explosive blasting were compared and analyzed. The results show that the ground vibration velocity calculation model based on the time-domain recurrence method can better evaluate the ground vibration characteristics induced by supercritical CO₂ phase change fracturing, and the vibration velocity decays exponentially with the fracture distance; the peak pressure of supercritical CO₂ phase change is only 1/3.36 of that of explosive blasting, but the action time is 100 times longer than that of explosive blasting, and the peak vibration velocity is 1/74~1/78 of that of equivalent explosive blasting.

1. Introduction

With the development of water transportation projects, the mainline of the Yangtze River needs to be de-reefed and broken to widen the scale of the channel, for which safe and efficient rock-breaking methods are needed. Traditional explosives have the advantages of economy and high efficiency [1,2]; however, their high risk and high pollution make them subject to strict controls and even bans [3,4]. In contrast, supercritical CO₂ phase change rock-breaking technology has received more and more attention in the field of rock-breaking due to the advantages of the spark-free fracturing process and pollution-free phase change products [5,6]. This technology is safer and more environmentally friendly than traditional explosive blasting, because the fracturing process is a physical one [7,8]; this is because the liquid CO₂ is instantly vaporized by heating the heat pipe, and the gas inside the pipe expands to form the impact pressure and splitting pressure to break the rock effectively.

The vibration generated by blasting Inevitably affects adjacent buildings and the surrounding environment at different levels, and blasting vibration effects have therefore become an important research topic [9]. In practical engineering, vibration velocity and main vibration frequency are two key indicators of the blasting vibration safety criterion, and some tests have shown that the vibration velocity induced by CO₂ phase change fracturing is only 20% of that of explosive blasting [10]. However, there are not many studies on the ground surface vibration law caused by CO₂ phase change fracturing, and a corresponding ground surface vibration velocity calculation model has not been established.

Vibration velocity calculation plays a key role in blast vibration control. Previous studies on vibration velocity are mainly based on field vibration monitoring, and the peak particle velocity (PPV) is obtained from the monitoring data and used to characterize the vibration velocity [11]. At present, this research is mainly focused on the field of explosive blasting, and there are few studies on the vibration velocity of supercritical CO₂ phase change fracturing. In order to obtain the PPV values, many empirical formulas for the vibration velocity due to explosive blasting have been proposed [12,13]; these are related to the explosive charge and geological conditions, and therefore have some limitations in use. With the rapid development of computer technology, the application of artificial intelligence has been extended to the field of vibration velocity prediction, and corresponding models such as the artificial neural network model, neuro-fuzzy model, and hidden Markov model have been applied to the study of blasting vibration velocity under different environments, one after another; certain application results have been achieved [14,15,16], and these have formed the direction of research concerning supercritical CO₂ blasting vibration velocity for the future. At the same time, stress wave propagation theory has also been gradually applied to the field of supercritical CO₂ blasting, and existing research methods of this theory mainly include the three basic conservation laws, the characteristic line method, and the displacement discontinuity method [17,18,19]. Li and Ma combined the momentum conservation theory and the displacement discontinuity method, and proposed the time-domain recursive analysis method (TDRM) [20], which avoids complex mathematical operations such as Fourier transform and effectively improves computational efficiency. Scholars have extended this method to analyze the propagation law of different types of stress waves in different rock masses, and fruitful research results have been achieved [21,22]. These reflect the effectiveness of this method, and also provide a good direction for future studies of supercritical CO₂ phase change blast vibration velocity.

As mentioned earlier, current methods for calculating the ground vibration velocity mostly use field measurements and empirical formulas, which have certain limitations. In this paper, from the essential mechanism of ground vibration induced by supercritical CO₂ phase change, we describe the process of supercritical CO₂ blasting stress waves on rock bodies, based on stress wave theory and the time-domain recurrence method, and establish a model to calculate the ground vibration velocity after supercritical CO₂ blasting, considering time change. Subsequently, LS-DYNA software was used to simulate the vibration velocity of supercritical CO₂ blasting, and the calculation results were compared with the calculation formula of surface vibration velocity, as well as with the equivalent explosive blasting. The results of the study may provide a new calculation method for the calculation of surface vibration velocity induced by supercritical CO₂ blasting.

2. Law of Surface Reflection Based on the Time-Domain Recursive Method

Before supercritical CO₂ blasting, the fission tube is typically placed inside a pre-drilled hole, and then the hole is sealed for blasting. Figure 1 illustrates the arrangement of the fission tube, with a length of h₂. The ground level is situated h₁ above the fission tube, while the front represents the critical surface, and the remaining direction corresponds to the infinite rock body. Upon igniting the fission tube for supercritical CO₂ blasting, a stress wave, specifically a P-wave, is generated in the rock body due to the impact pressure of the gas. The propagation process of the P-wave produces various effects on the ground surface, forming the basis for studying the velocity of the blasting vibration.

2.1. A Generalized Model of the Transmissivity Reflection of the Stress Wave on a Free Surface

Assuming that blasting induces a stress wave in the rock body, this wave propagates as a column surface wave. The fracturing tube’s blasting location is considered a point, from which the stress wave (referred to as the P wave) is emitted to the surrounding area. The P wave can approach the air and ground surfaces at arbitrary angles, and at each point of incidence, it undergoes reflection and generates corresponding transmitted P and S waves, as well as reflected P and S waves.

The distribution of the reflected and incident waves is symmetrical concerning the vertical line of the incident surface, as illustrated in Figure 2. This symmetry enables us to calculate the propagation direction of the reflected wave based on the angle of incidence. The propagation process of the incident P wave at the critical surface can be broadly divided into the following stages:

• The bursting of supercritical CO₂ causes the emission of a spherical P wave that reaches the critical surface.
  The P wave undergoes reflection, generating a corresponding surface wave.
• The reflected surface wave propagates to the surface and overlaps and superimposes with the original emitted wave, resulting in surface vibrations.

The reflection angles during the transmission and reflection of stress waves all obey Snell’s theorem, i.e.,:

$$\frac{\sin \beta }{\sin \alpha }=\frac{C_s}{C_p}=\sqrt{\frac{1-2v}{2\left(1-v\right)}}$$

(1)

where C_p and C_s are the P-wave and S-wave velocity, and ν is Poisson’s ratio of intact rock.

2.2. Ground Vibration Velocity Based on the Time-Domain Recursive Method

Aiming to describe the vibration effect caused by a plane P wave on the free surface, the surface wave superposition expression is based on the time-domain recursive method, which states that the particle vibration velocity is equal to the superposition of the body wave reaching this point, and to the second surface wave generated by the velocity difference of all adjacent particles in front of the point, as follows:

$$\left[\begin{array}{c}{V}_n\left(t\right)\\ V_{\tau }\left(t\right)\end{array}\right]=\left[\begin{array}{c}{V}_{nb}\left(t\right)\\ V_{\tau b}\left(t\right)\end{array}\right]+\underset{\Delta l\to 0}{\lim }{{\displaystyle \sum _{l_i=0}^L\left[\begin{array}{c}{V}_{ns}\left(t\right)\\ V_{\tau s}\left(t\right)\end{array}\right]}}_{l_i}$$

(2)

where V_n(t) and V_τ(t) are the tangential and normal particle vibration velocities at the measuring point. V_nb(t) and V_τb(t) are the body wave of the tangential and normal velocities at the measuring point. V_ns(t) and V_τs(t) are the tangential and normal particle vibration velocities of the surface wave, respectively, at the measuring point, which are caused by the difference in vibration velocities between adjacent particles.

According to the energy conservation of wavefront, the velocity of the body wave at the measuring point is affected by the geometric diffusion of wavefront material damping, and the calculation is as follows:

$$\left[\begin{array}{c}{V}_{nb}\left(t\right)\\ V_{\tau b}\left(t\right)\end{array}\right]=\left(C-B^{-1}AD\right)\cdot v_{Ip}\left(t\right)$$

(3)

$$v_{Ip}\left(t\right)={\left[\frac{r_1^2-r_0^2}{r_2^2-r_0^2}{\left(\frac{r_2}{r_1}\right)}^{1.5}\right]}^{0.5}\mathrm{e}\mathrm{x}\mathrm{p}\left(-\eta \cdot \Delta r\right)\cdot v_{Ip}^0\left(t-\frac{\Delta r}{c_p}\right)$$

(4)

The second surface wave velocity caused by the vibration velocity difference between two adjacent particles will attenuate with the distance to the measuring point, as shown in the calculation formula:

$${\left[\begin{array}{c}{V}_{ns}\left(t\right)\\ V_{\tau s}\left(t\right)\end{array}\right]}_{l_i}=\left[\begin{array}{c}{v}_{ns}\left(x_i,l_i/c_R\right)\cdot {\left(\frac{x_i}{L}\right)}^{0.5}\cdot \exp \left(-\eta l_i\right)\\ v_{\tau s}\left(x_i,l_i/c_R\right)\cdot {\left(\frac{x_i}{L}\right)}^{0.5}\cdot \exp \left(-\eta l_i\right)\end{array}\right]$$

(5)

$$\left[\begin{array}{c}{v}_{ns}\left(x_i,t\right)\\ v_{\tau s}\left(x_i,t\right)\end{array}\right]={\left(C-B^{-1}AD\right)}_{\left(i+1\right)}\cdot v_{Ip\left(i+1\right)}\left(t\right)-{\left(C-B^{-1}AD\right)}_{\left(i\right)}\cdot v_{Ip\left(i\right)}\left(t\right)$$

(6)

$$A=\left[\begin{array}{l}{Z}_p\cos 2\beta \\ Z_p\sin 2\beta \tan \beta /\tan \alpha \end{array}\right]$$

(7)

$$B=\left[\begin{array}{ll}{Z}_p\cos 2\beta & -Z_s\sin 2\beta \\ -Z_p\sin 2\beta \tan \beta /\tan \alpha & -Z_s\cos 2\beta \end{array}\right]$$

(8)

$$C=\left[\begin{array}{l}\cos \alpha \\ \sin \alpha \end{array}\right]$$

(9)

$$D=\left[\begin{array}{ll}\cos \alpha & \sin \beta \\ \sin \alpha & \cos \beta \end{array}\right]$$

(10)

where α and β represent the incident angle and reflection angle obtained using the Snell theorem, respectively. v_Ip(t) is the particle vibration velocity when the incident body wave arrives at the measuring point, and the velocity is mainly affected by the geometric diffusion of the wavefront and material damping.

3. Numerical Simulation Validation

3.1. Establishment of Numerical Model of Supercritical CO₂ Cracking

To further verify the reliability of the theoretical model, LS-DYNA is used to establish the simulation model of phase change cracking, assuming that the rock mass is a complete elastic body. The specific parameters are as follows in

Table 1. Rock material parameters.

Density (kg/m−3)	Elastic modulus (K) (GPa)	Shear modulus (G) (GPa)	P wave speed (cp) (m/s)	S wave speed (cs) (m/s)	Surface wave velocity (cr) (m/s)
2650	90	38.8	5828	3826	3445

:

The supercritical CO₂ cracking tubes are arranged in the manner shown in Figure 3. The buried depth of the fracturing tubes is 3 m, the resistance line from the free surface is 2 m, the model size is 4 × 6 × 0.01 m, and the fracturing hole diameter is 0.1 m. The non-reflective boundary is used to avoid the reflection of the stress wave on the calculated boundary, except on the free surface. We can calculate the ground vibration caused by CO₂ fracturing by applying a fracturing pressure curve to the fracturing hole. Measuring points are set up at the horizontal distance of 0.5, 1.0, 1.52, 2.0, and 2.5 m to monitor the horizontal and vertical closing velocities at each point.

In the modeling process, SOLID164 solid cells are used to simulate the rock mass; each cell has eight nodes, and the vibration effect of the rock mass can be fully displayed through the node displacement. After the trial analysis, the mesh size of the rock body near the rock burst hole is divided into 20 mm, the mesh size of the rock body far from the burst hole is divided into 100 mm, and the finite element model is shown in Figure 4. At the same time, in order to improve the calculation efficiency, we must fully consider the symmetry of the calculation model, and only let the 1/4 model be involved in the calculation. The calculation model is shown in Figure 5.

3.2. Analysis of Numerical Model Calculation Results

The numerical simulation results are shown in Figure 6. From the calculation, it can be seen that the wavefront excited by fracturing at the initial stage of fracturing diverges outward in the form of a concentric body wave at the center of the fracturing hole, and the vibration velocity of particles decays with the increase in transmission distance. When the wavefront is transmitted to the free surface, the vibration of the ground surface is caused by transmission and reflection.

A monitoring point was set at (0.2,0.2,0) coordinates, i.e., 0.282 m from the center of the fracturing hole, to monitor the vibration velocity of the surface plasma point after supercritical CO₂ blasting, and the time curve of the vibration velocity of the plasma point was obtained, as shown in Figure 7. As can be seen from the figure, the vibration velocity of the mass appeared two peak points, respectively, in 1.0 ms and 1.7 ms, the maximum vibration velocity of 3.56 cm/s. The second peak point is due to the vibration of the mass caused by the reflection of the stress wave at the critical surface and then returned to the measurement point.

3.3. Comparison of Numerical Simulation and Theoretical Calculation

The vibration velocity of each point on the free surface is calculated using the time-domain recursive method programmed with Matlab, and using numerical simulation. Additionally, the normal (Y) and tangential (X) particle vibration velocities at each measuring point are shown in Figure 8.

According to the calculation results, the difference table between numerical simulation and theoretical calculation results can be obtained as follows in

Table 2. Difference between numerical simulation and theoretical calculation results.

Blastpoint Distance (m)	x Direction	y Direction	Combined Speed
0.5	15.81%	14.51%	14.62%
1.0	11.16%	11.91%	11.70%
1.5	6.04%	7.30%	6.67%
2.0	5.36%	8.51%	6.67%
2.5	−4.52%	1.94%	−2.93%

:

It can be seen from Figure 8 that the horizontal and tangential velocity time history curves of each measuring point obtained using theoretical calculation and numerical simulation coincide, and the peak velocity difference between them is shown in Table 2; the maximum difference is 15.81%. It can be seen that the calculation formula based on the time-domain recursive method, used to achieve the superposition of the time domain through Matlab, can better calculate the surface vibration effect caused by supercritical CO₂ phase change cracking.

The variation rule of the peak velocity of each measured point with fracturing distance obtained via the theoretical model and numerical simulation is shown in Figure 9. From the figure, it can be seen that the peak velocity obtained using theoretical calculation and numerical simulation is consistent; the normal velocity decreases linearly and rapidly with distance, the tangential velocity increases exponentially with distance, and the combined velocity decreases linearly.

4. Comparison of Vibration Effects between Supercritical CO₂ Cracking and Dynamite Blasting

4.1. Equivalent Substitution of Supercritical CO₂ and Explosives

To compare the ground vibration velocity between phase change cracking and dynamite blasting, the energy released by supercritical CO₂ phase change cracking was calculated using the calculation model of the blasting energy of a compressed gas and water vapor container [23], and then converted into an equivalent explosive. The energy released by supercritical CO₂ phase change cracking can be calculated using the following formula:

$$E_g=\frac{PV}{k-1}[1-{(\frac{0.1013}{p})}^{\frac{k-1}{k}}]\times {10}^3$$

(11)

where E_g is the energy of high-pressure gas blasting, kJ. P is the absolute pressure of gas in the liquid storage tank, MPa. V is the volume of the liquid storage tank; k is the adiabatic index of gas, and the adiabatic index of CO₂ is 1.295.

The equivalent explosive TNT with the energy released by phase change cracking is calculated using the following formula:

$$W_{TNT}=\frac{E_g}{Q_{TNT}}$$

(12)

where W_TNT is the quality of chemical explosives. Q_TNT is the value of the energy of 1 kg TNT, which can be determined according to the selected explosive performance parameters.

Typical RDX-based explosives are selected as comparative chemical explosives. For chemical explosives, see the following

Table 3. Explosive parameters.

Ρ density (g/cm3)	D Burst speed (km/s)	A (pa)	B (pa)	R1	R2
1.65	8.19	937.29	16.279	5.2	1.0
ω	E0 (J/mm3)	γC	PC (Gpa)	QC (kJ/g)
0.43	8.5	3.0	27.67	5.152

[24]. The value of Q_TNT is 5151 kJ.

According to the results of the literature [25], the volume size of a type-85 fracturing tube is 1597.97 cm³, the liquid charge size is 1210 g, and the fracturing peak pressure is 242.09 mpa. The pressure–time curve of the type-85 fracturing tube in the process of fracturing the rock is shown in Figure 10, and it can be seen from the figure that the fracturing pressure is triangularly distributed, the peak fracturing pressure is about 242.09 MPa, and the pressure duration is 0.02 s.

Using the node merging method, the parameters of a 0.25 kg typical RDX-based explosive and the drilling parameters of the same rock are calculated using LS-DYNA, and then the pressure–time history curve of the wall detonation of RDX-based explosives can be obtained, as shown in Figure 11. As shown in the figure, the explosion pressure is also triangular, but the peak pressure is about 812.5 MPa, and the fracturing pressure duration is about 200 us. It can be seen that the peak pressure of explosive blasting is 3.36 times that of phase change cracking, but its action time is only 1/100 of that of phase change cracking.

4.2. Comparison of Vibration Velocities Caused by Supercritical CO₂ and Explosives

Normal (Y) and tangential (X) particle vibration velocities are obtained by calculating the vibration effect of the free surface under the action of equivalent explosive blasting with the same calculation method. The normal and tangential peak particle vibration velocities of the free surface are shown in Figure 12. The evolutionary trend of the normal and tangential peak particle vibration velocities is consistent with supercritical CO₂ phase change.

However, the ground vibration caused by supercritical CO₂ cracking is quite different from that caused by dynamite blasting after equivalent substitution. The following Figure 13 shows the attenuation curve of particles with increasing distance caused by supercritical CO₂ cracking and dynamite blasting.

The difference in particle vibration velocities between the two cracking modes is shown in

Table 4. Comparison of particle vibration velocity between two cracking modes.

Distance	Explosives Exploded (m/s)	Supercritical CO2 Cracking (m/s)	Ratio
0.5	2	0.0257	77.82
1.0	1.94	0.0249	77.91
1.5	1.85	0.0238	77.73
2.0	1.73	0.0224	77.23
2.5	1.6	0.0211	75.83
5	1.16	0.0155	74.84
10	0.75	0.0010	75.38
20	0.4	0.0054	74.63
30	0.25	0.0033	75.53
40	0.16	0.0022	73.73
50	0.11	0.0015	74.32

. It can be seen from the calculation results that the vibration velocities of supercritical CO₂ phase change cracking are only 1/74~1/78 of those of explosive blasting. Compared with the particle vibration caused by the explosive explosion, supercritical CO₂ cracking has an obvious vibration reduction effect when the energy of the two cracking modes is equal.

5. Conclusions

To investigate the effect of surface vibration induced by supercritical CO₂ phase change fracturing, this paper describes the process of supercritical CO₂ blasting stress waves on rock mass according to the stress wave theory and time-domain recurrence method, and establishes a model for calculating the surface vibration velocity after supercritical CO₂ blasting, considering time variation; the LS-DYNA dynamic finite element program is used to establish a three-dimensional numerical model to verify the reliability of the established theoretical calculation model. The reliability of the established theoretical model was verified. Based on the blasting energy calculation models of compressed gas and water vapor containers, the ground vibration velocities and decay laws of supercritical CO₂ phase change fracturing and equivalent explosive blasting were compared and analyzed.

(1)

The surface vibration induced by supercritical CO₂ phase change blasting is divided into three stages: blasting-induced column surface waves, surface waves generated by trans reflection, and superposition of reflected waves and original waves, which in turn cause surface vibration.

(2)

The vibration velocity calculation model based on the time-domain recurrence method can better evaluate the surface vibration effect caused by supercritical CO₂ phase change fracturing, and the maximum error between the theoretical model and the numerical simulation is 15.81%, demonstrating a good overall fit.

(3)

The surface vibration velocity under the effect of supercritical CO₂ phase fracturing decays exponentially with the increase in fracture distance, but the maximum vibration velocity is only 3.56 cm/s.

(4)

The peak pressure of supercritical CO₂ phase change is 1/3.36 times that of explosive blasting, the action time is 100 times that of explosive blasting, and the peak vibration velocity is only 1/74~1/78 of that of equivalent explosive blasting, which has a good vibration damping effect.

This paper focuses on supercritical CO₂ fracturing technology and its ground vibration rate after blasting. Compared with dynamite blasting, supercritical CO₂ fracturing technology is cheaper, and the blasting energy is more controllable; it also avoids the disadvantages of the loud noise and significant pollution produced by dynamite blasting. However, there is the phenomenon of “flying pipe” involved in its actual use, and there is a certain safety hazard, which is also a problem that needs to be solved in the future to improve the device. In addition, the vibration velocity characterization method for the use of this blasting technology on the ground surface has not been thoroughly studied, and most of monitoring is field-based; thus, at present, there is no perfect theoretical support for the process of its use. In this paper, a combination of theoretical methods and numerical simulations is used to study the vibration velocity of the ground surface during supercritical CO₂ phase change blasting, and certain results have been obtained.