Rock Breaking and Dynamic Response Characteristics of Carbon Dioxide Phase Transition Fracturing Considering the Gathering Energy E ﬀ ect

: Carbon dioxide phase transition fracturing has been widely used in rock mass excavation under complex environments, and its special rock breaking process shows obvious gathering energy effect. In this paper, the gathering energy effect of this technology is considered and then the impact reduction coefficient is defined and determined. Eventually, a combined method of field tests and numerical simulations is used to study the crack propagation characteristics and spatiotemporal changes of dynamic response. The results show that the cracks grow more and more slowly as time goes by; the peak displacement and peak point velocity in the primary impact direction are both greater than those in the secondary impact direction. The peak point velocity in different directions decreases as the distance from borehole increases and it decays more and more slowly. With the increase of distance from the borehole, the peak effective stress in the primary impact direction constantly decreases. However, it increases first and then decreases in the secondary impact direction. The results mentioned above can provide effective guidance for later experimental research and engineering.


Introduction
With the continuous progress of underground space technology, the requirements for construction safety and environmental safety are becoming stricter. Explosive blasting is limited to being used in insensitive areas because of the shortcomings of large exothermic reactions, large vibration and high overpressure. Carbon dioxide phase transition fracturing (CDPTF) is extensively used as an alternative to explosive blasting when performing the gas drainage in a coal mine or excavating rock mass in sensitive areas such as hospitals, residential areas, schools and gas stations [1][2][3][4][5][6]. This technology utilizes the phase energy difference of the carbon dioxide phase physical transformation process to break the rock. It has the advantages of no sparks and low vibration. In addition, it can absorb heat and suppress fire during the whole process. In short, CDPTF is a safe, green and reliable rock breaking technology [6][7][8].
At present, theoretical research of CDPTF lags behind the practical level. The carbon dioxide phase transition load is the focus of theoretical research. Until now, numerous breakthroughs have been made in the CDPTF loading characteristics by the theoretical derivation and experimental test, and a variety of methods have been proposed to obtain the fracturing load. Among them, the most widely used method is the TNT equivalent method with fewer calculation parameters. It is recognized that the CDPTFT energy is linked to the storage pipe volume, rupture disc strength and external pressure of the fracturing pipe in this method [9][10][11]. In addition, Sun et al. [12] tested the pressure pressure curve of gas explosion nozzle under the specific temperature and pressure, selected the MATLAB to fit the Jones-Wilkins-Lee (JWL) parameters in the corresponding case and used the JWL equation to reflect the CDPTFT load; Guo [13] combined the root-mean-square (RMS) analysis of the vibration data and energy analysis to quantify the total energy of rock breaking; Zhou et al. [14] and Ke et al. [15] measured the pressure in the fracturing pipe, analyzed the pressure change law and gave the total energy calculation method. Kang et al. [16] studied the gas pressure change in the whole fracturing process by the improved van der Waals equation and obtained the total fracturing energy; Xie et al. [17] utilized the changes in carbon dioxide entropy and enthalpy before and after the reaction, combined with the burst energy release formula of the superheated liquefied gas container and obtained the function between the shear pressure and the total energy.
As we mentioned above, due to the lack of understanding for the CDPTF characteristics and because the blasting is a mature technology, scholars converted the total energy of CDPTF into TNT equivalents to simulate the CDPTF [11]. As is known to all, the explosive pressure in the same borehole section is equal, and explosive cracks grow in all directions around the borehole. The cracks simulated by this TNT equivalent method also have such characteristics. In fact, due to the special mechanical structure of the fracturing pipe, the energy of CDPTF mainly gathers in one direction [18]. Zhang [19] verified this gathering energy effect of the CDPTF by using laboratory concrete tests and applied the loads only to the borehole wall directly corresponding to the gas outlet to conduct the numerical simulation. However, under the high-pressure impact, even if the load of the borehole wall not corresponding to the gas outlet is significantly smaller than that of the borehole wall corresponding to the gas outlet in actual engineering, the borehole wall not corresponding to the gas outlet is inevitably bearing the impact load during the CDPTF. So far, there is no one method that can perfectly reflect the gathering energy effect in actual engineering.
In this paper, we fully considered the load of the borehole wall which did not correspond to the gas outlet, defined the impact reduction coefficient, took the crack length as the comparison standard, and determined the most reasonable impact reduction coefficient by combining field tests and numerical simulation results. In addition, based on the most reasonable impact reduction coefficient, rock breaking and dynamic response characteristics of CDPTF were obtained through numerical simulations. The load application method based on the impact reduction coefficient provides a load method that can quantitatively reflect the rock breaking in the actual engineering. The crack propagation and dynamic response characteristics we obtained can provide a reference for the rock breaking mechanism of CDPTF.

Procedure of Carbon Dioxide Phase Transition Fracturing
The core device of the CDPTF is the carbon dioxide phase transition fracturing pipe, which is composed of a detonating head, a heater, a liquid storage pipe, a rupture disc, a discharge head and two sealing gaskets, as shown in Figure 1 [20]. Carbon dioxide is the main medium for CDPTF. As Figure 2 shows, when the temperature is higher than 31.4 °C and the pressure exceeds 73.8 bar, the interface between liquid carbon dioxide and gaseous carbon dioxide will suddenly disappear and all the carbon dioxide reaches a special supercritical state. Supercritical carbon dioxide has some properties of both gas and liquid phases. Its density is close to that of liquid and its diffusion coefficient is close to that of gas. It is because of this Carbon dioxide is the main medium for CDPTF. As Figure 2 shows, when the temperature is higher than 31.4 • C and the pressure exceeds 73.8 bar, the interface between liquid carbon dioxide and gaseous carbon dioxide will suddenly disappear and all the carbon dioxide reaches a special supercritical state. Supercritical carbon dioxide has some properties of both gas and liquid phases. Its density is close to Energies 2020, 13, 1336 3 of 16 that of liquid and its diffusion coefficient is close to that of gas. It is because of this nature that it can be quickly heated in the liquid storage pipe, and the pressure in the storage pipe can rise rapidly. The energy difference between supercritical carbon dioxide and gaseous carbon dioxide is utilized to break the rock. During the CDPTF, firstly the liquid carbon dioxide is filled into the liquid storage pipe, and then the filled pipe is sent into the pre-drilled borehole. The heater is energized to generate a large amount of heat immediately after the fracturing pipe is installed. After the heat is absorbed, liquid carbon dioxide turns into a supercritical state, and the internal pressure of the liquid storage pipe continues to rise. When the internal pressure exceeds the pre-set pressure of the rupture disc, the rupture disc is destroyed. At the same time, the supercritical carbon dioxide is instantly released and turns to high-pressure gas. The high-pressure gas brings an extremely strong impact immediately, which damages the rocks around the borehole, resulting in a large number of cracks. After the high-pressure gas enters the impact cracks, the cracks are affected by the gas wedge and continue to grow. The whole procedure of CDPTF is illustrated in Figure 3 [21][22][23].
Energies 2020, 13, x FOR PEER REVIEW 3 of 16 nature that it can be quickly heated in the liquid storage pipe, and the pressure in the storage pipe can rise rapidly. The energy difference between supercritical carbon dioxide and gaseous carbon dioxide is utilized to break the rock. During the CDPTF, firstly the liquid carbon dioxide is filled into the liquid storage pipe, and then the filled pipe is sent into the pre-drilled borehole. The heater is energized to generate a large amount of heat immediately after the fracturing pipe is installed. After the heat is absorbed, liquid carbon dioxide turns into a supercritical state, and the internal pressure of the liquid storage pipe continues to rise. When the internal pressure exceeds the pre-set pressure of the rupture disc, the rupture disc is destroyed. At the same time, the supercritical carbon dioxide is instantly released and turns to high-pressure gas. The high-pressure gas brings an extremely strong impact immediately, which damages the rocks around the borehole, resulting in a large number of cracks. After the high-pressure gas enters the impact cracks, the cracks are affected by the gas wedge and continue to grow. The whole procedure of CDPTF is illustrated in Figure 3 [21][22][23].

Gathering Energy Characteristics of Carbon Dioxide Phase Transition Fracturing
The energy release outlet controls the energy to gather in a specific direction when the rock mass is broken by CDPTF. Just after the energy is released, the high-energy carbon dioxide impacts the borehole, and the tangential tensile stress generates inside the rock mass. When the dynamic tensile stress exceeds dynamic tensile strength, the rock mass is damaged and a radial crack generates. The crack propagation direction is perpendicular to the tensile stress direction [19]. Not only does the gathering energy effect promote the formation of tensile cracks perpendicular to the energy release direction, but also it reduces the crushing towards the adjacent rock by the impact. Furthermore, the energy utilization rate is increased. Until now, CDPTF research considering the gathering energy nature that it can be quickly heated in the liquid storage pipe, and the pressure in the storage pipe can rise rapidly. The energy difference between supercritical carbon dioxide and gaseous carbon dioxide is utilized to break the rock. During the CDPTF, firstly the liquid carbon dioxide is filled into the liquid storage pipe, and then the filled pipe is sent into the pre-drilled borehole. The heater is energized to generate a large amount of heat immediately after the fracturing pipe is installed. After the heat is absorbed, liquid carbon dioxide turns into a supercritical state, and the internal pressure of the liquid storage pipe continues to rise. When the internal pressure exceeds the pre-set pressure of the rupture disc, the rupture disc is destroyed. At the same time, the supercritical carbon dioxide is instantly released and turns to high-pressure gas. The high-pressure gas brings an extremely strong impact immediately, which damages the rocks around the borehole, resulting in a large number of cracks. After the high-pressure gas enters the impact cracks, the cracks are affected by the gas wedge and continue to grow. The whole procedure of CDPTF is illustrated in Figure 3 [21][22][23].

Gathering Energy Characteristics of Carbon Dioxide Phase Transition Fracturing
The energy release outlet controls the energy to gather in a specific direction when the rock mass is broken by CDPTF. Just after the energy is released, the high-energy carbon dioxide impacts the borehole, and the tangential tensile stress generates inside the rock mass. When the dynamic tensile stress exceeds dynamic tensile strength, the rock mass is damaged and a radial crack generates. The crack propagation direction is perpendicular to the tensile stress direction [19]. Not only does the gathering energy effect promote the formation of tensile cracks perpendicular to the energy release direction, but also it reduces the crushing towards the adjacent rock by the impact. Furthermore, the energy utilization rate is increased. Until now, CDPTF research considering the gathering energy

Gathering Energy Characteristics of Carbon Dioxide Phase Transition Fracturing
The energy release outlet controls the energy to gather in a specific direction when the rock mass is broken by CDPTF. Just after the energy is released, the high-energy carbon dioxide impacts the borehole, and the tangential tensile stress generates inside the rock mass. When the dynamic tensile stress exceeds dynamic tensile strength, the rock mass is damaged and a radial crack generates. The crack propagation direction is perpendicular to the tensile stress direction [19]. Not only does the gathering energy effect promote the formation of tensile cracks perpendicular to the energy release direction, but also it reduces Energies 2020, 13, 1336 4 of 16 the crushing towards the adjacent rock by the impact. Furthermore, the energy utilization rate is increased. Until now, CDPTF research considering the gathering energy effect is still in the primary stage of the development. In order to quantify the gathering energy effect on the crack propagation, the impact reduction coefficient k is introduced to describe the CDPTF load in this paper.
It is assumed that the maximum radial compressive stress in the direction along the gas outlet (primary impact direction) during rock breaking is P m , and the maximum radial compressive stress in the direction perpendicular to the gas outlet (secondary impact direction) is: where P s is the peak load on the borehole in the secondary impact direction, MPa; k is the impact reduction coefficient, which is the ratio of the peak pressure in the secondary impact direction to the peak pressure in the primary impact direction. Equation (1) representss the load in the secondary impact direction, reflects the difference in load in the primary and secondary impact direction, and quantitatively represents the energy accumulation effect of CDPTF.
In the existing theory, when the compressive stress is applied to the borehole, the borehole wall whose direction is vertical to the compressive stress generates tensile stress correspondingly, which can be expressed as [23]: where T is the tension stress, MPa; µ is the Poisson's ratio of the rock mass; P 0 is the compressive stress on the borehole, MPa. When the tensile stress exceeds the dynamic tensile strength of the rock mass, the rock mass is damaged and radial tensile cracks generate. Since the peak compressive stress in the primary impact direction is greater than that of the secondary impact direction, the borehole wall in the secondary impact direction is subjected to greater tensile stress. Therefore, these cracks are mainly distributed along the secondary impact direction [19]. The peak pressure can be obtained by testing the mechanical properties of the rupture disc. In the elastic region of the adjacent rock mass, the radial peak pressure attenuation equation is [24]: where σ r is radial peak compressive stress in the rock mass, MPa; a is the borehole radius, m; r is the length of the line connecting the calculation point with the center of the borehole, m; α is the stress attenuation coefficient, α= 2 − µ/1 − µ.
In the process of shock wave propagation, the relationship between tangential stress and radial stress around the borehole is [25]: where σ θ is the tangential stress, MPa; b is the ratio of the S-wave velocity to the P-wave velocity of the rock. When the tangential stress is greater than the dynamic tensile strength of the rock, the rock is damaged by tension and the cracks grow forward under the action of shock waves.

Field Tests of Carbon Dioxide Phase Transition Fracturing
Since the existing research has not yet covered the distribution law of crack length in different directions to CDPTF, field tests and numerical simulations are used to determine the most reasonable impact reduction coefficient for CDPTF in this paper. In order to provide the experimental basis for numerical simulation, we conducted multiple single-hole fracturing tests to obtain the rock breaking effect of CDPTF in actual engineering. The field test was conducted in Nyingchi, Tibet. The lithology of the tested rock mass was feldspar quartz sandstone. The thickness of the sandstone layer was greater than 3 m. The rock mass was complete and unweathered. There were no unstable blocks and the structural plane was not developed in the field test area. The rock quality designation (RQD) value of the rock mass was 96%, which belongs to a good quality rock mass. According to the rock mass ratings (RMR) classification method, the quality of the rock mass was classified as the first level. Common type 73 fracturing pipes were used in the test; their outer diameter was 73 mm. The diameter of the energy release outlet was 4 mm, and the pre-set pressure of the rupture disc was 140 MPa. There were two symmetrical energy release outlets on the discharge head, as shown in Figure 4. Before the test, boreholes with a diameter of 100 mm and a depth of 3 m were drilled. After half of the borehole depth was blocked with yellow mud and crushed stone, the filled fracturing pipe was inserted into the borehole. Then, yellow mud and crushed stone were utilized to seal the borehole in order to prevent the pipe from moving outward during the fracturing process. It is worth explaining that the selected yellow mud was taken directly from the test site, and the crushed stone was about 1 cm in diameter. The proportion of crushed stones in the sealing material was about 70%. The HTSW-500D high-energy pulse detonator was selected to remotely control the heater in the fracturing pipe to generate heat and then the fracturing pipe was detonated to break the rock mass. Three single-hole CDPTF experiments were carried out so as to ensure that the crack propagation can be accurately obtained. The rock breaking effect of CDPTF is shown in Figure 5. After the test, we observed that the cracks are mainly distributed in a direction perpendicular to the direction of the gas outlet. In addition, the widths of the main cracks on the borehole walls were about 3-4 cm, and the main crack lengths of the three experiments were counted respectively as 1.59 m, 1.56 m, and 1.45 m.
Energies 2020, 13, x FOR PEER REVIEW 5 of 16 which belongs to a good quality rock mass. According to the rock mass ratings (RMR) classification method, the quality of the rock mass was classified as the first level. Common type 73 fracturing pipes were used in the test; their outer diameter was 73 mm. The diameter of the energy release outlet was 4 mm, and the pre-set pressure of the rupture disc was 140 MPa. There were two symmetrical energy release outlets on the discharge head, as shown in Figure 4. Before the test, boreholes with a diameter of 100 mm and a depth of 3 m were drilled. After half of the borehole depth was blocked with yellow mud and crushed stone, the filled fracturing pipe was inserted into the borehole. Then, yellow mud and crushed stone were utilized to seal the borehole in order to prevent the pipe from moving outward during the fracturing process. It is worth explaining that the selected yellow mud was taken directly from the test site, and the crushed stone was about 1 cm in diameter. The proportion of crushed stones in the sealing material was about 70%. The HTSW-500D high-energy pulse detonator was selected to remotely control the heater in the fracturing pipe to generate heat and then the fracturing pipe was detonated to break the rock mass. Three single-hole CDPTF experiments were carried out so as to ensure that the crack propagation can be accurately obtained. The rock breaking effect of CDPTF is shown in Figure 5. After the test, we observed that the cracks are mainly distributed in a direction perpendicular to the direction of the gas outlet. In addition, the widths of the main cracks on the borehole walls were about 3-4 cm, and the main crack lengths of the three experiments were counted respectively as 1.59 m, 1.56 m, and 1.45 m.

Acquisition of Rock Mass Mechanical Parameters
We collected on-site rock samples and performed laboratory mechanical tests after processing. First, the rock samples were ground and processed, and the thin section identification was taken under the polarizing microscope of China University of Geosciences. According to the results, the rock sample contains 20% feldspar, 20% quartz, 40% calcite, and 10% amphibole (as shown in Figure  6), which belongs to feldspar quartz sandstone. Then we took a whole rock sample, used the direct drilled method to get a lot of Φ100 mm × 50 mm cylindrical samples and Φ50 mm × 25 mm disk samples. The samples had no obvious defects on the surface and similar P-wave velocities were selected to conduct the uniaxial compression tests (as is shown in Figure 7) and Brazil split tests to which belongs to a good quality rock mass. According to the rock mass ratings (RMR) classification method, the quality of the rock mass was classified as the first level. Common type 73 fracturing pipes were used in the test; their outer diameter was 73 mm. The diameter of the energy release outlet was 4 mm, and the pre-set pressure of the rupture disc was 140 MPa. There were two symmetrical energy release outlets on the discharge head, as shown in Figure 4. Before the test, boreholes with a diameter of 100 mm and a depth of 3 m were drilled. After half of the borehole depth was blocked with yellow mud and crushed stone, the filled fracturing pipe was inserted into the borehole. Then, yellow mud and crushed stone were utilized to seal the borehole in order to prevent the pipe from moving outward during the fracturing process. It is worth explaining that the selected yellow mud was taken directly from the test site, and the crushed stone was about 1 cm in diameter. The proportion of crushed stones in the sealing material was about 70%. The HTSW-500D high-energy pulse detonator was selected to remotely control the heater in the fracturing pipe to generate heat and then the fracturing pipe was detonated to break the rock mass. Three single-hole CDPTF experiments were carried out so as to ensure that the crack propagation can be accurately obtained. The rock breaking effect of CDPTF is shown in Figure 5. After the test, we observed that the cracks are mainly distributed in a direction perpendicular to the direction of the gas outlet. In addition, the widths of the main cracks on the borehole walls were about 3-4 cm, and the main crack lengths of the three experiments were counted respectively as 1.59 m, 1.56 m, and 1.45 m.

Acquisition of Rock Mass Mechanical Parameters
We collected on-site rock samples and performed laboratory mechanical tests after processing. First, the rock samples were ground and processed, and the thin section identification was taken under the polarizing microscope of China University of Geosciences. According to the results, the rock sample contains 20% feldspar, 20% quartz, 40% calcite, and 10% amphibole (as shown in Figure  6), which belongs to feldspar quartz sandstone. Then we took a whole rock sample, used the direct drilled method to get a lot of Φ100 mm × 50 mm cylindrical samples and Φ50 mm × 25 mm disk samples. The samples had no obvious defects on the surface and similar P-wave velocities were selected to conduct the uniaxial compression tests (as is shown in Figure 7) and Brazil split tests to

Acquisition of Rock Mass Mechanical Parameters
We collected on-site rock samples and performed laboratory mechanical tests after processing. First, the rock samples were ground and processed, and the thin section identification was taken under the polarizing microscope of China University of Geosciences. According to the results, the rock sample contains 20% feldspar, 20% quartz, 40% calcite, and 10% amphibole (as shown in Figure 6), which Energies 2020, 13, 1336 6 of 16 belongs to feldspar quartz sandstone. Then we took a whole rock sample, used the direct drilled method to get a lot of Φ100 mm × 50 mm cylindrical samples and Φ50 mm × 25 mm disk samples. The samples had no obvious defects on the surface and similar P-wave velocities were selected to conduct the uniaxial compression tests (as is shown in Figure 7) and Brazil split tests to obtain the mechanical parameters. Moreover, the physical parameters of the rock were also tested. All results are presented in Table 1 below.
Energies 2020, 13, x FOR PEER REVIEW 6 of 16 obtain the mechanical parameters. Moreover, the physical parameters of the rock were also tested. All results are presented in Table 1 below.

Numerical Simulation of Carbon Dioxide Phase Transition Fracturing
In 2018, Lei [18] proposed that the loads at different positions in the same cross section are not equal in the CDPTF. Shortly afterward Zhang et al. [19] damaged the concrete by CDPTF and performed numerical simulations, but he only applied the load to the borehole wall in the primary impact direction, without considering the load in the secondary impact direction. In fact, since the borehole and the fracturing pipe are not coupled, the borehole wall in the secondary impact direction will inevitably be impacted. It is not certain for the relationship between the impact load in the primary impact direction and that in the secondary impact direction. In order to determine the most reasonable impact reduction coefficient for this test, the non-linear dynamics software LS-DYNA was used to investigate the crack propagation. It is worth mentioning that LS-DYNA is a famous finite element analysis program, and it is also the recognized originator of explicit integral calculation programs and the best explicit analysis software. The Lagrange algorithm and the explicit solution are most widely used in LS-DYNA. This software is suitable for high-speed collisions, explosions and material forming. The crack growth length is the main parameter for evaluating the effect of CDPTF, so we use the crack length as the basis for the comparison between numerical simulations and field tests. The most reasonable impact reduction coefficient was obtained by comparing the crack length  obtain the mechanical parameters. Moreover, the physical parameters of the rock were also tested. All results are presented in Table 1 below.

Numerical Simulation of Carbon Dioxide Phase Transition Fracturing
In 2018, Lei [18] proposed that the loads at different positions in the same cross section are not equal in the CDPTF. Shortly afterward Zhang et al. [19] damaged the concrete by CDPTF and performed numerical simulations, but he only applied the load to the borehole wall in the primary impact direction, without considering the load in the secondary impact direction. In fact, since the borehole and the fracturing pipe are not coupled, the borehole wall in the secondary impact direction will inevitably be impacted. It is not certain for the relationship between the impact load in the primary impact direction and that in the secondary impact direction. In order to determine the most reasonable impact reduction coefficient for this test, the non-linear dynamics software LS-DYNA was used to investigate the crack propagation. It is worth mentioning that LS-DYNA is a famous finite element analysis program, and it is also the recognized originator of explicit integral calculation programs and the best explicit analysis software. The Lagrange algorithm and the explicit solution are most widely used in LS-DYNA. This software is suitable for high-speed collisions, explosions and material forming. The crack growth length is the main parameter for evaluating the effect of CDPTF, so we use the crack length as the basis for the comparison between numerical simulations and field tests. The most reasonable impact reduction coefficient was obtained by comparing the crack length

Numerical Simulation of Carbon Dioxide Phase Transition Fracturing
In 2018, Lei [18] proposed that the loads at different positions in the same cross section are not equal in the CDPTF. Shortly afterward Zhang et al. [19] damaged the concrete by CDPTF and performed numerical simulations, but he only applied the load to the borehole wall in the primary impact direction, without considering the load in the secondary impact direction. In fact, since the borehole and the fracturing pipe are not coupled, the borehole wall in the secondary impact direction will inevitably be impacted. It is not certain for the relationship between the impact load in the primary impact direction and that in the secondary impact direction. In order to determine the most reasonable impact reduction coefficient for this test, the non-linear dynamics software LS-DYNA was used to investigate the crack propagation. It is worth mentioning that LS-DYNA is a famous finite element analysis program, and it is also the recognized originator of explicit integral calculation programs and the best explicit analysis software. The Lagrange algorithm and the explicit solution are most widely used in LS-DYNA. This software is suitable for high-speed collisions, explosions and material forming. The crack growth length is the main parameter for evaluating the effect of CDPTF, so we use the Energies 2020, 13, 1336 7 of 16 crack length as the basis for the comparison between numerical simulations and field tests. The most reasonable impact reduction coefficient was obtained by comparing the crack length of numerical simulations whose impact reduction coefficients were between 0.5 and 1 with the field test results.

Establishment of the Numerical Model
A 1/4 numerical model with a side length of 4.0 m × 4.0 m was established. The blast hole was located at the lower left of the model with a radius of 5 cm. Because the whole model contained only one material-rock, the Lagrangian algorithm and Soild164 element were selected. In addition, it can be seen from Section 3.1 that the tested rock mass had good integrity and belonged to the first-level rock mass, so the rock mass in the model was a continuous medium material. A method combining free mesh and mapped mesh was employed in this model. The free mesh was used in the vicinity of the blast hole to avoid the mesh distortion caused by the mapped mesh from controlling the crack propagation direction. Since this 1/4 model simulated the crack propagation in infinite rock mass caused by CDPTF, constraints were imposed on the symmetry axis. A nonreflecting boundary condition was applied to the model in order to avoid stress wave reflection at the boundary. The numerical calculation model is shown in Figure 8.
Energies 2020, 13, x FOR PEER REVIEW 7 of 16 of numerical simulations whose impact reduction coefficients were between 0.5 and 1 with the field test results.

Establishment of the Numerical Model
A 1/4 numerical model with a side length of 4.0 m × 4.0 m was established. The blast hole was located at the lower left of the model with a radius of 5 cm. Because the whole model contained only one material-rock, the Lagrangian algorithm and Soild164 element were selected. In addition, it can be seen from Section 3.1 that the tested rock mass had good integrity and belonged to the first-level rock mass, so the rock mass in the model was a continuous medium material. A method combining free mesh and mapped mesh was employed in this model. The free mesh was used in the vicinity of the blast hole to avoid the mesh distortion caused by the mapped mesh from controlling the crack propagation direction. Since this 1/4 model simulated the crack propagation in infinite rock mass caused by CDPTF, constraints were imposed on the symmetry axis. A nonreflecting boundary condition was applied to the model in order to avoid stress wave reflection at the boundary. The numerical calculation model is shown in Figure 8.

Selection of Load and Model Parameter
The CDPTF load is special, and its pressure curve is between that of blasting and hydraulic fracturing, and it is usually converted to the explosive load in the TNT equivalent method. There are obvious differences between the numerical simulation results and the actual rock breaking effect. According to the existing laboratory tests, the pressure increase rate in the gas ejection direction is about 200 GPa/s, and the gas wedge duration is about 5 ms while the high-energy gas impacts the borehole [18,19,26]. Moreover, the pressure increase rate in the secondary impact direction is significantly smaller than that in the primary impact direction. There is no obvious difference for the decrease section of load between the primary impact direction and the secondary impact direction, because the main load manifestation of two directions in the descending section is the static wedge effect of high-pressure gas. A triangular equivalent loading method commonly used in engineering blasting was used in this model [27], and the loading curve is shown in Figure 9. It is assumed that when the high-energy gas bursts out of the energy release outlet, the loading area of the borehole in the primary impact direction corresponds to the area of the energy release outlets [28]. In Figure 10, the red arrows and the blue arrows represent the impact load of the primary impact direction and the secondary impact direction, respectively.

Selection of Load and Model Parameter
The CDPTF load is special, and its pressure curve is between that of blasting and hydraulic fracturing, and it is usually converted to the explosive load in the TNT equivalent method. There are obvious differences between the numerical simulation results and the actual rock breaking effect. According to the existing laboratory tests, the pressure increase rate in the gas ejection direction is about 200 GPa/s, and the gas wedge duration is about 5 ms while the high-energy gas impacts the borehole [18,19,26]. Moreover, the pressure increase rate in the secondary impact direction is significantly smaller than that in the primary impact direction. There is no obvious difference for the decrease section of load between the primary impact direction and the secondary impact direction, because the main load manifestation of two directions in the descending section is the static wedge effect of high-pressure gas. A triangular equivalent loading method commonly used in engineering blasting was used in this model [27], and the loading curve is shown in Figure 9. It is assumed that when the high-energy gas bursts out of the energy release outlet, the loading area of the borehole in the primary impact direction corresponds to the area of the energy release outlets [28]. In Figure 10, the red arrows and the blue arrows represent the impact load of the primary impact direction and the secondary impact direction, respectively.
The most used constitutive model in LS-DYNA is the MAT_PLASTIC_KINEMATIC model. It can describe isotropic hardening materials and follow-hardened plastic materials. The Cowper-Symonds model used in this constitutive model considers the strain rate effect. The MAT_PLAST-IC_KINEMATIC Energies 2020, 13, 1336 8 of 16 model is applicable to beams, shells and solid elements [29]. When this constitutive model is used for numerical simulations, the calculation efficiency is also very high, so this model is also selected to simulate rocks in this numerical simulation. The model parameters are shown in Table 2 below.  The most used constitutive model in LS-DYNA is the MAT_PLASTIC_KINEMATIC model. It can describe isotropic hardening materials and follow-hardened plastic materials. The Cowper-Symonds model used in this constitutive model considers the strain rate effect. The MAT_PLAST-IC_KINEMATIC model is applicable to beams, shells and solid elements [29]. When this constitutive model is used for numerical simulations, the calculation efficiency is also very high, so this model is also selected to simulate rocks in this numerical simulation. The model parameters are shown in Table  2 below. In addition, the keyword MAT_ADD_EROSION was added to delete the failed units, which was convenient to obtain the crack growth length. The CDPTF had the crushed zones, which is similar to the blasting. However, it was too small and the rock mass was mainly damaged by tensile stress. Thus, the tensile failure criterion named mneps in MAT_ADD_EROSION was used.   The most used constitutive model in LS-DYNA is the MAT_PLASTIC_KINEMATIC model. It can describe isotropic hardening materials and follow-hardened plastic materials. The Cowper-Symonds model used in this constitutive model considers the strain rate effect. The MAT_PLAST-IC_KINEMATIC model is applicable to beams, shells and solid elements [29]. When this constitutive model is used for numerical simulations, the calculation efficiency is also very high, so this model is also selected to simulate rocks in this numerical simulation. The model parameters are shown in Table  2 below. In addition, the keyword MAT_ADD_EROSION was added to delete the failed units, which was convenient to obtain the crack growth length. The CDPTF had the crushed zones, which is similar to the blasting. However, it was too small and the rock mass was mainly damaged by tensile stress. Thus, the tensile failure criterion named mneps in MAT_ADD_EROSION was used.  In addition, the keyword MAT_ADD_EROSION was added to delete the failed units, which was convenient to obtain the crack growth length. The CDPTF had the crushed zones, which is similar to the blasting. However, it was too small and the rock mass was mainly damaged by tensile stress. Thus, the tensile failure criterion named mneps in MAT_ADD_EROSION was used.

Determination of Impact Reduction Coefficient
A set of numerical models to impact reduction coefficients between 0.5 and 1 was calculated and the crack propagation length of each model was counted. The numerical simulation results of some models are shown in Figure 11. Because of the different impact reduction coefficients, the crack growth Energies 2020, 13, 1336 9 of 16 length of each model was significantly different. However, the difference in the angle between the crack and the secondary impact direction of all models was very small, ranging from 13 to 16 degrees. When the impact reduction coefficient was less than 0.9, the crack growth length increased with the increase of impact reduction coefficient. All the statistical data of the crack length are shown in Figure 12. As the impact reduction coefficient exceeded 0.9, the impact difference between the two impact directions was small. The tensile stress in the secondary impact direction was not enough to damage the rock, and the crack could not generate.
Energies 2020, 13, x FOR PEER REVIEW 9 of 16 A set of numerical models to impact reduction coefficients between 0.5 and 1 was calculated and the crack propagation length of each model was counted. The numerical simulation results of some models are shown in Figure 11. Because of the different impact reduction coefficients, the crack growth length of each model was significantly different. However, the difference in the angle between the crack and the secondary impact direction of all models was very small, ranging from 13 to 16 degrees. When the impact reduction coefficient was less than 0.9, the crack growth length increased with the increase of impact reduction coefficient. All the statistical data of the crack length are shown in Figure 12. As the impact reduction coefficient exceeded 0.9, the impact difference between the two impact directions was small. The tensile stress in the secondary impact direction was not enough to damage the rock, and the crack could not generate.   Energies 2020, 13, x FOR PEER REVIEW 9 of 16 A set of numerical models to impact reduction coefficients between 0.5 and 1 was calculated and the crack propagation length of each model was counted. The numerical simulation results of some models are shown in Figure 11. Because of the different impact reduction coefficients, the crack growth length of each model was significantly different. However, the difference in the angle between the crack and the secondary impact direction of all models was very small, ranging from 13 to 16 degrees. When the impact reduction coefficient was less than 0.9, the crack growth length increased with the increase of impact reduction coefficient. All the statistical data of the crack length are shown in Figure 12. As the impact reduction coefficient exceeded 0.9, the impact difference between the two impact directions was small. The tensile stress in the secondary impact direction was not enough to damage the rock, and the crack could not generate.   Fitting the data points in Figure 12, the relationship between the impact reduction coefficient and the crack growth length is obtained as follows: where y is the crack length, cm; x is the impact reduction coefficient. After knowing the crack length, combined with Equation (5), we can infer the most reasonable impact reduction coefficient in this case. In addition, the correlation coefficient of Equation (5) is 0.9707, so the most reasonable impact reduction coefficient calculated according to this equation has high accuracy and low uncertainty.

Crack Propagation Characteristics
According to the field test results and Equation (5), it can be recognized that when the impact reduction coefficient is 0.7, the numerical simulation can accurately reflect the field test results of CDPTF. Therefore, the model whose impact reduction coefficient was 0.7 was selected to analyze the crack propagation characteristics. The crack propagation in this condition is shown in Figure 13. From Figure 13, we can find that cracks of CDPTF were mainly X-shaped along the secondary impact direction in good quality rock mass, these cracks were generally flat. Moreover, no cracks generate in the primary impact direction. The results are consistent with those of Zhang's experiments. The boundary conditions of these numerical simulations were determined according to the field test, and the model parameters are obtained according to the laboratory test. The results were almost identical to the field test, which proves that our model is effective. The CDPTF fracture zone had a dominant direction, which was different from the explosive fracture zone.
Energies 2020, 13, x FOR PEER REVIEW 10 of 16 Fitting the data points in Figure 12, the relationship between the impact reduction coefficient and the crack growth length is obtained as follows: = 268.6 − 36.06 ( 2 = 0.9707) Where y is the crack length, cm; x is the impact reduction coefficient. After knowing the crack length, combined with Equation (5), we can infer the most reasonable impact reduction coefficient in this case. In addition, the correlation coefficient of Equation (5) is 0.9707, so the most reasonable impact reduction coefficient calculated according to this equation has high accuracy and low uncertainty.

Crack Propagation Characteristics
According to the field test results and Equation (5), it can be recognized that when the impact reduction coefficient is 0.7, the numerical simulation can accurately reflect the field test results of CDPTF. Therefore, the model whose impact reduction coefficient was 0.7 was selected to analyze the crack propagation characteristics. The crack propagation in this condition is shown in Figure 13. From Figure 13, we can find that cracks of CDPTF were mainly X-shaped along the secondary impact direction in good quality rock mass, these cracks were generally flat. Moreover, no cracks generate in the primary impact direction. The results are consistent with those of Zhang's experiments. The boundary conditions of these numerical simulations were determined according to the field test, and the model parameters are obtained according to the laboratory test. The results were almost identical to the field test, which proves that our model is effective. The CDPTF fracture zone had a dominant direction, which was different from the explosive fracture zone. Crack propagation is an important index for evaluating the rock breaking effect. Therefore, the crack length of CDPTF at each moment was calculated, and a scatter plot of the crack length with time is drawn in Figure 14. The crack tip continuously accumulated fracture energy so that the crack stopped moving forward during 0.8-0.84 ms, 1.28-1.4 ms and 1.6-2.52 ms. The longer the crack was, the more time it took for the crack to gather fracture energy, and the more difficult it was for a crack to re-initiate. Fitting the data in Figure 14, we can get the function between the crack length and time as follows: Where y is the crack length, cm; t is the time of CDPTF in μs.
As is shown in Figure 14, the crack length generally increased more and more slowly, indicating that the crack propagation velocity decreased continuously with time. Crack propagation is an important index for evaluating the rock breaking effect. Therefore, the crack length of CDPTF at each moment was calculated, and a scatter plot of the crack length with time is drawn in Figure 14. The crack tip continuously accumulated fracture energy so that the crack stopped moving forward during 0.8-0.84 ms, 1.28-1.4 ms and 1.6-2.52 ms. The longer the crack was, the more time it took for the crack to gather fracture energy, and the more difficult it was for a crack to re-initiate. Fitting the data in Figure 14, we can get the function between the crack length and time as follows: y = 54.213 ln(t − 193.566) − 275.21R 2 = 0.96363 (6) where y is the crack length, cm; t is the time of CDPTF in µ s. As is shown in Figure 14, the crack length generally increased more and more slowly, indicating that the crack propagation velocity decreased continuously with time.
The quantitative analysis of crack propagation speed was first carried out by Mott in metal materials in 1948. He believed that the crack propagation velocity was smaller than the sound velocity inside the object, and then gave an expression for the crack propagation velocity [30]: where V is the crack propagation velocity, m/s; V p is the P-wave velocity of the material, m/s; a 0 is the initial crack length, m; a is the instantaneous crack length, m. According to the Mott's crack propagation velocity equation, the maximum ratio of the peak crack propagation velocity to the P-wave velocity of the material is 0.38. The quantitative analysis of crack propagation speed was first carried out by Mott in metal materials in 1948. He believed that the crack propagation velocity was smaller than the sound velocity inside the object, and then gave an expression for the crack propagation velocity [30]: where V is the crack propagation velocity, m/s; Vp is the P-wave velocity of the material, m/s; a0 is the initial crack length, m; a is the instantaneous crack length, m. According to the Mott's crack propagation velocity equation, the maximum ratio of the peak crack propagation velocity to the Pwave velocity of the material is 0.38.
In 1957, Stroh proposed that when the propagation speed is close to the speed of sound, the peak propagation velocity should be the Rayleigh wave velocity, and the expression of the crack propagation velocity can be written as [31]: Where Vr is the wave velocity of the Rayleigh wave, m/s. In addition, there is the following relationship between the Rayleigh wave velocity and the Swave velocity [32]: Where Vs is the S-wave velocity, m/s; μ is the Poisson's ratio of the rock mass. At the same time, the following equation can reflect the relationship between P-wave velocity and S-wave velocity [33]: From Equations (8)-(10), we can draw a conclusion that the theoretical peak crack propagation velocity of feldspar quartz sandstone in the field test is Where Vmax is the peak crack propagation velocity, m/s. Since the crack propagation velocity of CDPTF is extremely fast and the duration is extremely short, it is assumed that the ratio of the crack growth length to the time between the two time points In 1957, Stroh proposed that when the propagation speed is close to the speed of sound, the peak propagation velocity should be the Rayleigh wave velocity, and the expression of the crack propagation velocity can be written as [31]: where V r is the wave velocity of the Rayleigh wave, m/s. In addition, there is the following relationship between the Rayleigh wave velocity and the S-wave velocity [32]: where V s is the S-wave velocity, m/s; µ is the Poisson's ratio of the rock mass. At the same time, the following equation can reflect the relationship between P-wave velocity and S-wave velocity [33]: From Equations (8)-(10), we can draw a conclusion that the theoretical peak crack propagation velocity of feldspar quartz sandstone in the field test is where V max is the peak crack propagation velocity, m/s. Since the crack propagation velocity of CDPTF is extremely fast and the duration is extremely short, it is assumed that the ratio of the crack growth length to the time between the two time points is the crack propagation velocity at the intermediate moment. The crack propagation velocity reaches the maximum value of 1894.25 m/s at 0.9 ms in the numerical simulation. From Table 1, the P-wave velocity of this feldspar quartz sandstone measured in the laboratory test was 4197 m/s. The peak crack propagation velocity was 45.13% of the P-wave velocity, which is slightly less than the theoretical peak velocity of 2316.5 m/s. Furthermore, the average crack propagation velocity was 639.68 m/s in the numerical test, which is 15.2% of the P-wave velocity. The law of crack propagation velocity obtained by the numerical simulation is consistent with the law obtained by other scholars in the laboratory tests [34,35], and the rationality of numerical simulation was also verified from the side.

Characteristics of Displacement
In order to investigate the displacement change rule in different directions during carbon dioxide phase transition fracturing, particle displacements of 0.5 m, 1 m, 1.5 m, 2 m, 2.5 m, 3 m and 3.5 m from the center of the borehole to the monitoring point in the primary and secondary impact directions were monitored. The stress-time curves of each particle were obtained and shown in Figure 15. There was a significant difference in displacement between the monitoring points in the primary and the secondary impact direction, but the displacement direction was the same as the impact direction. Although the peak displacements in the primary impact direction were greater than those in the secondary impact direction, they showed significant time delays in those two directions. In addition, they decreased with the increase of particle-borehole distance, and their decline rate gradually decreased. The peak displacement in the secondary impact direction showed different characteristics; it increased first and then decreased with the increase of the particle-borehole distance. Moreover, the peak displacements with the particle-borehole distance not exceeding 1.5 meters were smaller than that with the particle-borehole distance exceeding 2.0 meters. The reason is that the stress distribution around the monitoring points was affected, and it influenced the peak displacements in the secondary impact direction. The displacement-time curves of the primary impact direction had no obvious fluctuations, but the curves of the other direction had small fluctuations. It is because the rock mass was subjected to complex tensile stress and compression stress in the secondary impact direction, which influenced the displacement along this direction.
Energies 2020, 13, x FOR PEER REVIEW 12 of 16 theoretical peak velocity of 2316.5 m/s. Furthermore, the average crack propagation velocity was 639.68 m/s in the numerical test, which is 15.2% of the P-wave velocity. The law of crack propagation velocity obtained by the numerical simulation is consistent with the law obtained by other scholars in the laboratory tests [34,35], and the rationality of numerical simulation was also verified from the side.

Characteristics of Displacement
In order to investigate the displacement change rule in different directions during carbon dioxide phase transition fracturing, particle displacements of 0.5 m, 1 m, 1.5 m, 2 m, 2.5 m, 3 m and 3.5 m from the center of the borehole to the monitoring point in the primary and secondary impact directions were monitored. The stress-time curves of each particle were obtained and shown in Figure  15. There was a significant difference in displacement between the monitoring points in the primary and the secondary impact direction, but the displacement direction was the same as the impact direction. Although the peak displacements in the primary impact direction were greater than those in the secondary impact direction, they showed significant time delays in those two directions. In addition, they decreased with the increase of particle-borehole distance, and their decline rate gradually decreased. The peak displacement in the secondary impact direction showed different characteristics; it increased first and then decreased with the increase of the particle-borehole distance. Moreover, the peak displacements with the particle-borehole distance not exceeding 1.5 meters were smaller than that with the particle-borehole distance exceeding 2.0 meters. The reason is that the stress distribution around the monitoring points was affected, and it influenced the peak displacements in the secondary impact direction. The displacement-time curves of the primary impact direction had no obvious fluctuations, but the curves of the other direction had small fluctuations. It is because the rock mass was subjected to complex tensile stress and compression stress in the secondary impact direction, which influenced the displacement along this direction.

Characteristics of Vibration Velocity
Rock mass vibration caused by the rapid impact is usually expressed by the vibration velocity of particles. The vibration velocity of particles is a continuous reciprocating movement process along with the particles' equilibrium position caused by rapid impact. Vibration velocity is a reference between vibration intensity and seismic intensity and it is a safety criterion for damage to surrounding rock mass by the impact. Clarifying the vibration characteristics of CDPTF is of great significance for controlling the adverse effects. In order to explore the vibration velocity change in different directions, the same monitoring points like these in the previous section are selected to obtain the velocity-time curves in two directions, as shown in Figure 16. In general, the peak particle velocities (PPVs) in the primary impact direction are much larger than those in the secondary impact

Characteristics of Vibration Velocity
Rock mass vibration caused by the rapid impact is usually expressed by the vibration velocity of particles. The vibration velocity of particles is a continuous reciprocating movement process along with the particles' equilibrium position caused by rapid impact. Vibration velocity is a reference between vibration intensity and seismic intensity and it is a safety criterion for damage to surrounding rock mass by the impact. Clarifying the vibration characteristics of CDPTF is of great significance for controlling the adverse effects. In order to explore the vibration velocity change in different directions, the same monitoring points like these in the previous section are selected to obtain the velocity-time curves in two directions, as shown in Figure 16. In general, the peak particle velocities (PPVs) in the primary impact direction are much larger than those in the secondary impact direction. Since the PPV in the primary impact direction decreased at a rate of 26.16%, 20.85%, 16.81%, 15.98%, 11.84% and 10.44% each 0.5 m, the attenuation rate of PPV decreased in this direction with the increase of particle-borehole distance, indicating that the PPV attenuation trend gradually slowed down. In the secondary impact direction, the PPV decreased rapidly by 48.56% between 0.5 and 1 m, and the decline rate gradually decreased after a distance of 1 m, they were 20.29%, 6.9%, 4.9%, 2.9% and 1.6%, Energies 2020, 13, 1336 13 of 16 respectively. Although there were significant differences between the PPVs in the primary direction and those in the other direction, the PPVs showed obvious delays in two different directions. Moreover, the peak particle velocity direction was consistent with the impact direction in every monitoring point towards both directions.
direction. Since the PPV in the primary impact direction decreased at a rate of 26.16%, 20.85%, 16.81%, 15.98%, 11.84% and 10.44% each 0.5 m, the attenuation rate of PPV decreased in this direction with the increase of particle-borehole distance, indicating that the PPV attenuation trend gradually slowed down. In the secondary impact direction, the PPV decreased rapidly by 48.56% between 0.5 and 1 m, and the decline rate gradually decreased after a distance of 1 m, they were 20.29%, 6.9%, 4.9%, 2.9% and 1.6%, respectively. Although there were significant differences between the PPVs in the primary direction and those in the other direction, the PPVs showed obvious delays in two different directions. Moreover, the peak particle velocity direction was consistent with the impact direction in every monitoring point towards both directions.
(a) (b) Figure 16. Particle velocity curves of monitoring points at different directions:(a) particle velocitytime curves in the primary impact direction; (b) particle velocity-time curves in the secondary impact direction.

Characteristics of Effective Stress
Effective stress is a reliable mechanical criterion for judging whether the material fails. Determining the effective stress distribution characteristics in the rock mass during the rock breaking process can help to clarify the mechanism of CDPTF. In this monitoring, we select the same monitoring points like those in Section 5.2 to obtain the effective stress in the two impact directions, and the stress-time curves are shown in Figure 17. In general, the peak effective stresses in the primary impact direction are all greater than those at the corresponding positions in the secondary impact direction. The peak effective stress decreased when the unit-borehole distance increased in the primary impact direction, and its decline rate gradually decreased. Differently, as the unitborehole distance increased, the peak effective stress in the secondary impact direction first increased and then decreased, and the peak effective stress at 1.5 m was the largest. Although the change rule of effective stress was quite different in both primary and secondary directions, the peak effective stress did not reach the yield strength to damage the rock mass in the two directions during CDPTF. Furthermore, the peak effective stress showed a significant time delay in the two directions.

Characteristics of Effective Stress
Effective stress is a reliable mechanical criterion for judging whether the material fails. Determining the effective stress distribution characteristics in the rock mass during the rock breaking process can help to clarify the mechanism of CDPTF. In this monitoring, we select the same monitoring points like those in Section 5.2 to obtain the effective stress in the two impact directions, and the stress-time curves are shown in Figure 17. In general, the peak effective stresses in the primary impact direction are all greater than those at the corresponding positions in the secondary impact direction. The peak effective stress decreased when the unit-borehole distance increased in the primary impact direction, and its decline rate gradually decreased. Differently, as the unit-borehole distance increased, the peak effective stress in the secondary impact direction first increased and then decreased, and the peak effective stress at 1.5 m was the largest. Although the change rule of effective stress was quite different in both primary and secondary directions, the peak effective stress did not reach the yield strength to damage the rock mass in the two directions during CDPTF. Furthermore, the peak effective stress showed a significant time delay in the two directions.

Conclusion
CDPTF is a safe and reliable rock breaking technology. According to the energy accumulation characteristics of CDPTF, we established a mechanical model based on gathering energy load, carried out field tests and numerical simulations. The following conclusions are reached: (1) There is a clear gathering energy effect in the CDPTF. The impact reduction coefficient we

Conclusion
CDPTF is a safe and reliable rock breaking technology. According to the energy accumulation characteristics of CDPTF, we established a mechanical model based on gathering energy load, carried out field tests and numerical simulations. The following conclusions are reached: (1) There is a clear gathering energy effect in the CDPTF. The impact reduction coefficient we defined can effectively characterize the force difference between the primary impact direction and the secondary impact direction during the CDPTF. The mechanical model based on the impact reduction coefficient can be very well applied to the CDPTF numerical simulation in actual engineering.
(2) The cracks generated by CDPTF has obvious directionality, generally X-shaped in a rock mass of good quality. Crack propagation is fast and the entire crack propagation process takes a short time in CDPTF. The crack propagation velocity gradually decreases with time. The average value and the peak value of the longitudinal wave are 15.2% and 45.13%, respectively. Carbon dioxide phase transition fracturing shows intermittent expansion, and the length of cracks grows more and more slowly with time.
(3) There are obvious differences in the dynamic response characteristics in the primary and secondary directions. Although the displacements in the primary and secondary impact direction are both small, the displacements in the primary impact direction are greater than those of the other direction. The displacements of the secondary impact direction increase and then decrease with the increasing particle-borehole distance. In both directions, the peak point velocity decreases and decays more and more slowly as the particle-borehole distance increases. The PPV of the primary impact direction is greater than that of the other direction, and it decays slower than that of the secondary impact direction near the blast hole. For the peak effective stress, it decreases in the primary impact direction with the increase of the unit-borehole distance, and it increases first and then decreases in the other direction. Although the peak effective stress is quite different in both primary and secondary directions, the peak effective stress will not cause the rock damage in both two directions.