Optimization of the Air Deck Blasting Parameters on the Basis of the Holmquist-Johnson-Cook Constitutive Model

The present study considers the so-called air deck blasting, one of the most commonly used techniques for the improvement of blasting efficiency in mining applications. In particular, it aims to improve the operating conditions of large-scale equipment, increase the efficiency of the slope enlarging process, and reduce the mining cost. These objectives are implemented through a two-fold approach where, first, a program for slope enlarging based on the middle air-deck charge blasting-loosening technology is proposed, and second, the physical mechanism responsible for the stress wave attenuation in the rock is analyzed in the framework of a Holmquist–Johnson– Cook constitutive model. Field test results verify that the proposed approach is highly efficient and economical when used in medium-hard rock blasting and thus provides a basis for the practical implementation of slope enlargement for high benches in Open-Pit Mines.


Introduction
Given that rocks are non-uniform materials, their fracture process depends on the strain rate, and their air-deck charge has a complicated mechanism and flexible setting [1][2][3][4][5][6]. Therefore, studies on the change process of the stress field around the blast hole that used the dynamic crack propagation model have theoretical significance and an engineering value. The numerical simulation method has been utilized to research the factors that insluence dynamic crack propagation model. These factors include the loading rate of the explosive load, distance from the blast hole to the free boundary, and size of the air deck in the blast hole [7][8][9][10][11][12][13][14][15]. The stress wave or initial velocity is adopted as an input, and analyzed step by step in accordance with the time step. The influence of mass and acceleration on the mechanical balance is considered in each step [16]. Moreover, the elastic dynamic finite element program has been used for stress analysis. The maximum tensile stress criterion and the Mohr-Coulomb criterion have also been employed to determine whether an element is damaged and obtain the dynamic crack propagation process of the rock material under the explosive load of the air-deck charge [17].
Typically, metal mines have a bench height of 10-12 m. There is a problem in Open-Pit Mine is that when second-phase slope enlargment uses the conventional step-by-step, near-slope method, the temporary slope becomes steep, the slope enlargment site becomes narrow, and the operation safety of large-scale equipment is difficult to guarantee. Therefore, a technical solution for the slope enlargment of a 24 m high bench is introduced in this study. The proposed solution addresses the technical problem of high-efficiency and low-cost blasting of a 24 m high bench in the case of a deep vertical hole with a large-diameter [18][19][20][21][22].
A large air-deck charge structure and a large density coefficient are usually utilized to arrange blast holes. However, the blast hole has no charge in the middle, making the explosion energy insufficient and the production of large rocks easy [23]. The Open-Pit has undergone the 24 m high-level blasting tests in its early stages. Although several useful results have been achieved, promoting the use of such tests in mines remains problematic. Specifically, the slope enlargment of the 24 m-high bench still needs to address three technical problems [24]. First, the perforation depth of a domestic ordinary roller drill cannot work for the large diameter drilling of the 24 m-high benches directly. Second, the slope of the bench is limited. With the increase in step height, the resistance line of the chassis on the front row of blast holes is doubled. Lastly, the unit consumption of explosives become much higher than ordinary [25][26][27][28].
The upper air-deck charge has been partially studied in the early stage, and it has been promoted and applied in some working conditions, and good results have been achieved [29]. This study uses LS-DYNA to investigate the stress distribution generated by columnar charge with an air-deck charge structure on the middle. A fluid-structure coupling simulation is conducted on a 24 m-high bench blasting under this condition to analyze the stress of the blasted rock mass and provide a reference for on-site high-bench blasting tests.

Holmquist-Johnson-Cook Theoretical Model
Holmquist, Johnson, and Cook proposed the Holmquist-Johnson-Cook (HJC) damage constitutive model to describe the nonlinear deformation and fracture characteristics of engineering and geological materials under large strain and high variability. This model considers high hydrostatic pressure, high strain rate, and the material damage effect [30][31][32][33][34][35].

Strength Model
The HJC strength model is described by the characteristic equivalent stress, as shown in Fig. 1. It is expressed in this study as follows: where r Ã ¼ r f c is the characteristic equivalent stress (r Ã S max , in which S max is the maximum characteristic equivalent stress that the material can achieve). r is the actual equivalent stress, f c is the quasi-static uniaxial compressive strength, and D is the damage variable.
is the characterized strain rate, _ e is the actual strain rate and _ e 0 ¼ 1:0 s À1 is the reference strain rate. C is the strain rate effect parameter. A, B, N, and S max are collectively referred to as limit surface parameters. A is the characteristic viscous strength coefficient, B is the characteristic pressure hardening coefficient, and N is the pressure hardening coefficient.

Damage Model
The damage in the HJC model is formed by the accumulation of plastic strain. Plastic strain has two parts: shear deformation plastic strain (equivalent plastic strain) and volume compressive plastic strain, as shown in Fig. 2. The damage model is expressed as follows: where De p and Dl p are the equivalent plastic strain increment and equivalent plastic volume strain increment of the element in a calculation cycle, respectively e f p and l f p are the equivalent plastic strain and equivalent plastic volume strain at the current integration step. T* is the maximum characteristic tensile pressure that the material can withstand.
where T is the tensile strength of the material. D 1 , D 2 is the material damage parameter. EF MIN is the minimum plastic strain when the material fails, and it is used to control the brittle cracking of the material caused by the tensile stress wave.

State Equation
The state equation of the HJC model is shown in Fig. 3. It is used to describe the relationship between hydrostatic pressure and volumetric strain. It is divided into three stages: elastic compression, compaction deformation, and post-compaction deformation. The first stage (P < P c ) involves the linear elastic zone. Hydrostatic pressure has a linear relationship with the volumetric strain, and the relationship is expressed as follows: The second stage (P c ≤ P ≤ P l ) involves the plastic transition zone. At this time, the voids in the material are gradually compressed to produce plastic deformation. The loading and unloading equations are as follows: The third stage (P > P l ) involves the high-compression zone. At this stage, the material has been completely broken. The loading and unloading equations are as follows: where the meanings of the symbols are listed in Table 1.
K 3 the pressure constants l l the volume strain at the limit of compaction l p the volume deformation

Basic Physical Parameters
A cylindrical specimen with a diameter of 40 mm and height of 80 mm was used to measure density via the volumetric method. Under the same water-bearing condition, each group of test pieces was divided into three pieces, and the density of the test pieces was measured. The average density of each group of samples was obtained as the rock density. The density ρ of dolomite was 2,950 kg/m 3 . Afterward, the specimen was subjected to a uniaxial compression test. On the basis of the test results, the uniaxial compressive strength (UCS) of the specimen was determined to be 159.65 MPa, elastic modulus E was 24.16 GPa, and Poisson's ratio v was 0.238. The final failure form of the specimen is shown in Fig. 4.
The rock specimen had a diameter of 40 mm, height of 20 mm, and a disk shape. It was placed horizontally between the pressure plates of the press. A hard steel wire mat with a diameter of 2 mm was also placed between the upper and lower bearing plates. The cushion strip was perpendicular to the symmetry axis of the test piece. The cushion strip changed the applied pressure into a linear load and made the test piece produce a tensile force perpendicular to the direction of the upper and lower loads. The test piece was damaged in tension along the radial direction (Fig. 5). On the basis of the test result, the tensile strength (TS) of dolomite was 7.84 MPa.
Elastic modulus, shear modulus, and other measures were calculated by using the formula: E is the elastic modulus and r E is the ultimate elastic strength of sandstone. e E1 is the axial strain corresponding to the ultimate elastic strength of sandstone, and e E2 is the radial strain corresponding to the ultimate elastic strength of sandstone. G and K are the shear modulus and bulk modulus, respectively. T* is the maximum characteristic tensile pressure to bear. T is the tensile strength of the material, whose value is obtained via a splitting experiment. f c is the quasi-static compressive strength, which can be obtained through a uniaxial compression test. The basic mechanical parameters of the rock are shown in Table 2. The cohesion force and internal friction can be obtained by the trixial compresion test (Fig. 6).

Strain Rate Effect Parameter
The SHPB experimental system was used to measure the uniaxial dynamic strength of sandstone under different strain rates. The result was combined with the measured quasi-static uniaxial compressive strength to obtain the dimensionless equivalent stress of the sandstone specimen at different strain rates as follows:

Limit Surface Parameters A is the characteristic viscous strength coefficient in
Without considering the damage and rate effects, the strong expression of the HJC model becomes On the basis of plastic yield surface theory, if damage and rate effects are not considered, then the HJC model and the cap model would pass through two points of pure shear and uniaxial compression on the compression meridian plane, as shown in Fig. 7. The following derivation establishes the relationship between the cap model parameters and the characteristic viscosity strength parameters.
The cap model has two stress tensor invariants: the first stress tensor invariant I 1 and the second stress . S ij is the stress deflection, which consists of a dilatancy surface and cap surface. The cap surface can control the yield caused by compression and the dilatancy of the material under shear.
B and N are the characteristic pressure hardening coefficient and pressure hardening coefficient, respectively. A universal testing machine was used to perform a triaxial compression experiment, which r 1 is the axial compression and r 2 ¼ r 3 is the confining pressure. Under different confining pressure conditions, the values of different hydrostatic pressures P and the corresponding maximum differential pressure Dr were measured. Then, the obtained data were normalized to derive a series of points ðP Ã ; r Ã Þ. The data points were fitted to determine the values of B and N.
3.4 Pressure, Damage, and Software Parameters P c and u c are the hydrostatic pressure and volume strain at the elastic limit; P c ¼ f c 3 , and u c ¼ P c K . Based on the Hugoniot relationship expressed by the empirical constant of impact compression, P l , u l , K 1 , K 2 , and K 3 were obtained by fitting. The damage parameter is expressed as D 1 ¼ 0:01 , where a value of 1 is assigned to material damage parameter D 2 . EF MIN is the minimum plastic strain when the material is broken. It is used to control the brittle cracking of the material caused by the tensile stress wave and was given a value of 0.01 in this study. _ e and f s are the software recommended parameters. The final parameters of the HJC constitutive model are shown in Table 3.

Numerical Simulation Research on the Air-Deck Charge of the 2m-High Bench Blasting
The high bench blasting was analyzed using the HJC constitutive model. The explosive and air results employed a multi-material ALE algorithm, and the rock used the LAGRANGE algorithm. The fluid-solid coupling algorithm reflects the stress transfer and attenuation between explosives, air, and rocks.

Model Establishment
The basic principle is to combine the explosive detonation research results with the explosive gas state equation to describe the dynamic expansion of the entire explosion cavity. LS-DYNA can directly simulate the explosion process of high-energy explosives. After igniting the explosive, the explosive load is generated to act on the surrounding medium and the pressure at one point in the explosive source at any time is expressed as where P is the explosion pressure; F is the chemical energy release rate of the explosive; D is the explosive velocity; t and t 1 are the current and detonation times of one point in the explosive. A e max and v e are the maximum cross-sectional area and the explosive unit volume.
With regard to explosive model selection, *MAT_HIGH_EXPLOSIVE_BURN high-performance explosive material was used to simulate the detonation of explosives, several of the parameters are shown in Table 4. The JWL equation of state was used to simulate the entire detonation process. In the numerical simulation of the explosion field, due to the large pressure fluctuation range of the detonation product, a state equation that is suitable for all pressure ranges is difficult to find. In 1965, Kury et al. [36] proposed the JWL equation based on the work of Jones and Wilkins. Its general form is where P eos is the pressure. V is the relative volume; E 0 is the initial specific internal energy; and A, B, R 1 , R 2 , and x are five independent physical constants describing the JWL equation. The parameter values of ammonium explosives are shown in Table 4.  The parameters of the air model (*MAT_NULL) and the state equation (*EOS_LINEAR_ POLYNOMIAL) are shown in Table 5. Air density was taken as 1.225 kg/m 3 , and the pressure was simulated using a linear polynomial equation of state (*EOS-LINEAR-POLYNOMIAL). The equation is expressed as where C 0 -C 6 are constants, l is the specific volume, and E is the ratio of internal energy to the initial volume. The values of the coefficients of the above-mentioned equation are as follows:

Numerical Simulation Results and Analysis of the Air Deck in the Middle Part
The upper and lower explosives were detonated at the same time. After the detonation, the stress wave spread around in the form of a spherical wave. The rocks around the explosion point were destroyed first, resulting in a cavity. The load stress on the blast hole wall at the blasting source reached the maximum value of 1.3 Gpa, as shown in Fig. 8. At 6 ms, the bottom of the columnar charge completed the detonation and met the upper part of the charge in middle air-deck, and the mutual superposition effect was obvious. With the propagation of the detonation wave, the columnar charge reaction was completed at 13 ms, and the stress wave expanded in an ellipsoidal shape. At the same time, the detonation wave generated by the explosion of the ammonium explosive propagated upward and downward along the blast hole, which characterized the detonation transfer process of the ammonium explosives. At 18 ms, the detonation wave was transmitted to the bottom and top of the blast hole, and the wavefront of the explosion stress wave was centered on the blast hole and expanded outward in an ellipsoid shape, as shown in Fig. 8d.   Figure 8: Continued At 25 ms, the second blast hole bagan to detonate, and the initial process of detonation was similar to that of the first blast hole. At 38 ms, the reaction of the second blast hole was completed, and the stress of the first blast hole attenuated to zero; they acted together on the rock between the two blast holes. As the stress spread, the stresses of the two blast holes in the front row was transmitted to all parts of the bench and superimposed, thereby destroying the rock and creating a free surface for the holes in the back row. As shown in Fig. 9.  The two stress cloud diagrams of the front row of holes show that the position where the stresses were superimposed occupied most of the slope of the entire bench, indicating that most of the rock at the slope of the bench was destroyed with a good crushing effect.

Field Blasting Test on a 24 m-High Bench
The design of the high bench blasting used vertical drilling. Given the large density factor, the hole spacing was 18 m, and the row spacing was 6 m. The method of hole-by-hole initiation was adopted, and the interval time difference was 25 ms. The shape of the muck pile after blasting is shown in Fig. 10. The blasting effect was better, and the degree of looseness was higher.

Conclusion
(1) The analysis results obtained from the effective stress-time diagram of the air-deck charge structure showed that when the air-deck was in the hole, the effective stress peak value of the stress monitoring point on the slope surface was greater than the rock tensile strength. The blasting effect of the charge structure with the air deck in the hole produced good results.
(2) The comparison of the stress cloud diagrams and the effective stress of the rock at the observation point indicated that the middle air-deck charge exerted a good effect on different blast holes.
(3) The comparison and analysis results of the three blast holes showed that when the air-desk was in the hole, the gas pressure was balanced, and the action time was prolonged. This condition reduced the use of explosives and improved the bottom rock shattering effect.
(4) The application of the air-deck charge structure in the hole improved the quality of loose blasting, reduced the bulk rate, and decreased the appearance of roots. The cost of bench blasting was reduced and achieved good results in Open-Pit Mine.
Acknowledgement: The author would like to thank the Barun Mining Co., Ltd. of Baotou Steel Group, China.
Funding Statement: The authors received no specific funding for this study.
Conflicts of Interest: The authors declare that they have no conflicts of interest to report regarding the present study.