#### Abstract

Rock damage and vibration attenuation are the basis of blasting design, which will affect rock breaking results and the safety of structures (buildings). In this paper, the effect of blast-hole arrangement, delay time, and decoupling charge on the rock damage and the vibration attenuation in multihole blasting were numerically investigated. Through dynamic analysis software ANSYS/LS-DYNA, the JOHNSON_HOLMQUIST_CONCRETE (JHC) rock model and fluid-solid coupling method were used to establish single-hole and multihole rock blasting models. Based on the analysis of single-hole rock damage and vibration, firstly, the effect of the above three factors on the rock damage characteristics of multihole blasting was analyzed. Then, the extracted peak particle velocity (*PPV*) data of multihole blasting were fitted to the United States Bureau of Mines (USBM) equation to obtain the blasting vibration attenuation parameters (the site constant K and *α*). It is found that the blast-hole arrangement, delay time, and decoupling charge have a much smaller effect on *α* than K. At the same time, a delay-time-dependent model of PPV correction coefficient was proposed, and its rationality was verified by field data. The results obtained in this paper have reference significance for optimizing the blasting design and improving the blasting effect.

#### 1. Introduction

As an important construction method, engineering blasting is widely used in production and construction activities such as mining engineering, hydraulic engineering, tunnel engineering, and nuclear power engineering [1]. During the blasting process, a small part of the explosion energy is used to break the rock, causing the fragmentation and cracking of rock mass in the near zone, that is, the rock mass damage. The range of the near zone of the blasting source is affected by the lithology, explosive performance, and blasting methods, and it is estimated to be approximately 2–33/100 times the charge radius [2]. The near-zone damage characteristics caused by the explosion have an important influence on the subsequent grinding, crushing, or slagging efficiency and energy utilization. The remaining part of the explosion energy will be dissipated in the form of shock waves, thermal energy, and blasting seismic waves, accompanied by adverse effects such as vibrations, shock waves, fly rocks, toxic gas, and noise [3]. Among such effects, blast-induced vibration is considered the most adverse effect. The blast-induced vibration can reduce the stability and bearing capacity of the rock mass, resulting in damage or even destruction of buildings (structures). In the study of blast-induced vibration effects, usually, the peak particle velocity (*PPV*) is used as the vibration variable to evaluate the intensity of blast-induced vibration and control blasting safety [4–6]. The attenuation law of blasting vibration is one of the main bases for blasting design. Researchers from various countries have established some different empirical equations (7)–(11), as listed in Table 1.

Due to the anisotropy, nonlinearity, variability of mechanical properties of rock materials, and the coexistence of solid, liquid, and gas phases, the response of rock mass under explosive load is very complex [12]. Many scholars and engineers have conducted long-term research on rock blasting damage and blast-induced vibration. Based on the tensile expression method of the classic blasting damage model, Hu et al. [13, 14] constructed a new tensile-compressive damage model, which can quantitatively calculate the damage. The model was linked to LS-DYNA to analyze the excavation damage effects of smooth blasting and presplit blasting on high slopes. Numerical simulation results show that the blasting damage of smooth blasting has significant cumulative characteristics. Based on field tests and numerical simulations, Li et al. [15] studied the explosive load characteristics, the blast-induced damage characteristics for rock mass, and the PPV attenuation law and established the relationship between PPV at 30 m away from the charge hole and the damage depth of rock mass. A rock damage control method is proposed to restrict the rock damage depth under blasting excavation according to the safety threshold of vibration velocity for rock mass. Khandelwal and Singh [16] conducted multiple sets of field blasting tests and monitored the ground vibration values and then studied the artificial neural network method (ANN) to predict PPV. By comparing with the prediction results of PPV empirical formulas in various countries, it proved that ANN has a better prediction level.

The research on single-hole blasting law is representative, but the actual engineering blasting is mainly multihole blasting. Therefore, the study of the multihole blasting effect is of great significance to engineering practice. In multihole blasting, the proper hole arrangement will make the blasting energy evenly distributed, improve the crushing result of the rock, reduce the large lump rate, and increase the looseness. Based on field survey data, Zhou et al. [17] used numerical simulation and physical model test methods to simulate and optimize the four typical cut hole arrangements (9-hole, single-spiral-hole, double-spiral-hole, 2-empty-hole-diamond) of shaft formation by one deep-hole blasting, and the relevant blasting parameters obtained were successfully applied to field tests. Through numerical simulation of rock slope presplit blasting under in-situ stress, Wang et al. [18] compared the blasting results of rock mass with the square and plum-blossom-shape arrangement of blast-holes. The results show that the plum-blossom-shape arrangement can improve the blasting quality under prestress conditions. Delayed blasting has an important effect on rock damage and blast-induced vibration. With the gradual popularization of electronic detonators, the precise control of delayed blasting has attracted more attention from researchers [19–21]. Based on field data, Chen et al. [22] studied the influence of millisecond delay time, detonating sequence, and properties of the seed signal on the vibration reduction effect of millisecond blasting. Based on ANSYS/LS-DANY finite element software, Shao et al. [23] established a three-dimensional tunnel blasting model. The effects of different delay times on the dynamic response of tunnel lining are quantitatively studied, and combined with the *PPV* safety criterion, a reasonable millisecond delay time is obtained. Qiu et al. [24] theoretically analyzed the formation mechanism of blasting craters in the near-field and the mechanism of vibration reduction in the far-field under the condition of millisecond blasting. A comparison of field tests was carried out, and it was found that only when the delay time is less than the formation of a new free surface, it is possible to form a normal blasting crater. Compared with simultaneous blasting, delayed blasting can effectively reduce the *PPV* in the near-field. Using LS-DYNA, Hashemi and Katsabanis [25] simulated delayed blasting in the rock mass and explored the influence of initiation time on damage generation and morphological distribution. It is found that the delay time provides favorable conditions for the growth of cracks in the rock mass adjacent to the blast hole. Johansson and Ouchterlony [26] conducted small-scale tests on short-delay blasting. The results show that when the delay time is in the time range of the stress wave interactions, there is no distinct difference in the rock fragmentation compared with no shock wave interactions. Tatsuya et al. [27] simulated the relationship between *PPV* at certain monitoring points and delay time through computer programs and selected the best delay time. A design concept based on combined delayed blasting is proposed, which can effectively reduce the *PPV* of important monitoring points. The decoupling charge is a charge structure often used in controlled blasting. Blasting with decoupling charge can greatly reduce underexcavation, overexcavation, and side pumice, improve construction efficiency, and have an important effect on rock damage and blast-induced vibration. Through numerical simulation, Shao et al. [28] studied the law of rock crack propagation in smooth blasting. The calculation formula of charge decoupling coefficient in smooth blasting is proposed, which can reduce the damage and improve the stability of surrounding rock in smooth blasting. Qi [29] used software ANSYS to simulate rock blasting with decoupling charges and analyzed the corresponding stress response values and nephogram under different radial decoupling coefficients. It is found that the rock mass blasting stress response is the largest in the case of a coupled charge and gradually decreases with the increase of the decoupling coefficient, indicating that the air medium can buffer the effect of detonation waves. Wang et al. [30] also reached the same conclusions. Through theoretical analysis and blasting model tests of decoupling charge blasting, Lou et al. [31] studied the formation and propagation of shock waves, as well as the impact pressure characteristics of blast-hole walls. Results show that when the decoupling coefficient is 1.5–3.5, the impact pressure is evenly distributed, the impact time is longer, and the blast-induced vibration is weaker. With the increase of the decoupling coefficient, the initial impact pressure shows a trend of first increasing and then decreasing.

Existing studies have used the methods of theoretical analysis, field test, numerical simulation, and model test to study rock blasting damage and blast-induced vibration. Among them, numerical simulation is a common method for studying blasting, which can easily control blasting conditions such as the number, size, arrangement, and initiation time of blast-holes, and has the advantages that other methods cannot replace. Among the numerical methods, the large-scale general explicit dynamic finite element analysis software LS-DYNA is widely used in the dynamic response of rock engineering under impact loads. LS-DYNA has two methods for blasting analysis: the Lagrange method and the multimaterial fluid-solid coupling method. In the Lagrange method, the grid and the analysis structure are integrated, and the finite element nodes are material points, and clear material interfaces can be obtained. However, serious distortions will occur during the blasting, leading to calculation errors. In the multimaterial fluid-solid coupling method, explosives and other fluid materials (air, water, etc.) use Euler or ALE (Arbitrary Lagrange-Euler) algorithm, and solid materials use the Lagrange algorithm, and the two interact through fluid-solid coupling. In numerical modeling, grid division is an important step. The quantity, uniformity, and division quality of model grids have a great influence on calculation accuracy, scale, and stability. Under the condition that the computer hardware performance can support the calculation, it is necessary to ensure the regular shape of the units, and the parts of concern and with large deformation should be appropriately refined.

The outline of the work given in this paper is as follows. First, we establish single-hole and multihole rock blasting models by software ANSYS/LS-DYNA in Section 2. Then, we analyze the single-hole blasting damage and vibration effect in Section 3, followed by the analysis of multihole rock damage in Section 4 and the analysis of multihole vibration attenuation characteristics in Section 5. Finally, we conclude our paper in Section 6.

#### 2. Models and Materials

The ANSYS/LS-DYNA numerical simulation software was used to simulate the cylindrical charge blasting. To avoid calculation errors caused by the distortion of the model grid, the multimaterial fluid-solid coupling method was adopted. That is, the explosive and air materials were assigned to the ALE algorithm, and the rock material was assigned to the Lagrange algorithm, and then the two algorithms were coupled.

##### 2.1. Model Design

According to the different analysis focus, the limitation of the model grids number and the high-quality grid’s requirements, when studying the rock damage in the near zone and the blasting vibration effect in middle and far zone, the analysis models of the near zone and middle and far zone are established, respectively. The size of the near-zone damage analysis model is 2.5 m × 2 m × 6 m (*XY* × *Z* direction). To better observe the damaged area, dense meshing is carried out, and the element size is set small. The size of the middle and far zone analysis model is 30 m × 4 m × 6 m (*X* × *Y* × *Z* direction), and the grid is sparser than that in the near area. At the same time, 5 vibration measurement points are set at different horizontal distances (*R*) from the blast-hole: *R* = 5 m, 7.6 m, 11.6 m, 17.7 m, and 27 m. In all numerical models, the charge radius is 0.035 m, the blast-hole radius is 0.045 m (adjusted according to the decoupling coefficients in the analysis of decoupling charge blasting). The depth of the blast-hole is 5 m, of which the stemming length is 2 m, and the charge height is 3 m. The detonation points are all set as the geometric center of the explosive. Nonreflecting boundary conditions are applied to the bottom and surroundings of the models to avoid the influence of reflected waves. The upper surface of the model is free. Section *a* is perpendicular to the *Z*-axis, and section *b* is perpendicular to the *Y*-axis, both passing through the center point of the explosive. The single-hole analysis model considering near zone damage is shown in Figure 1, and the single-hole analysis model considering middle and far zone vibration is shown in Figure 2.

**(a)**

**(b)**

**(c)**

**(a)**

**(b)**

**(c)**

When building the multihole blasting model, based on the single-hole model, only the number of blast-holes, the hole arrangements, the delay times, and the decoupling coefficients were changed. Numerical simulations were then carried out, and the results of single-hole blasting rock damage and vibration were taken as comparison standards to analyze the influence of hole arrangement, delay time, and decoupling charges on the blasting effect.

##### 2.2. Material Models

###### 2.2.1. Explosive

The explosive material is defined by the MAT_HIGH_EXPLOSIVE_BURN keyword and the Jones-Wilkins-Lee (JWL) state equation are used to define the detonation state, and the keyword is JWL_EOS. This state equation expresses the relationship between the detonation pressure *P*_{e} and the relative volume *V*_{e} of the detonation product and the internal energy per unit volume *E*_{e}:where *A*_{e}, , , , are constants related to the properties of explosives, *V*_{e0} is the initial relative volume, and *E*_{e0} is the initial internal energy per unit volume. The related parameters of the explosive are listed in Table 2, where VoD represents the detonation velocity, and PCJ represents the Chapman–Jouget pressure.

###### 2.2.2. Air

The air material adopts the MAT_NULL model, which can better simulate the fluid properties. The corresponding state equation is defined by the JWL_LINEAR_POLYNOMIAL keyword, as shown below:where *P*_{a} represents the detonation pressure, *V*_{a} represents the relative volume, *V*_{a0} represents the initial relative volume, *E*_{a} represents the internal energy per unit volume, *E*_{a0} represents the initial internal energy per unit volume, *C*_{0}, *C*_{1}, *C*_{2}, *C*_{3}, *C*_{4}, *C*_{5} and *C*_{6} are coefficients, and *μ* represents the relative density. The related parameters of the air are listed in Table 3.

###### 2.2.3. Rock

The rock material adopts the MAT_JOHNSON_HOLMQUIST_CONCRETE (JHC) constitutive model, which can describe large strains, high strain rates, and high pressures. JHC constitutive model has many parameters. Among them, the density *ρ* and uniaxial compressive strength *f*_{c} are based on laboratory test results. Shear modulus *G* is a combination of elastic modulus *E* and Poisson’s ratio *ν* and is calculated according to the following formula:

The rock Brazilian tensile strength (*T*) is calculated from *f*_{c} according to the following formula [33]:

The limit surface parameters (*A*, *B*, and *N*) have a great influence on the blast-induced vibration, which is determined by the basic rock parameters according to the following method [34, 35].

The JHC normalized equivalent stress is expressed as a function of pressure, strain rate, and damage:where *A* is the normalized cohesive strength, *B* is the normalized pressure hardening coefficient, *N* is the pressure hardening exponent, *C* is the strain rate coefficient, *D* is the damage parameter, is the normalized hydrostatic pressure, and is the dimensionless strain rate. , and are defined as follows:where *σ* is the realistic equivalent stress, *P* is the realistic hydrostatic pressure, is the realistic strain rate, and is the reference strain rate.

Without considering the effects of damage and strain rate, (5) is abbreviated as follows:

Based on plastic theory, both JHC and Mohr-Coulomb models pass through the following two points on the compression meridian plane: uniaxial compression and pure shear . The following relationship can be obtained:

Based on the Mohr–Coulomb model, the maximum principal stress has the following relationship with the confining pressure :

According to the Mohr–Coulomb model, the maximum principal stress has the following relationship with the confining pressure , the cohesion *c*, and the internal friction angle *φ*:*σ* and *P* are expressed as follows:

Substituting the laboratory data of triaxial compression into (7), the other two limit surface parameters of the JHC model (*B* = 2.298 and *N* = 1.0344) were obtained.

##### 2.3. Definition of Rock Damage

The damage variable *D* of the JHC model (0 ≤ *D* ≤ 1, *D* = 0 means intact rock, *D* = 1 means completely broken rock) is accumulated from the equivalent plastic strain and the plastic volume strain , expressed as follows:where is the normalized maximum Brazilian tensile strength (*T* is the realistic maximum Brazilian tensile strength), *D*_{1} and *D*_{2} are material constants.

The damage strength *DS* is defined as follows:where *SF*max is the normalized maximum strength.

In addition, the relevant parameters of the JHC model also include the plastic strain before fracture (*EF*min), the crushing pressure (*p*_{c}), the crushing volume strain (*u*_{c}), the compaction pressure *p*_{l}, the compaction volume strain *u*_{l}, and three material constants *K*_{1}, *K*_{2} and *K*_{3}. These parameters are less sensitive to the calculation results, and the initial parameter values of the JHC model can be used. The JHC model is used to simulate the slightly weathered granite and the related parameters of the JHC model are listed in Table 4.

#### 3. Simulation of the Single-Hole Blasting

The model considering near zone damage (Figure 1) and the model considering middle and far zone vibration (Figure 2) of the single hole were calculated, and the rock damage and vibration attenuation results were obtained.

##### 3.1. Rock Damage

The damage nephogram in the near zone of the single-hole blasting is shown in Figure 3, and different colors represent different damage levels, that is, different *D* values. It can be seen from the damage nephogram of section *a* that the damaged area is centered on the charge and extends outward along the radial direction of the blast hole. The damage level changes from heavier to lighter and finally dissipates in a scattering shape. From the damage nephogram of section *b*, it can be seen that the damaged area extends outward from the blast hole. The damage level of the initiation point is smaller than the upper and lower parts of the explosive, and the boundary of the damage zone is smoother than section *a*.

**(a)**

**(b)**

According to the research of Hu et al. [13] on the safety threshold of blasting damage, the critical damage value of retained rock mass is *D* = 0.19, so the area of *D* = 0.19–1 is taken as the damaged area. Calculated from the damage contour of the single-hole blasting model, the damaged area of section *a* is about 35 times of the blast-hole cross section, and the damaged area of section *b* is about 8 times of the blast-hole longitudinal section.

##### 3.2. Blast-Induced Vibration

The *PPV*s of the five measurement points of single-hole blasting are shown in Figure 4(a), the horizontal radial direction (*X*-direction) is the largest, the vertical direction (*Z*-direction) is second, and the horizontal tangential direction (*Y-*direction) is the smallest. The horizontal radial velocity-time curve at *R* = 5 m is shown in Figure 4(b). It can be seen that the vibration duration is about 50 ms, the part with large vibration amplitudes is basically within 20 ms.

**(a)**

**(b)**

Based on the USBM equation (Table 1), the single-hole numerical calculation data is subjected to regression analysis, and the horizontal radial attenuation parameters of the single-hole blasting are obtained. The fitting result is shown in Figure 5, and the attenuation equation is as follows:

#### 4. Multihole Blasting Rock Damage

Based on the single-hole analysis model considering near zone damage (Figure 1), only the number of blast-holes, the blast-hole arrangements, the delay times, and the decoupling coefficient are changed to establish the multihole blasting models and conduct numerical simulations.

##### 4.1. Hole Arrangement

The hole arrangements often used in multihole blasting can be summarized into two types: rectangle and triangle. The rectangle and triangle arrangements represented by 4 blast-holes are shown in Figure 6.

**(a)**

**(b)**

Figure 7 shows the near-zone damage nephograms of blasting models with different hole arrangements under the same explosive quantity and the same initiation time. It can be seen that, compared with rectangular holes, the size of rock mass fragments under the triangular hole arrangement is more uniform and smaller. The rectangular hole arrangement has a greater possibility of producing larger-sized fragments, which is not conducive to rock mass breaking and engineering production.

**(a)**

**(b)**

Calculating the damaged area in the range of *D* = 0.19–1 under the rectangular hole arrangement is 1.543 m^{2}, accounting for 34.3% of the total cross-sectional area. The damaged area of the triangular hole arrangement is 1.581 m^{2}, accounting for 35.1% of the total cross-sectional area. The damaged area of the triangular hole arrangement is 1.024 times that of the rectangle. That is, the breaking capacity of the triangle on the rock mass is slightly larger than that of the rectangle.

##### 4.2. Delay Time

Millisecond blasting is a controlled blasting technique often used in precise blastings such as rock foundation excavation, mineral mining, and demolition of structures. In millisecond blasting, multiple blast-holes in a group of charges are detonated in sequence according to a set time interval. It is an effective method to improve rock fragmentation and reduce blast-induced vibration. The analysis is simplified to a group of two blast-holes. As shown in Figure 8, explosive No. 1 detonates first, and explosive No. 2 detonates after an interval of Δ*t*.

Considering that the vibration duration is about 50 ms in single-hole blasting, the rock damage conditions when the delay time Δ*t* is 0 ms, 1 ms, 5 ms, 10 ms, 15 ms, 20 ms, and 50 ms are simulated, respectively. Among them, Δ*t* = 0 ms means that No. 1 and No. 2 explosives detonate at the same time. Figure 9 shows the damage nephogram when Δ*t* is 0 ms, 5 ms, and 50 ms.

**(a)**

**(b)**

**(c)**

Define the damage area ratio *λ* as follows:where *j* = 1, 2 represent No. 1 and No. 2 explosive; *i* = 1, 5, 10, 15, 20, 50 represent the delay time; represents the damaged area of No. *J* explosive when Δ*t* = *i*, and represents the damaged area of explosive No. *J* when detonated simultaneously; is the damage area ratio of explosive No. *j*.

Figure 10 shows the damage area ratio of No. 1 and No. 2 explosives with different delay times. Compared with simultaneous detonation, in the case of delayed detonation: (1) The damaged area of explosive No. 1 which detonated first was slightly affected by Δ*t*; (2) Delayed detonation provides time for crack propagation, that is, the damaged area of explosive No. 2 after detonation increases with the increase of Δ*t*. When Δ*t* = 50 ms, it is close to the duration of the vibration waveform. At this time, *λ*_{2} = 1.406, which is an increase of about 40% compared with Δ*t* = 0 ms. It can be considered that *λ* can increase by about 40% at most compared with simultaneous detonation. According to Stagg [36], the possible reason is as follows: when detonating at the same time, what is produced between explosives No. 1 and No. 2 is the superposition of stress waves, which breaks the rock. During the delay blasting, the No. 1 explosive that detonated first produces detonation gas, which did not dissipate within the short-delay time. The stress wave generated after the detonation of the No. 2 explosive interacts with the detonation gas produced by the No. 1 explosive, resulting in greater crushing power. That is, delayed blasting promotes the interaction between the stress wave induced by the latter blast-hole and the detonation gas generated by the previous blast-hole, and the interaction between the stress wave and the explosive gas can promote rock breaking result more than the superposition of the stress wave.

##### 4.3. Decoupling Charge

The charge structure is one of the factors that have a decisive influence on the blasting effect. In controlled blasting, decoupling charges are usually used. The decoupling coefficient is an important control parameter for the decoupling charge blasting, which is divided into the radial decoupling coefficient and the axial decoupling coefficient. The radial decoupling coefficient refers to the ratio of the blast-hole diameter to the charge diameter, and the axial decoupling coefficient refers to the ratio of the blast-hole length to the charge length in the direction of the blast-hole length. The construction operability of radial decoupling charges is better than that of axial decoupling charges. Therefore, radial decoupling charges are often used in engineering practice. In the radial decoupling charge blasting, the air layer between the explosive and the blast-hole wall can reduce the shock wave pressure generated by the explosion, thereby reducing excessive damage to the rock mass and achieving better crushing results.

The radial decoupling coefficient *K*_{r} is as follows:where *d*_{b} and *d*_{e} are the diameter of the blast-hole and the diameter of the charge, respectively.

The model is simplified as two blast-holes, and the decoupling coefficient *K*_{r} is set to 1, 1.5, 2, 2.5, 3, and 5, respectively. *K*_{r} = 1 is the decoupling charge under the limit condition, which is equivalent to the coupling charge.

The damage nephogram of blasting models with *K*_{r} = 1, 2, and 5 is shown in Figure 11. According to the rock damage results, with the increase of *K*_{r} value, the damaged area around and the connecting part of the blast holes gradually decreases. When *K*_{r} = 1 (coupling charge), the explosive detonation wave acts on the rock directly, and the energy transfer rate is the highest. When *K*_{r} > 1 and increases, it means that the air layer between the explosive and the rock is getting thicker, and the energy transfer rate, blast-hole pressure, and rock damage area are all reduced. Smooth blasting and presplit blasting use this principle to achieve the expected rock-breaking result.

**(a)**

**(b)**

**(c)**

#### 5. Multihole Blast-Induced Vibration Attenuation

Based on the single-hole analysis model considering middle and far zone vibration (Figure 2), only the number of blast-holes, the blast-hole arrangements, the delay times, and the decoupling coefficient are changed to establish the multihole blasting models considering middle and far zone vibration and conduct numerical simulations.

##### 5.1. Hole Arrangement

The *PPV*s of the measurement points in the middle and far zone under different hole arrangements are shown in Figure 12(a). There is almost no difference in vertical vibration *PPV*s between rectangular and triangular holes. The main differences lie in horizontal radial and horizontal tangential vibrations: (1) In rectangular hole blasting, the horizontal radial *PPV* is larger than that of the triangular hole, and the difference is about 30 cm/s. (2) In triangular holes, the horizontal tangential *PPV* is larger than that of the rectangular hole, indicating that the triangular hole arrangement is more likely to stimulate the generation of shear waves. Therefore, in the blasting practice, the blast-hole arrangement method should be selected according to the situation of protected items. For structures or instruments that are sensitive to radial vibrations, a triangular blast-hole should be selected, and for tangential vibrations, a rectangular blast-hole can be considered. (3) At a far distance(≥30 m), the effect of the blast-hole arrangement on vibration tends to be weakened.

**(a)**

**(b)**

According to the USBM equation, the regression curves and optimum equations of horizontal radial blasting vibrations under different hole arrangements are shown in Figure 12(b). It can be seen that the site constants *K* and *α* of the rectangular-hole blasting are 60.062 and 1.7886, which are larger than those of the triangular-hole blasting (*K* = 45.699, *α* = 1.5883).

##### 5.2. Delay Time

Calculate the vibration PPV when the delay time Δ*t* is 0 ms, 1 ms, 5 ms, 10 ms, 15 ms, 20 ms, and 50 ms. Figure 13 takes the *X*-direction blasting vibration waveforms at *R* = 5 m of Δ*t* = 0 ms, 15 ms, and 50 ms as examples. The results show that with the increase of Δ*t*, two peaks gradually appear. When Δ*t* = 50 ms, the two wavebands are independent because 50 ms is close to the duration of a single-hole waveform.

**(a)**

**(b)**

**(c)**

Figure 14 shows the variation of horizontal radial PPVs with different delay times. It can be seen that the delayed blasts produce distinctly smaller PPVs than does the simultaneous blasting: (1) When Δ*t* = 0 ms, that is, simultaneous blasting, the PPV value is the largest, and the vibration superimposition degree is the highest. The reason is that millisecond blasting can reduce the charge quantity per delay, which allows the energy of the explosion to be released in stages. In the case of simultaneous blasting, the charge quantity is the total amount of 2 holes, while in the delayed blasting, the charge quantity per delay is the amount of 1 hole. The blast-induced vibration increases with the increase of the charge quantity. (2) When Δ*t* ≠ 0, with the increase of Δ*t*, the degree of vibration superposition decreases, and the *PPV* value decreases. When Δ*t* = 1 ms, PPV is close to Δ*t* = 0, and the vibration superimposition degree is greater. When Δ*t* = 15 ms and 20 ms, the *PPV* is close to the single-hole *PPV*, and the vibration superposition is less. When Δ*t* = 50 ms, there is no difference between the PPV value and single-hole PPV value. At this time, the interval time is the same as the waveform duration, and there is almost no vibration superimposition. (3) Compared with simultaneous blasting, millisecond blasting can effectively reduce PPV values, especially in nearby zones: At measurement points, *R* = 5 m, 7.6 m, 11.6 m, 17.7 m, and 27 m, the distribution ranges of PPVs under different delay times are 42–66 cm/s, 21–29 cm/s, 11–15 cm, 6–8 cm/s, and 3–4 cm/s, the difference is small in the far distance and large in the near. It indicates that the delay time has a small influence on the PPVs of long-distance measurement points but has a large influence on short-distance measurement points. Setting a large delay time is better for vibration control, but it will increase the overall vibration duration, which is unreasonable for engineering practice [23]. Therefore, the delay time needs to be selected reasonably in engineering practice, and it should not be too large or too small.

The horizontal radial PPVs of two holes blasting with different delay intervals and the regression curves and optimum equations based on the USBM equation are plotted in Figure 15. When Δ*t* = 1 ms, 5 ms, 10 ms, 15 ms, 20 ms, 50 ms, the coefficients of determination are 0.99897, 0.99899, 0.9994, 0.99987, 0.99991, 0.99915. When Δ*t* = 0, it means simultaneous blasting, with the largest charge quantity *W* and the smallest *K* value. When Δ*t* ≠ 0, the charge quantity *W* is half of that when Δ*t* = 0. As Δ*t* increases, the *K* value tends to decrease. The constant *α* is in a stable range without large fluctuations.

**(a)**

**(b)**

###### 5.2.1. Correction Coefficient

To quantify the relationship between the change of blast-induced vibration PPVs and the delay times, define *ß* as the correction coefficient of the radial vibration *PPV* when the delay time is Δ*t*:where PPV_{single} is the radial vibration PPV of the single-hole blasting.

The empirical equations of blast-induced vibration with charge and propagation distance as independent variables have been well-known in the field of blasting engineering and have been widely used. Take the empirical USBM equation as an example, when Δ*t* is constant:where *K* and *α* are the site coefficients of single-hole blasting, *K*_{Δt} and *α*_{Δt} are the site coefficients of millisecond blasting when the delay time is Δ*t*. The single-hole blasting charge quantity is the same as the maximum charge quantity of delayed blasting, both are *W*. Substitute into equation (20):

It can be seen that *ß* and *R* have an inverse proportional function relationship. Let , , then equation (22) can be written as follows:

If the duration of the single-hole blasting waveform is *T*, let Δ*t* = *nT*; in this case, *T* = 50 ms, then *n* = Δ*t*/50. When the distance *R* is constant, it can be seen from Figure 16 that as *n* increases, *ß* is approximately inversely proportional to *n*. Assuming that when Δ*t* > 0, *ß* is in inverse proportion to *R* and *n*:where *b*_{1}, *b*_{2}, *b*_{3} are constants.

Take the logarithm of both ends of the equation (20) to get the following:

MATLAB was used to fit equation (21) based on the least-square method. *b*_{1} = 1.0666, *b*_{2} = 0.04046, and *b*_{3} = 0.0855 were obtained.

Obtaining equation (20) is as follows:

The regression coefficient of determination is 0.8866, and the sum of squared error is 0.04431.

###### 5.2.2. Equation Verification

According to the field tests conducted by Chen et al. [22] when studying the effect of millisecond blasting to reduce vibration, the signal collected from one detonation was taken as the single-hole vibration waveform, and the signal collected from two detonations with a certain delay time was taken as the double-hole vibration waveform. The nonel detonator was used for detonation, and millisecond times of 50 ms and 100 ms were selected, which are the delay times of the two sets of tests. Two measurement points were arranged with a distance of 18 m and 30 m. The single-hole charge quantity and the double-hole maximum charge quantity are both 10 kg. In the above numerical simulation of blasting, the single-hole charge quantity and the double-hole maximum charge quantity are both 11.6 kg. The value of the two charge quantities is close. The *PPV* results of the field tests are listed in Table 5. The waveform duration *T* of the measured time history curve in single-hole blasting is 200 ms.

According to the data in Table 5, the rationality of equation (22) is verified. The distance *R*, the delay time Δ*t*, and the waveform duration *T* of each measurement point were substituted into the equation (22) to obtain the calculated value of *ß*. The measured single-hole PPV values and the calculated values of *ß* were substituted into equation (16), and the predicted values of the double-hole PPV were obtained. The relevant calculation results are listed in Table 6. The relative errors between the predicted values and the measured values are not very small, which is inevitable. Each site has its unique lithology characteristics, and the types and weights of explosives, charging methods, delay accuracy are also different. The correction coefficient *ß* is affected by these factors, but it can still meet the needs of the project.

##### 5.3. Decoupling Charge

Change of the horizontal radial *PPV*s in the middle and far zone of blasting models with *K*_{r} = 1, 1.5, 2, 2.5, 3, and 5 is shown in Figure 17(a): (1) With the increase of *K*_{r} value, the PPVs decrease greatly, indicating that the decoupling charge can reduce the vibration effect of the explosion. (2) As the propagation distance increases, the influence of *K*_{r} on PPV gradually decreases, indicating that the vibration reduction effect of the decoupling charge will be affected by the distance. That is, the vibration reduction effect is large near and small at a distance.

**(a)**

**(b)**

The regression curves and optimum equations based on the USBM equation of horizontal radial blasting vibrations under different decoupling coefficients are shown in Figure 17(b). The site constants *K* gradually decreases with the increase of *K*_{r}, and the rate of decrease gradually slows down. The value of *α* is still stable.

##### 5.4. Comparison of Vibration Attenuation Law with Different Blast-Holes

Figure 18 is a comparison of horizontal radial *PPV*s at the 5 m measurement point with one hole, two holes (two holes_{1} refers to delay blasting, two holes_{2} refers to decoupling blasting), and four holes. The general rule is the blasting cases with more blast-holes produced larger PPVs than did the cases with fewer blast-holes. When the decoupling coefficient of two-hole blasting is large, the PPV will drop to less than that of single-hole blasting.

Based on the site constants *K*_{single} and *α*_{single} of single-hole blasting, the changes of *K*/*K*_{single} and *α*/*α*_{single} with the number of blast-holes under delayed blasting, decoupling charge blasting, and different hole arrangements are shown in Figure 19. The variation range of *K* is within 0.450–1.507. When *W*^{1/2}/*R* are the same (the single-hole blasting and the delayed blasting when Δ*t* ≠ 0), the site constant *K* increases as the number of blast-holes increases. When simultaneous blasting (the single-hole blasting, the delayed blasting when Δ*t* = 0, decoupling blasting, rectangular and triangular hole blasting), as the number of blast-holes increases, *W*^{1/2}/*R* increases, and the site constant *K* generally shows a decreasing trend. The reason is that as the number of blast-holes increases, the increase in charge quantity is greater than the increase in *PPV*. The site constants *K* for the cases of 20 ms and 50 ms are almost the same as those of the single-hole case since the vibration waveforms from the 20 ms and 50 ms blasting can be separated within the interval time. The variation range of *α* is within 0.932–1.118, and the number of blast-holes has a much smaller influence on *α* than *K*.

**(a)**

**(b)**

#### 6. Conclusions

Based on the ANSYS/LS-DYNA software, single-hole and multihole blasting simulation models have been established. For the damage in the near zone and the vibration effect in the middle and far zone of multihole blasting, the influence of the hole arrangement method, the delay time, and the decoupling charge are studied. Based on the USBM equation, the *PPV* data are fitted to obtain the blasting vibration attenuation parameters. A calculation equation of *PPV* correction coefficient with delay time as the independent variable is proposed, and the rationality of the equation is verified by field data. The blasting vibration parameters and *PPV*s of different numbers of blast-holes are compared. The following conclusions are obtained.(1)In single-hole blasting, in the cross section, the rock damaged area is about 35 times that of the blast-hole, and in the longitudinal section, it is about 8 times. The blasting vibration parameters of the single-hole blasting are *K* = 77.91, *α* = 1.6.(2)In multihole blasting, compared with rectangular holes, the size of rock mass fragments under the triangular hole arrangement is more uniform and smaller, which is beneficial to engineering production. The horizontal radial *PPV*s of rectangular holes are larger than that of triangular holes, but the horizontal tangential *PPV*s are smaller than that of triangular holes. Therefore, in blasting practice, the blast-hole arrangement should be selected according to the vibration-sensitive direction of the protected items. The site constants *K* and *α* of the rectangular-hole blasting are larger than those of the triangular-hole blasting.(3)In the two-hole millisecond blasting, the damaged area of the first blast-hole is not affected by the delay time, and that of the second blast-hole increases with the increase of delay time, but the speed of increase gradually slows down, increasing by about 40% at most. Compared with simultaneous blasting, delayed blasting reduces the charge quantity, thereby effectively reducing *PPV*, especially in nearby zones. The *PPV* correction coefficient is inversely proportional to the delay time. When Δ*t* ≠ 0, the *K* value tends to decrease as Δ*t* increases, and the constant *α* is in a stable range.(4)With the increase of the decoupling coefficient *K*_{r}, the rock damage area gradually decreases, and the horizontal radial and vertical *PPV*s decrease greatly, indicating that the decoupling charge can effectively reduce the vibration of the explosion. The site constants *K* gradually decreases with the increase of *K*_{r}, and the value of *α* is still stable.(5)The general rule is the blasting cases with more blast-holes produced larger *PPV*s than did the cases with fewer blast-holes. When *W*^{1/2}/*R* are the same, the site constant *K* increases as the number of blast-holes increases. When simultaneous blasting, the site constant *K* generally shows a decreasing trend as the number of blast-holes increases. The number of blast-holes has a much smaller influence on *α* than *K.*

It should be mentioned that the simulation in this paper is based on blasting with idealized conditions, and there is a certain gap between the complex blasting conditions in engineering practice. Thus, readers should be cautious of any extrapolation from the results above. Further investigations are needed to study the rock damage and blasting vibration under more complex blasting conditions through theoretical analysis, numerical methods, and field tests in the future.

#### Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

#### Conflicts of Interest

The authors declare that they have no conflicts of interest.

#### Acknowledgments

The research was financially supported by the National Key R&D Program of China (Grant no. 2020YFA0711802) and the National Nature Science Foundation of China (Grant nos. 51439008 and 51779248). The authors gratefully acknowledge all the support of this work.