Introduction

With the rapid development of engineering infrastructure construction at home and abroad, engineering blasting technology is used in more and more engineering fields, such as highway tunnel excavation, mine excavation, nuclear island bedrock excavation, slope excavation, and building demolition1,2,3. Explosive blasting facilitates engineering construction and may generate blasting hazards such as blasting flying stones, building cracking, slope landslides, and personnel injuries4,5,6. Aiming at blasting engineering, a large number of scholars have done a lot of theoretical research and field blasting experiments. Extensive research has been carried out on rock blasting damage models, monitoring means of blasting vibration effect, propagation characteristics of blasting vibration waves in various terrains, and analysis of blasting vibration monitoring signals characteristics7,8,9,10. The formulas for rock damage range characteristics and blasting vibration velocity have been established, and the safety standards for rock damage and buildings (structures) caused by blasting vibration have been formulated, achieving extensive and profound theoretical basis and engineering practice results11]– [12.

Although explosive blasting is commonly utilized, it brings environmental concerns such as flying rocks, vibration, dust, and noise13]– [14. In contrast, CO2 fracturing technology presents a novel solution that generates lower vibration intensity than traditional blasting methods. CO2 fracturing technology offers substantial environmental advantages15. Compared to the high-risk nature of explosive blasting, the liquid carbon dioxide phase transition fracturing technology (see Fig. 1) is a new technology that can be used as a supplementary method for on-site blasting that is not suitable for explosive blasting. The liquid carbon dioxide phase transition fracturing technology is an advanced, safe, and effective rock-breaking technique16,17,18.

Fig. 1
figure 1

The schematic and physical diagram of CO2 blasting device.

When liquid carbon dioxide phase transition fracturing blasting explodes in rock and soil media, a certain amount of energy is generated. Some of these energies are used to break the rock and soil media, while others will form blasting vibration waves19,20. A blasting vibration wave is a random wave; due to the change of charge mode, blasting medium, and terrain, the blasting vibration effect shows different propagation laws. It is of certain value to take corresponding research measures to predict and control the blasting vibration effect and formulate reasonable safety control standards. To achieve accurate prediction and control of the blasting vibration effect, it is necessary first to identify the propagation laws of blasting vibration waves in various rock and soil media, terrain conditions, joint characteristics, and charge structures. Many factors affect the blasting vibration effect, among which terrain factors have a great impact21,22. When the blasting vibration wave propagates over undulating terrain, the law of blasting vibration wave will be different from that of a flat field; it can be seen that the terrain has an essential influence on the propagation of the blasting vibration wave. At present, liquid carbon dioxide phase transition fracturing blasting differs from explosive blasting; in the application of liquid carbon dioxide phase transition fracturing blasting in rock excavation, there is still a lack of research on the vibration effect of terrain on liquid carbon dioxide phase transition fracturing blasting. It is significant to study the propagation law of vibration waves of liquid carbon dioxide phase transition fracturing blasting under different terrain conditions.

In terms of the influence of terrain on blasting vibration waves for deep hole blasting in open-pit mines, to study the laws of “elevation effect” and “whipping effect”, through on-site blasting experiments, Chen et al.23 found that the blasting vibration velocity generally shows a decreasing trend, with a high altitude amplification effect and obvious directionality; At the same elevation, the phenomenon of “whipping effect” is evident. Bai et al.24 conducted bench blasting tests under different charge parameters; by analyzing each measuring point’s vibration velocity amplitude and waveform characteristics, the existence law and triggering reasons of elevation amplification effect during the propagation of blasting seismic waves were studied. The intensity of the elevation amplification effect is not positively correlated with the charge amount, and the generation of the elevation amplification effect is related to the bulging movement of the rock and soil surface in the blasting center area. For the high rock slope blasting, Fu et al.25 found a significant vibration amplification effect at the outer edge and middle point of the horseway based on finite element analysis software LS-DYNA. Considering excavation damage, the amplification degree in the near area is significantly enhanced, and the amplification effect of blasting vibration is more evident at the slope steps with interlayers. Regarding the stepped terrain, Wu et al.26 analyzed the elevation amplification effect of slope blasting vibration velocity based on numerical simulation and vibration monitoring experiments. The amplification effect of peak particle velocity only occurs at a specific distance and height from the explosion source; the magnification does not increase monotonically with the height of the steps. Yu et al.27 analyzed the propagation law of blasting vibration along the slope of deep rock foundation pits and studied the elevation amplification effect at different positions. Within a certain range of height, the peak particle velocity and main frequency of measuring point in the edge slope of the foundation pit gradually increased with the increase of height, then decreased beyond a certain height. The accuracy of the vibration velocity prediction formula, considering the elevation effect, is better than that of the traditional Sadovsky empirical formula. When there are significant terrain undulations, the terrain will significantly affect the propagation law of blasting vibration waves. During the blasting construction process, it is inevitable to cause blasting vibration and adversely affect the safety of surrounding structures. Due to the complex interaction between explosive seismic waves and concave terrain, artificial or natural concave terrain can amplify or reduce the blasting vibration. It is necessary to analyze complex interactions and evaluate the effect of amplifying or reducing blasting vibration in concave terrain. Li et al.28 conducted on-site experiments and obtained the velocity of ground particles; it was found that the blasting vibration effect before the concave terrain is amplified, while the blasting vibration effect at the bottom and behind the concave terrain is reduced. The boundary integral equation method was used to obtain the theoretical formula for the motion of the entire field. The effectiveness of this theoretical method was verified by comparing it with the results of finite element calculations and field experiments. The complex interactions between those seismic waves and a rectangular concave terrain were described in detail. Tang et al.29 used the dimensional analysis method to analyze the physical parameters related to blasting vibration and obtained a blasting vibration formula to consider the elevation amplification effect. The blasting vibration formula corresponding to the dimensional analysis method can correctly reflect the amplification effect of positive elevation difference and is highly accurate; the error in predicting blasting vibration velocity is relatively small, only 10–14%.

By sorting out relevant research and summarizing the theoretical analysis, indoor and field experiments, and engineering application, there are extensive research on liquid carbon dioxide phase transition fracturing blasting and the influence of terrain on blasting vibration waves30,31,32,33. The research on liquid carbon dioxide phase transition fracturing blasting mainly focuses on the fracturing law of rock masses, the performance of fracturing pipes, rock breaking during subway tunnel excavation, and coal seam permeability fracturing. The influence of terrain on blasting vibration waves mainly targets the vibration waves caused by explosive blasting and earthquakes34,35. Liquid carbon dioxide phase transition fracturing blasting is different from explosive blasting. In the application process of rock mass excavation caused by liquid carbon dioxide phase transition fracturing blasting, there is currently a lack of systematic research on the influence of terrain on liquid carbon dioxide phase transition fracturing blasting waves. In view of the shortcomings of the above research, this paper conducts a study on the propagation law of liquid carbon dioxide phase transition fracturing blasting vibration waves under different terrain conditions, the influence laws of flat terrain, concave terrain, convex terrain, and slope terrain on the peak particle velocity of liquid carbon dioxide phase transition fracturing blasting are obtained, which is significant for controlling the safety of rock masses and structures.

Descriptions for field tests and PPV measurements

The nuclear power plant site where the CO2 blasting test is located in Cangnan County, Zhejiang Province, China, and the plant site is mainly hilly terrain. During the construction of the site area, mountain excavation, and site leveling above the ground design elevation shall be carried out first. Then, the negative excavation blasting of the foundation pit of the main structure shall be carried out, mainly in the form of layered steps. The test site for this study is the site leveling layer or negative excavation layer near the design elevation. The physical and mechanical parameters of the site rock mass are shown in Table 1. The rock is mainly moderately weathered and slightly weathered tuff, with a density of 2640 kg/m, an elastic modulus of 73.2GPA, a shear modulus of 29.2GPA, a Poisson’s ratio of 0.23, a uniaxial compressive strength of 107.5 MPa, a tensile strength of 5.5 MPa, a cohesion of 7.8 MPa, and an internal friction angle of 56.3°.

Table 1 The physical and mechanical parameters of the site rock mass.

In order to study the blasting vibration characteristics in the process of rock breaking by liquid carbon dioxide phase change, the field monitoring of blasting vibration signal was conducted. Referring to the vibration monitoring scheme of the traditional explosive blasting, a vibration measurement scheme with multiple monitoring points arranged on the vertical line of CO2 blasting holes was proposed. The monitoring points of blasting vibration on site were placed on the steps. The CO2 blasting vibration monitoring point layout diagram and test site are shown in Fig. 2. Three vector sensors that can monitor the horizontal radial, horizontal tangential, and vertical directions are arranged at each measuring point to monitor the vibration characteristics of each measuring point during the bench excavation of CO2 blasting.

The carbon dioxide fracturing tube of 100# type is used in the test, and its parameters are shown in Table 2. The parameters of hole network arranged on the site are shown in Table 3. Five blast holes are detonated in a single row at a time, and the borehole diameter and depth are 100 mm and 3.0 m, respectively. The minimum resistance line is 2 m, and the hole spacing is 3.0 m.

Fig. 2
figure 2

The CO2 blasting vibration monitoring point layout diagram and test site.

Table 2 The parameters for of 100# type carbon dioxide fracturing tube.
Table 3 The parameters of hole network arranged on the site.

Descriptions for test scheme considering terrain influence

The main research objective is to analyze the influence of terrain on the vibration wave propagation law of liquid carbon dioxide fracturing blasting. Based on the use of liquid carbon dioxide phase transition fracturing technology in engineering blasting, the numerical simulation of terrain propagation law on vibration wave of liquid carbon dioxide phase change fracturing blasting was carried out.

The numerical model of flat, concave, convex, and slope terrain was established in the UDEC program. The typical liquid carbon dioxide fracturing blasting vibration stress wave in the form of a Gaussian function stress wave (Fig. 3) was used as the input load41. The peak particle velocity of liquid carbon dioxide fracturing blasting vibration wave under various terrains was analyzed, and the influence of concave terrain width, concave terrain depth, blast-center distance, convex terrain height, slope terrain horizontal distance, and slope terrain height on the propagation law of peak particle velocity was obtained.

Applicability of UDEC software to blasting analysis

UDEC (Universal Distinction Element Code) is a two-dimensional numerical analysis program based on the theory of the discrete element method, which is specifically designed for handling discontinuous media36,37. The effectiveness of analysis results have been recognized by a large number of engineering workers. UDEC not only effectively solves practical engineering problems, but also has excellent accuracy and stability, making it one of the important tools in the field of geotechnical engineering.

UDEC can simulate the changes in multiple parameters such as displacement, stress, strain, and temperature of continuous or discontinuous media under static and dynamic forces. Discontinuous media can be regarded as discrete blocks, each set as rigid or deformable. To simulate the physical and mechanical properties of various objects, UDEC software has built-in multiple constitutive models, including the isotropic elastic model, transverse isotropic elastic model, D-P plastic model, M-C plastic model, and double yield model.

Numerous domestic and foreign researchers have extensively studied the propagation laws of vibration waves caused by explosive blasting in different terrains and joint distributions using UDEC38,39,40. By comparing with the results of on-site explosive blasting tests, it has been found that UDEC can effectively simulate the propagation of blasting vibration waves with excellent accuracy.

Determination on load of liquid carbon dioxide phase transition

In the process of carbon dioxide phase transition, supercritical carbon dioxide gas forms a high-pressure impact load around the drill hole. Under the shock wave and high-pressure gas, the rock mass causing the drill hole wall is subjected to both radial compressive stress and tangential tensile stress.

When the radial compressive stress of rock mass is much larger than its dynamic compressive strength, shock wave and high-pressure gas will cause a crushing ring phenomenon in the rock mass around the drill hole.

Zeng41 and Wang42 established a stress wave propagation model and analyzed the generation and propagation of stress waves in rock caused by the expansion of liquid carbon dioxide phase transition. When the liquid carbon dioxide high-pressure gas collides with the hole wall, the rock stress wave will be aroused in the borehole wall, and the following formula can accurately calculate the corresponding initial pressure P value.

$$P=n{P_1}$$
(1)
$${P_1}={P_m}{\left( {\frac{{{{\text{d}}_g}}}{{{d_0}}}} \right)^6}$$
(2)

where P is the initial pressure of rock stress wave caused by the impact of liquid carbon dioxide high-pressure gas on the hole wall, MPa; n is the pressure amplification factor after the impact, and the value of n ranges from 0 to 20. P1 is the CO2 high-pressure gas pressure before the collision hole wall, MPa; Pm is the pressure threshold under the limit failure state set by the energy release plates of the fracturing tube, MPa; dg and d0 are the outer diameter of fracturing tube and the diameter of fracturing drilling hole, respectively, mm.

As the shock wave formed by the fracturing of liquid carbon dioxide propagates towards the surrounding rock mass of the drill hole, the energy will attenuate. For the rock mass at a distance of r from the center of drill hole, the radial compressive stress can be expressed as:

$${\sigma _r}=P{\left( {\frac{{{{\text{r}}_{\text{b}}}}}{r}} \right)^\alpha }=P{\left( {\bar {r}} \right)^{ - \alpha }}$$
(3)

where \({\sigma _r}\) indicates radial compressive stress, MPa; \({r_b}\) represents the radius of borehole, m; r represents the distance from the center of the borehole, m; \(\bar {r}\) represents the specific distance: \(\bar {r}={r \mathord{\left/ {\vphantom {r {{r_b}}}} \right. \kern-0pt} {{r_b}}}\), \(\alpha\) represents the stress wave attenuation index, where the positive sign corresponds to the crushing zone and the negative sign corresponds to the crack zone; \(\mu {}_{{\text{d}}}\) represents the dynamic Poisson ratio of the rock, ranging from 0.25 to 0.40.

For the liquid carbon dioxide phase transition fracturing blasting load, the same as explosive blasting, the impact load pressure of blasting hole wall, the edge of the blasting hole crushing zone, and the edge of the near zone at a certain distance from the blasting hole can be calculated according to the Eqs. (1 ~ 3). The impact load pressure is applied at the corresponding position.

According to Zeng43, the carbon dioxide phase transition fracturing tube’s peak pressure range was between 169.1 and 260.2 MPa, with a positive pressure duration of 21.1ms to 74.8ms. A Gaussian curve can be used as the pressure curve for liquid carbon dioxide phase transition fracturing blasting. Taking the 100# model fracturing tube shown in Table 4 as an example, the peak pressure is taken as 225 MPa, the total loading time is taken as 25ms, the time required for zero pressure to reach the peak pressure is taken as 12.5ms, the drilling hole diameter is 100 mm, and the outer diameter of the fracturing tube is 95 mm. Using Eqs. (1 ~ 3), the peak pressure at a distance of r = 0.50 m is 23.5 MPa, the pressure curve of carbon dioxide fracturing load at r = 0.50 m is shown in Fig. 3.

Table 4 Load pressure parameters of liquid carbon dioxide phase transition fracturing blasting.
Fig. 3
figure 3

Pressure curve of carbon dioxide fracturing load at r = 0.50 m.

Boundary conditions and damping setting of numerical model

The excavation site of CO2 blasting during the site leveling in the Cangnan nuclear power plant was selected as the engineering background. In order to simulate the propagation of stress wave in semi-infinite space, a numerical calculation model with a certain size truncated boundary is established. For dynamic calculation, when stress waves propagate to the truncated boundary, the stress waves will be reflected into the interior of the mode due to the existence of the calculation model boundary. In the real case, the stress wave will propagate outwards indefinitely and will not reflect into the model. Therefore, due to the setting of the truncation boundary, the calculated result will deviate from the actual value.

In order to prevent the reflection of stress waves on the boundary, UDEC usually sets a viscous boundary at the boundary of the numerical model, which is equivalent to setting a normal and tangential damper at the boundary of the model to provide normal and tangential resistance. After the stress wave is transmitted to the truncated boundary, the energy is absorbed and will not be reflected to the model at the truncated boundary surface44,45. Considering symmetry, an axisymmetric boundary is set on the left boundary of the model. To prevent blasting vibration waves from returning to the model’s interior at the edge of the numerical model, the viscous boundary is set on the right and bottom of the numerical model.

The purpose of setting damping in the dynamic response analysis of rock masses is to simulate the real energy attenuation process in rocks under dynamic loads. In geotechnical media, the natural decay of this movement and energy has a hysteresis effect, which has a frequency correlation. In UDEC dynamic analysis, Rayleigh damping is set in the model, which can achieve a good simulation of real field vibration wave attenuation. To adopt Rayleigh damping, it is necessary to determine the critical damping ratio ξcr and the vibration center frequency f of the rock mass. This article takes the critical damping ratio ξcr = 0.05 and the vibration center frequency f = 100 Hz by repeatedly calculating and referring to the research results of the literature44,45.

Numerical modeling and scheme of CO 2 blasting in flat terrain

Based on the experimental site conditions, a numerical model was established using the discrete element software UDEC, with a height of 25 m and a length of 50 m. Horizontal constraints were applied to both the left and right sides of the model, and horizontal and vertical constraints were applied at the bottom of the model. In addition, axisymmetric boundaries are set on the left boundary of the model, and viscous boundaries are applied to the right and bottom of the model. The liquid carbon dioxide fracturing blasting hole is located at the left of the model, and the fracturing load is applied at the position where the real high-pressure carbon dioxide fracturing pipe sprays. The schematic diagram of the flat terrain calculation model is shown in Fig. 4. An ideal elastoplastic model for rock material types is adopted. In addition, the Mohr-Coulomb strength criterion is adopted as the yield criterion. The grid division diagram for the numerical calculation of flat terrain is shown in Fig. 5; the calculation model is divided into grids using a step-by-step refinement method.

Fig. 4
figure 4

The schematic diagram of the flat terrain calculation model of CO2 blasting.

Fig. 5
figure 5

The grid division diagram for numerical calculation of flat terrain in CO2 blasting.

The influencing factor is set as blast-center distance to study the effect of blast-center distance on the vibration signal of CO2 blasting under flat terrain. A parametric study was carried out by using a single-factor control variable method. For the vibration monitoring of liquid carbon dioxide cracking blasting, peak particle velocity monitoring points are set at 3 m, 5 m, 7 m, 10 m, 12.5 m, 15 m, 17.5 m, 20 m, 25 m, 30 m, and 45 m away from the blasting center, respectively, to obtain the horizontal radial peak particle velocity and vertical peak particle velocity. The changes of horizontal radial peak particle velocity and vertical peak particle velocity with different blast-center distances under flat terrain are obtained.

Numerical modeling and scheme of CO 2 blasting in concave terrain

In recent decades, numerous studies have examined the impact of different terrains, such as upward step and ditches, on the propagation of blast waves. These studies have utilized field investigations, theoretical analysis, and numerical simulations to understand how these terrains affect blast wave propagation. It has been observed that upward step, or elevated terrains, tend to amplify blast vibrations, while damping ditches serve as vibration absorbers.

When the blasting vibration wave encounters concave terrain during the propagation of rock media, it is equivalent to the waves entering from one medium to another (rock mass entering the air). As a result, the waves undergo reflection and refraction on the concave surface, causing energy loss of the waves. Meanwhile, diffraction occurs at the edge of the concave surface, and the waves continue to propagate in the form of volume waves and surface waves beneath the concave surface. Therefore, it can be considered that concave landforms play a shielding role, cutting off the propagation of some waves, and the blasting vibration waves undergo attenuation effects.

To study the influence of concave terrain, such as damping ditch, on the vibration signal of CO2 blasting, analyze whether the vibration wave of liquid carbon dioxide phase transition fracturing blasting after passing through concave terrain is consistent with that of a flat site. The concave terrain calculation model site was established, a measuring point was placed on the edge of the concave back explosion source, and the horizontal radial velocity and vertical velocity were monitored simultaneously. The schematic diagram of field vibration monitoring points for concave terrain is shown in Fig. 6, and the grid division diagram for numerical calculation of concave terrain is shown in Fig. 7. Under the concave terrain, the peak particle velocity at the different concave terrain width, concave terrain depth, and blast-center distance is recorded each time.

According to the research findings of Tang et al.21 and Xia et al.22, the attenuation coefficient η of peak particle velocity for concave terrain is defined as:

$$\eta ={\text{ }}\frac{{PP{V_{flat{\text{ }}terrain}}-PP{V_{concave{\text{ }}terrain}}}}{{PP{V_{flat{\text{ }}terrain}}}}$$
(4)

where attenuation coefficient η is an amplitude change factor of CO2 blast vibration that has been introduced to reflect the PPV variation caused by concave terrain, \(PP{V_{flat{\text{ }}terrain}}\) is the calculated peak particle velocity of flat terrain, which can be calculated by data regression analysis of Eq. (7) and Eq. (8), \(PP{V_{concave{\text{ }}terrain}}\) is the monitored peak particle velocity of measuring point on the edge of the concave back explosion source.

Fig. 6
figure 6

The schematic diagram of the concave terrain calculation model of CO2 blasting.

Fig. 7
figure 7

The grid division diagram for numerical calculation of concave terrain in CO2 blasting.

To study the influence of concave terrain, such as damping ditch, on the vibration signal of CO2 blasting, the influencing factors are set as concave terrain width, concave terrain depth, and blasting center distance. A parametric study was carried out using a factor control variable method. The width of concave terrain is taken as 3 m, 4 m, 5 m, and 6 m, respectively; the depth of concave terrain is taken as 3 m, 4 m, 5 m, and 6 m, respectively; the blast-center distances are taken as 10 m, 15 m, 20 m, and 30 m, respectively. The liquid carbon dioxide phase transition fracturing blasting test schemes under different concave terrains are shown in Table 5. The horizontal radial peak particle velocity and vertical peak particle velocity under different concave terrain can be monitored, and the variation law of the attenuation coefficient of horizontal radial peak particle velocity and vertical peak particle velocity with different concave terrain width, concave terrain depth and blasting center distance can be obtained.

Table 5 The CO2 blasting test schemes under different concave terrains.

Numerical modeling and scheme of CO2 blasting in convex terrain

In the process of nuclear power construction, the rock foundation of the main structure such as nuclear island, conventional island and combined pumping station should be excavated by CO2 blasting. To speed up the progress of the project, the negative excavation of the engineering foundation is carried out in several places, and several convex terrains are formed in the foundation pit. In order to understand the influence of topography on the propagation characteristics of CO2 blasting vibration wave and obtain the attenuation law of CO2 blasting vibration wave, several numerical simulation studies were carried out using the convex terrain formed during negative excavation.

To study the influence of convex terrain on the vibration signal of CO2 blasting, analyze whether the vibration wave of liquid carbon dioxide phase transition fracturing blasting after passing through the convex terrain is consistent with that of a flat site, a convex terrain calculation model site is established, and two measuring points are arranged on the edge of the lower and upper corners of the convex platform, simultaneously monitor the horizontal radial velocity and vertical velocity. The layout diagram of field vibration monitoring points for convex terrain is shown in Fig. 8, and the grid division diagram for numerical calculation of convex terrain is shown in Fig. 9. Under convex terrain, the peak particle velocity of each measuring point at different convex terrain height and detonation distance is recorded each time.

Fig. 8
figure 8

The schematic diagram of the convex terrain calculation model of CO2 blasting.

Fig. 9
figure 9

The grid division diagram for numerical calculation of convex terrain in CO2 blasting.

According to the research findings of Tang et al.21 and Xia et al.22, the amplification factor β of peak particle velocity for convex terrain is defined as:

$$\beta ={\text{ }}\frac{{PP{V_{upper{\text{ }}corner}}}}{{PP{V_{lower{\text{ }}corner}}}}$$
(5)

where amplification factor β is an amplitude change factor of CO2 blast vibration that has been introduced to reflect the PPV variation caused by convex terrain, \(PP{V_{upper{\text{ }}corner}}\)and \(PP{V_{lower{\text{ }}corner}}\) is the monitored vibration velocity on the edge of the upper and the lower corner placed on the convex bench, respectively.

To study the influence of convex terrain on the vibration signal of fracturing blasting, the influencing factors are set as the convex terrain height and the blasting center distance. A parametric study was carried out by using a single-factor control variable method. The heights of convex terrain are taken as 1 m, 3 m, 6 m, 9 m, 12 m, 15 m, and 20 m, respectively. The blasting center distances are taken as 10 m, 12.5 m, 15 m, 17.5 m, 20 m, 25 m, and 30 m, respectively. The test scheme of liquid carbon dioxide phase transition fracturing blasting under different convex terrains is shown in Table 6. The horizontal radial peak particle velocity and vertical peak particle velocity under different convex terrain can be monitored, and the variation law of the amplifying coefficient of horizontal radial peak particle velocity and vertical peak particle velocity with different convex terrain height and blasting center distance can be obtained.

Table 6 The CO2 blasting test schemes under different convex terrains.

Numerical modeling and scheme of CO2 blasting in slope terrain

The schematic diagram of the slope terrain calculation model site is shown in Fig. 10. To study the influence of slope terrain on the vibration signal of fracturing blasting, the slope numerical model with a slope height of 75 m and different slope angles of 15°, 30°, 45°, and 60° was established. A numerical simulation of liquid carbon dioxide phase transition fracturing blasting vibration under different slope gradients was conducted to analyze the influence of different slopes and elevations on the peak particle velocity of CO2 blasting.

Fig. 10
figure 10

The schematic diagram of the slope terrain calculation model of CO2 blasting.

(1) The variation law of peak particle velocity on slope surface in the same slope.

In order to analyze the attenuation law of vibration velocity at different monitoring points on the slope surface in the same slope, it is necessary to arrange several measuring points on the outer edge of the slope surface to monitor the horizontal peak particle velocity and radial peak particle velocity. Using the discrete element method numerical simulation software UDEC, the slope model based on a slope project is taken as an example to establish a slope with a slope height of 75 m and a slope angle of 45°. The grid division diagram of slope topography numerical calculation is shown in Fig. 11.

Fig. 11
figure 11

The grid division diagram of slope terrain with a slope height of 75 m and a slope angle of 45°.

The fracturing hole is directly located at the foot of the slope, and the stress wave propagates in the slope body. A series of monitoring points were selected along the slope surface for statistical analysis to achieve the variation law of peak particle velocity along the slope surface and the elevation effect. The test scheme of the variation law of peak particle velocity on the same slope surface is shown in Table 7. In the test scheme, the horizontal distance and elevation difference between the numerical monitoring point and the blasting source are 1 m, 5 m, 7 m, 10 m, 15 m, 20 m, 25 m, 30 m, 35 m, and 40 m, respectively.

Table 7 The CO2 blasting test scheme of peak particle velocity on the same slope surface.

(2) The variation law of peak particle velocity inside the same slope body.

To analyze the variation law of peak particle velocity within the same slope body, multiple measuring points need to be arranged at the same horizontal distance and different elevations within the same slope body to monitor the peak particle velocity. As shown in Fig. 12, using the discrete element method numerical simulation software UDEC, a slope model was established with a height of 75 m and a slope angle of 30 °with the background of a slope project. A series of monitoring points are selected at the same level and different elevations to analyze the elevation effect of blasting vibration velocity inside the same slope.

Fig. 12
figure 12

The grid division diagram of slope terrain with a slope height of 75 m and a slope angle of 30°.

At horizontal distances of 10 m, 20 m, 30 m, and 40 m from the fracturing hole, four monitoring points were selected in sequence at the slope surface from top to bottom with the vertical distance of 0 m, 1 m, 2 m, and 3 m for monitoring the peak particle velocity. The test scheme of the variation law of peak particle velocity within the same slope body is shown in Table 8.

Table 8 The CO2 blasting test scheme of peak particle velocity within the same slope body.

(3) The variation law of peak particle velocity on slope surface in the different slope.

To analyze peak particle velocity on the slope surface in different slopes, it is necessary to arrange multiple measuring points that need to be placed on the slope surface at the same horizontal distance with different slopes to monitor the peak particle velocity. As shown in Fig. 13, Using the discrete element method numerical simulation software UDEC, multiple groups slope numerical models with a slope height of 75 m and slope angle of 15°, 30°, 45°, and 60° were established.

Fig. 13
figure 13

The grid division diagram of slope terrain with a slope height of 75 m and different slope angle.

The peak particle velocity monitoring point is set on the slope surface at the same horizontal distance from the fracturing hole. Table 9 shows the test scheme for the variation law of peak particle velocity on slope surfaces with different slopes.

Table 9 The CO2 blasting test scheme of peak particle velocity on slope surface in the different slope.

Results and discussions

Results and discussions on CO 2 blasting in flat terrain

The field test results of peak particle velocity under CO2 blasting are shown in Table 10. The numerical simulation results of peak particle velocity of liquid carbon dioxide phase transition fracturing blasting on flat terrain are shown in Table 11. The CO2 fracturing tests of rock mass at the excavation site of the foundation pit were carried out at the Cangnan nuclear power plant, and a vibration velocity analysis was conducted. By comparing the numerical simulation results in this paper with the field test results, the attenuation curves of peak particle velocity in horizontal radial and vertical directions of CO2 blasting can be obtained, as shown in Fig. 14.

Table 10 The field test results of peak particle velocity under CO2 blasting.
Table 11 The CO2 blasting numerical simulation results under the flat terrain.

As seen from Tables 10 and 11, and Fig. 14, the peak particle velocity in both horizontal radial and vertical directions exhibits a rapid decay pattern with increasing distance from the explosion center. For the blasting near field and far field from the fracturing blasting center, the numerical simulation results are in good agreement with the field test results in both horizontal radial and vertical directions. Through comparative analysis, the field test results are consistent with the numerical simulation results, which indicate the reliability of the research results.

Fig. 14
figure 14

The attenuation curves of peak particle velocity of CO2 blasting under the flat terrain.

In general, the relationship between the peak particle velocity PPV of rock mass and the explosive quantity Q and the detonation distance R is given by the Sadovsky formula46,47:

$$PPV={\text{ }}K{\left( {\frac{{\sqrt[3]{Q}}}{R}} \right)^\alpha }$$
(6)

where K is the coefficient related to the medium and the blasting conditions, α is attenuation coefficient related to the terrain and geological conditions between the blasting point and the monitoring point, Q is the maximum explosive weight per delay (kg).

The Sadovsky formula is in the form of a power function. To analyze the attenuation law of blasting vibration waves in flat terrain, the power function fitting was performed based on the attenuation curves of peak particle velocity of liquid carbon dioxide phase transition fracturing blasting in the horizontal radial and vertical directions. K, α in the horizontal radial and vertical directions were obtained by regression calculation from the peak particle velocity of each measuring point for field test results and numerical simulation results.

$$PP{V_X}={\text{ }}30.65{R^{-0.77}}\qquad{R^2}=0.92$$
(7)
$$PP{V_Z}={\text{ }}13.88{R^{-0.51}}\qquad{R^2}=0.97$$
(8)

In the Sadovski empirical formula, the coefficient K is directly proportional to PPV. α represents the attenuation rate of peak particle velocity, and the larger the value of α, the faster the attenuation rate of peak particle velocity with the increase of blast-center distance.

By comparing and analyzing the K and α values of the regression formulas for horizontal radial and vertical attenuation, the K and α values in the horizontal radial direction are 30.65 and 0.77, respectively, while those in the vertical direction are 13.88 and 0.51, respectively. Since the K in the horizontal radial direction is greater than that in the vertical direction, and the α in the horizontal radial direction is greater than that in the vertical direction, this indicates that although the peak particle velocity in the horizontal radial direction in the near field of explosion source is higher, as the blast-center distance increases, the peak particle velocity in the vertical direction will be greater than that in the horizontal radial direction. Therefore, a single directional component of blasting vibration velocity should not be used as a safety criterion for blasting vibration. When monitoring the vibration signal of liquid carbon dioxide phase transition fracturing blasting, the peak particle velocity in the horizontal and vertical directions should be monitored simultaneously each time, and the maximum value should be used as the safety criterion for blasting vibration.

The fitting correlation coefficients R2 of Eq. (7) and Eq. (8) are 0.92 and 0.97, respectively, indicating that the blasting vibration velocity regression analysis has an excellent correlation; a prediction formula can be established for predicting the blasting vibration velocity. In practical engineering applications, the safety factor should be appropriately increased to ensure absolute safety of the surrounding environment.

Results and discussions on CO 2 blasting in concave terrain

The numerical simulation results of peak particle velocity of liquid carbon dioxide phase transition fracturing blasting on different concave terrain are shown in Table 12. The attenuation coefficient h of the peak particle velocity in the horizontal radial and vertical directions defined by Eq. (4) is plotted as a function of concave terrain width b, concave terrain depth h, and blasting center distance R, as shown in Figs. 15, 16 and 17, respectively.

Table 12 The CO2 blasting numerical simulation results under the concave terrain.
Fig. 15
figure 15

The attenuation coefficient η of peak particle velocity varies with the concave terrain width.

Fig. 16
figure 16

The attenuation coefficient η of peak particle velocity varies with the concave terrain depth.

From Table 12; Fig. 15, it can be seen that for concave terrain with a depth of 4 m and a blast-center distance of 10 m, as the width of the concave terrain increases from 3 m to 6 m, the attenuation coefficient of the horizontal radial peak particle velocity increases from 37.6 to 49.6%; The attenuation coefficient of the vertical peak particle velocity increased from 9.6 to 22.1%; Increasing the width of concave terrain can significantly reduce the peak particle velocity in both horizontal radial and vertical directions. The same variation in concave terrain width has a greater impact on the attenuation coefficient of peak particle velocity in horizontal radial direction than vertical direction. It can be seen that the damping effect of concave groove increases with the increase of concave terrain width.

Fig. 17
figure 17

The attenuation coefficient η of peak particle velocity varies with the blasting center distance.

From Table 12; Fig. 16, it can be seen that for concave terrain with a width of 4 m and a blast-center distance of 10 m, as the concave terrain depth increases from 3 m to 6 m, the attenuation coefficient of the horizontal radial peak particle velocity increases from 27.6 to 69.6%; the attenuation coefficient of the vertical peak particle velocity increases from 7.3 to 26.1%; Increasing the depth of concave terrain can significantly reduce the peak particle velocity in both horizontal radial and vertical directions. Changing the depth of concave terrain has a greater impact on the attenuation coefficient of peak particle velocity in the horizontal radial direction than in the vertical direction. It can be seen that the damping effect of the concave groove increases with the increase of the concave terrain depth.

From Table 12; Fig. 17, it can be seen that for concave terrain with a width of 4 m and a depth of 4 m, as the blasting center distance increases from 10 m to 30 m, the attenuation coefficient of the horizontal radial peak particle velocity decreases from 41.1 to 8.9%, the attenuation coefficient of the vertical peak particle velocity decreased from 14.1 to 7.4%. The attenuation coefficients of horizontal radial and vertical peak particle velocity show a decreasing trend with the increase of blasting center distance. The effect of varying the same blasting center distance on the attenuation coefficient of horizontal radial peak particle velocity is greater than that on the attenuation coefficient of vertical peak particle velocity. It can be seen that the damping effect of concave terrain decreases with increasing blasting center distance, and the damping effect is more obvious in the zone near the explosion source.

The attenuation coefficients of horizontal radial and vertical peak particle velocity can be fitted with the following Eqs. (9 ~ 14) as a function of concave terrain width, concave terrain depth, and blasting center distance:

$$\eta_x= 0.040b+0.254\qquad R^{2}=0.99$$
(9)
$$\eta_z= 0.042b-0.027\qquad R^{2}=0.99$$
(10)
$$\eta_x= 0.139h-0.144\qquad R^{2}=0.99$$
(11)
$$\eta_z= 0.063h-0.113\qquad R^{2}=0.99$$
(12)
$$\eta_x= 9.91R^{-1.385}\qquad R^{2}=0.99$$
(13)
$$\eta_z= 0.483R^{-0.572}\qquad R^{2}=0.89$$
(14)

The fitting correlation coefficients R2 of Eqs. (9 ~ 14) are 0.99, 0.99, 0.99, 0.99, 0.99, and 0.89, respectively. This indicates that the peak particle velocity attenuation coefficient of the concave terrain used for regression analysis is correlated with the concave terrain width, depth, and blasting center distance. The prediction formulas can be established to predict the peak particle velocity attenuation coefficient of liquid carbon dioxide phase transition fracturing blasting.

Results and discussions on CO 2 blasting in convex terrain

Table 13 shows the numerical simulation results of the peak particle velocity of liquid carbon dioxide fracturing blasting on different convex terrains. The amplification coefficient β of peak particle velocity in the horizontal radial and vertical directions defined by Eq. (5) is plotted as a function of the convex terrain height h and the blasting center distance R, as shown in Figs. 18 and 19, respectively.

Table 13 The CO2 blasting numerical simulation results under the convex terrain.
Fig. 18
figure 18

The amplification coefficient β of peak particle velocity varies with the convex terrain height.

Fig. 19
figure 19

The amplification coefficient β of peak particle velocity varies with the blasting center distance.

From Table 13; Fig. 18, it can be seen that for a blast-center distance of 10 m, as the convex terrain height increases from 1 m to 6 m, the amplification coefficient of the peak particle velocity in the horizontal radial direction increases from 101 to 158%; Subsequently, as the convex terrain height increases from 6 m to 20 m, the amplification coefficient of the horizontal radial peak particle velocity gradually decreases from 158 to 120%. For the vertical direction, as the convex terrain height increases from 1 m to 12 m, the amplification coefficient of the peak particle velocity in the vertical direction increases from 101 to 209%; then, as the convex terrain height increases from 12 m to 20 m, the amplification coefficient of vertical peak particle velocity gradually decreased from 209 to 137%. It can be seen that the amplification coefficient of peak particle velocity in the horizontal radial and vertical directions does not monotonically increase with the increase of the elevation difference of the convex terrain. As the elevation difference of the convex terrain increases, the amplification coefficient of the peak particle velocity in the horizontal radial and vertical directions shows a characteristic of first increasing and then gradually decreasing. When the height range of convex terrain is 1 ~ 6 m, the amplification coefficient in the horizontal radial direction is slightly larger than that in the vertical direction; when the height range of convex terrain is 6 ~ 20 m, the vertical amplification coefficient is slightly larger than the horizontal radial amplification coefficient.

From Table 13; Fig. 19, it can be seen that for a convex terrain height of 9 m, as the blasting center distance increases from 10 m to 30 m, the amplification coefficient of peak particle velocity in the horizontal radial direction gradually decreases from 155 to 106%; the amplification coefficient of peak particle velocity in the vertical direction gradually decreases from 190 to 157%. It can be seen that as the blasting center distance increases, the amplification coefficient β of peak particle velocities in both horizontal and vertical directions exhibits a monotonically decreasing characteristic. When the distance range of the explosion center is 10 ~ 30 m, the amplification coefficient of peak particle velocity in the vertical direction is always greater than that in the horizontal direction.

The variation curves of peak particle velocity amplification coefficient in the horizontal radial and vertical direction with the convex terrain height and the blasting center distance can be fitted to the following Eqs. (15 ~ 18):

$$\beta_x=- 0.005h^2+0.109h+0.988\qquad R^{2}=0.88$$
(15)
$$\beta_z=-0.008h^{2}+0.179h+0.733\qquad R^{2}=0.72$$
(16)
$$\beta_x=-0.024R+1.793\qquad R^{2}=0.99$$
(17)
$$\beta_x=-0.024R+1.793\qquad R^{2}=0.96$$
(18)

The fitting correlation coefficients R2 of Eqs. (15 ~ 18) are 0.88, 0.72, 0.99, and 0.96, respectively, indicating that the amplification coefficient of peak particle velocity of convex terrain used for regression analysis correlates with the convex terrain height and the blasting center distance. The prediction formula can be used to predict the amplification coefficient of the peak particle velocity of liquid carbon dioxide phase transition fracturing blasting.

Results and discussions on CO 2 blasting in slope terrain

(1) The variation law of peak particle velocity on slope surface in the same slope.

The numerical results of the peak particle velocity variation law on the same slope surface are shown in Table 14. The variation curves of horizontal radial peak particle velocity and vertical peak particle velocity along with horizontal distance R and elevation difference H are drawn by using data analysis software origin, as shown in Fig. 20.

Table 14 The numerical results of peak particle velocity on the same slope surface.
Fig. 20
figure 20

The variation curve and fitting surface of PPV with horizontal distance and elevation difference.

Table 14; Fig. 20 show that for slope terrain, the horizontal radial peak particle velocity and vertical peak particle velocity generated by liquid carbon dioxide fracturing blasting decrease with the increase of horizontal distance and elevation difference. When the horizontal distance and elevation difference are both 1 m, the horizontal radial peak particle velocity and vertical peak particle velocity are 3.05 cm/s and 4.08 cm/s, respectively; When the horizontal distance and elevation difference are both 10 m, the horizontal radial peak particle velocity and vertical peak particle velocity are 1.84 cm/s and 1.46 cm/s, respectively, with a reduction of 39.7% and 63.2%, respectively; When the horizontal distance and elevation difference are both 30 m, the horizontal radial peak particle velocity and vertical peak particle velocity are 0.73 cm/s and 0.95 cm/s, respectively, with a reduction of 76.1% and 76.7%, respectively. It can be seen that with the increase of the horizontal distance and elevation difference, the horizontal radial and vertical peak particle velocity shows a decreasing trend and has a significant attenuation characteristic along the slope surface.

Referring to the Sadovsky empirical Eq. (6) for nonlinear fitting, the attenuation formulas for the horizontal radial peak particle velocity and vertical peak particle velocity with horizontal distance R and elevation difference H can be obtained as follows:

$$PP{V_x}={\text{ }}3.29 \cdot {R^{-0.057}} \cdot {H^{-0.285}}{R^2}=0.87$$
(19)
$$PP{V_z}={\text{ }}4.19 \cdot {R^{ - 0.215}} \cdot {H^{ - 0.215}}{R^2}=0.{\text{97}}$$
(20)

The fitting correlation coefficients R2 of Eqs. (19 ~ 20) are 0.87 and 0.97, respectively; it shows that the blasting vibration velocity of the monitoring data used for regression analysis has a good correlation and can intuitively reflect the influence of horizontal distance and height difference on the peak particle velocity. A prediction formula can be established to predict the peak particle velocity of blasting. In practical engineering applications, the safety factor should be appropriately increased to ensure the absolute safety of the surrounding environment.

(2) The variation law of peak particle velocity inside the same slope body.

Table 15 shows the numerical results of the peak particle velocity variation law inside the same slope body. Figure 21 shows the variation curve of the peak particle velocity elevation effect for interior points on the same slope under the liquid carbon dioxide phase transition fracturing blasting.

Table 15 The numerical results of peak particle velocity inside the same slope body.
Fig. 21
figure 21

The variation curve of the peak particle velocity elevation effect for interior points on the same slope.

According to Table 15; Fig. 21, for a slope terrain with a height of 75 m and a slope angle of 30 °; For the horizontal radial direction, at a horizontal distance of 10 m, as the elevation difference between the monitoring point and blasting source decreases from 5.77 m to 2.77 m, the distance between the monitoring point and the slope surface varies from 0 to 3 m, the peak particle velocity gradually decreases from 2.87 cm/s to 2.59 cm/s, showing the characteristic that the peak particle velocity increases with the increase of elevation, the peak particle velocity of monitoring point near the slope surface is larger. When the horizontal distance is 20 m, 30 m, and 40 m, as the elevation difference between monitoring point and blasting source decreases from 11.54 m to 8.54 m, from 17.32 m to 14.32 m, and from 23.09 m to 20.09 m, the distance between the monitoring point and the slope surface varies from 0 to 3 m, the peak particle velocity gradually increases from 1.72 cm/s to 1.92 cm/s, from 0.85 cm/s to 1.09 cm/s, and from 0.59 cm/s to 0.67 cm/s, showing that the peak particle velocity decreases with the increase of elevation, the peak particle velocity of monitoring point near the slope surface is smaller. For the vertical direction, at horizontal distances of 10 m, 20 m, 30 m, and 40 m, as the elevation difference between the monitoring point and blasting source decreases from 5.77 m to 2.77 m, from 11.54 m to 8.54 m, from 17.32 m to 14.32 m, and from 23.09 m to 20.09 m, the distance from the monitoring point to the slope surface varies from 0 to 3 m, the peak particle velocity gradually decreases from 1.88 cm/s to 1.71 cm/s, from 0.97 cm/s to 0.86 cm/s, from 0.86 cm/s to 0.77 cm/s, and from 0.83 cm/s to 0.68 cm/s, respectively; showing that the peak particle velocity increases with elevation and the peak particle velocity of monitoring point near the slope surface is larger.

It can be seen that at the same horizontal distance within the slope body, with the increase of elevation, although the distance between the monitoring points and the blasting center increases, the peak particle velocity shows inconsistent trends, indicating the influence of elevation on the peak particle velocity of slope blasting. For the horizontal radial direction, the peak particle velocity in the near field of CO2 blasting shows an amplification effect with the increase of elevation, and the peak particle velocity in the far field of CO2 blasting shows a decreasing effect with the increase of elevation. For vertical direction, the peak particle velocity in the near field and far field of CO2 blasting shows an amplification effect with the increase of elevation.

(3) The variation law of peak particle velocity on slope surface in the different slope.

Table 16 shows the numerical results of the peak particle velocity variation law on the slope surface at different slopes. The variation curve of the peak particle velocity on the slope surface at the same horizontal distance of different slopes is shown in Fig. 22.

Table 16 The numerical results of peak particle velocity on slope surface in the different slope.
Fig. 22
figure 22

The variation curve of the peak particle velocity on slope surface in the different slope.

According to Table 16; Fig. 22, For horizontal radial direction, for horizontal distances of 5 m, 10 m, 20 m, 30 m, and 40 m, as the elevation difference between the monitoring point and blasting source increases from 1.33 m to 8.66 m, from 2.67 m to 17.3 m, from 5.35 m to 34.6 m, from 8.03 m to 51.9 m, and from 10.7 m to 69.2 m, the slope gradient increases from 15° to 60°, the peak particle velocity of the slope surface shows a decreasing trend; the peak particle velocity decreased from 5.80 cm/s to 1.45 cm/s, from 4.12 cm/s to 1.13 cm/s, from 2.44 cm/s to 0.71 cm/s, from 1.93 cm/s to 0.43 cm/s, and from 1.71 cm/s to 0.24 cm/s, respectively. In the vertical direction, for horizontal distances of 5 m, 10 m, and 40 m, as the elevation difference between monitoring point and blasting source increases from 1.33 m to 8.66 m, from 2.67 m to 17.3 m, and from 10.7 m to 69.2 m, the slope gradient increases from 15° to 60°, the peak particle velocity of the slope surface shows a decreasing trend; the peak particle velocities decreased from 4.68 cm/s to 1.72 cm/s, from 1.93 cm/s to 1.13 cm/s, and from 0.88 cm/s to 0.44 cm/s, respectively; for horizontal distances of 20 m and 30 m, as the elevation difference between monitoring point and blasting source increases from 5.35 m to 34.6 m, and from 8.03 m to 51.9 m, the slope gradient increases from 15° to 60°, the peak particle velocity of the slope surface shows a trend of first decreasing, then increasing, and then decreasing again; When the slope gradient is 15°, 30°, 45°, and 60°, the peak particle velocities at a horizontal distance of 20 m are 1.05 cm/s, 0.97 cm/s, 1.12 cm/s, and 1.02 cm/s, respectively; the peak particle velocities at a horizontal distance of 30 m are 0.99 cm/s, 0.86 cm/s, 0.95 cm/s, and 0.79 cm/s, respectively.

Conclusions

In this study, in response to the influence of terrain on the propagation law of vibration waves in liquid carbon dioxide phase change fracturing blasting, the numerical models of flat terrain, concave terrain, convex terrain and slope terrain were established by using UDEC. The propagation law of vibration waves in liquid carbon dioxide phase change fracturing blasting in different terrains was obtained. The main conclusions are as follows:

  1. 1.

    For flat terrain, the peak particle velocity in both horizontal radial direction and vertical direction exhibit a rapid decay pattern with increasing distance from the explosion center. For the blasting near field and far field from the fracturing blasting center, the field test results are consistent with the numerical simulation results, which indicate the reliability of the research results.

  2. 2.

    For concave terrain, increasing the width and depth of the concave terrain can significantly reduce the horizontal radial and vertical peak particle velocity. As the blast-center distance increases, the attenuation coefficients of the horizontal radial and the vertical peak particle velocity both show a decreasing trend. For convex terrain, as the elevation difference in convex terrain increases, the amplification coefficient of peak particle velocity in horizontal radial and vertical directions shows a characteristic of first increasing and then gradually decreasing. As blast-center distance increases, the amplification coefficient of peak particle velocity in both horizontal and vertical directions exhibits a monotonically decreasing characteristic.

  3. 3.

    For slope terrain, in terms of the variation law of peak particle velocity on slope surface in the same slope; as the horizontal distance and elevation difference increase, the horizontal radial and vertical peak particle velocity significantly decreases along the slope surface. The attenuation formulas for the horizontal radial and vertical peak particle velocity with respect to the horizontal distance and elevation difference are\(PP{V_X}={\text{ }}3.29 \cdot {R^{-0.057}} \cdot {H^{-0.285}}\) and \(PP{V_Z}={\text{ }}4.19 \cdot {R^{ - 0.215}} \cdot {H^{ - 0.215}}\), respectively.

  4. 4.

    For slope terrain, in terms of the variation law of peak particle velocity inside the same slope body; at the same horizontal distance in the slope, the elevation has an influence on the peak particle velocity of slope blasting. For horizontal radial direction, the peak particle velocity in the near field of CO2 blasting shows an amplification trend with the increase of elevation, showing that the peak particle velocity of monitoring point near the slope surface is larger; the peak particle velocity in the far field of CO2 blasting shows a decrease trend with the increase of elevation, showing that the peak particle velocity of monitoring point near the slope surface is smaller. For vertical direction, the peak particle velocity shows an amplification trend in the near field and far field of CO2 blasting with the increase of elevation, showing that the peak particle velocity of monitoring point near the slope surface is larger.

  5. 5.

    For slope terrain, in terms of the variation law of peak particle velocity on slope surface in the different slope; in the horizontal radial direction, for horizontal distances of 5 m, 10 m, 20 m, 30 m, and 40 m, as the slope gradient increases from 15 ° to 60 °, the peak vibration velocity of the slope surface shows a decreasing trend; in the vertical direction, for horizontal distances of 5 m, 10 m, and 40 m, as the slope gradient increases from 15 ° to 60 °, the peak vibration velocity of the slope surface shows a decreasing trend; for horizontal distances of 20 m and 30 m, as the slope gradient increases from 15 ° to 60 °, the peak vibration velocity of the slope surface shows a trend of first decreasing, then increasing, and then decreasing again.

Although this article has successfully used the discrete element software UDEC to establish numerical models for carbon dioxide phase change fracturing blasting in flat terrain, concave terrain, convex terrain, and slope terrain, and obtained the propagation law of liquid carbon dioxide phase change induced blasting vibration waves in different terrains, some research results have been achieved. However, due to the complexity of the carbon dioxide fracturing blasting process and numerous influencing factors, the research work in this article is still in the exploratory stage and needs further in-depth exploration. The shortcomings and further work of the paper are: (1) Further conduct research on the damage and vibration effects of carbon dioxide fracturing rock masses by using various means such as acoustic emission monitoring, CT scanning, acoustic velocity measurement, and vibration monitoring. (2) Carbon dioxide fracturing blasting is caused by stress waves and high-pressure gas combined. In numerical simulations, load simplification is generally only applied to stress waves without considering the effect of high-pressure gas. To better understand rock mass’s damage and failure mechanism, the ABAQUS and ANSYS LS-DYNA3D can be used for dynamic load and fluid-structure coupling analysis of stress waves, high-pressure gas, and rock masses.