A Method of Calculating Critical Depth of Burial of Explosive Charges to Generate Bulging and Cratering in

For underground explosions, a thin to medium thickness layer near the cavity of an explosion can be considered a theoretical shell structure. Detonation products transmit the effective energy of explosives to this shell which can expand thus leading to irreversible deformation of the surrounding medium. Based on mass conservation, incompressible conditions, and boundary conditions, the possible kinematic velocity fields in the plastic zone are established. Based on limit equilibrium theory, this work built equations of material resistance corresponding to different possible kinematic velocity fields. Combined with initial conditions and boundary conditions, equations of motion and material resistance are solved, respectively. It is found that critical depth of burial is positively related to a dimensionless impact factor, which reflects the characteristics of the explosives and the surrounding medium. Finally, an example is given, which suggests that this method is capable of calculating the critical depth of burial and the calculated results are consistent with empirical results.


Introduction
Underground explosions have been studied for different purposes [1].They can be classified into contained explosions (or camouflet explosions), bulging explosions, and cratering explosions.Contained explosions occur when the depth of burial of the explosive is sufficiently deep [2] and the model best suited to such conditions is an explosion in an infinite medium [3].Bulging explosions occur when the depth of burial of the explosive is within a certain range [4] and are best modelled as explosions in a semi-infinite medium.Cratering explosions may occur in various ways when the depth of burial of the explosive is relatively shallow and they are best modelled as explosions in a semi-infinite medium [5].A widely accepted zonal model divides the medium near the cavity into a grinding zone, a radial cracking zone, and an elastic zone [6], which can be the basis for the study of bulging and cratering explosions.Owing to their significant military and engineering application, bulging and cratering explosions have arisen widespread concern for a long time.
Actually, bulging and cratering explosions are more complex subjects which should take the influence of the free facet into consideration in dynamic models.
At present, the cratering mechanism of blasting can be divided into three kinds of effects: compressional stress wave effects, tensile reflected wave effects, and gas pressure effects.Since crushing and plastic deformation are omnipresent in the medium surrounding the explosive source, compressional stress wave effects predominate.After Hino [7] proposed the concept of blastability coefficient, reflected tensile wave effects have been widely accepted and increasingly investigated by other researchers.To clarify the role of gas pressure in the fragmentation of an underground blast, Kutter and Fairhurst [8] designed a system of model tests in which the combustion products were simulated by pressurised oil.They concluded that the gas pressure played an important role in blasting.Hagan agreed with this view and coined the term "pneumatic wedging" to describe previous observations [9,10].Meanwhile, many other researchers, such as Dally et al. [11], Dally et al. [12], and Harries [13], verified the reasonability of the aforementioned view in different ways.However, these research findings exhibit fundamental mechanisms of explosive cratering and cannot be used to solve most practical problems.Therefore, researchers have to formulate different hypotheses or conduct many tests to obtain any more practical computational formulae.
Nowadays, the widely used theories for calculating the explosive charge or the depth of burial of cratering explosions include Livingston's crater theory [14,15], Boreskov's formula [16], Langefors' formula [17,18], Vlasov's formula [19,20], and Pokrovskii's formula [21,22].However, these theories cannot be used to calculate critical depth of burial of explosive charges to generate bulging and cratering in rock.Practically, it is important and pressing to investigate the mechanism or principle behind bulging explosions.A guide to bulging blasts in roads when earth-penetrating bombs are used on airport runway or highway pavements is needed.
As shown in Figure 1 [23], "n" denotes throwing index, which is a ratio of the radius of the crater (B) to the depth of burial of explosive (h).For contained explosions,  = 0, shown in Figure 1(a); for bulging explosions, n < 0.75, as shown in Figure 1(b); for cratering explosions, n > 0.75, as shown in Figure 1(c).A critical depth of burial, which forms the dividing line between a contained explosion and a bulging explosion, can be called the first critical depth of burial (FCDOB).Similarly, a critical depth of burial which forms the dividing line between a bulging explosion and a cratering explosion can be called the second critical depth of burial (SCDOB).When the depth of burial is between these two critical depths, a bulging explosion tends to occur.This work focuses on two types of critical depth of burial, which can be used to predict the type of explosion and provide criteria for the occurrence of a bulging explosion.

The Dynamic Model: An Introduction
As is known, scabbing or perforation will occur when a target is subjected to blast loads [24,25] (Figure 2).If the target is a semi-infinite medium and the explosive is detonated in the medium, scabbing or perforation will occur on the free facet and irreversible deformation happens near the blast hole; therefore, the progress is similar to that of a shallow buried explosion (Figure 3).Furthermore, the progress of a bulging explosion is similar to that of scabbing and the progress of a cratering explosion is similar to that of perforation.Therefore, a method [26,27] to calculate scabbing and perforation can be used to investigate the category of an underground explosion; however, the density of the medium, and the effective energy of the explosives, should be considered when establishing any dynamic model of a shallow buried explosion.
Before compressive waves reach the free surface, the processes in a shallow explosion are similar to those in a fully contained explosion.Hence, the initial conditions in a shallow-burial explosion are the same as the state parameters when coupled charges are detonated in an infinite media.
Figure 4 shows the proven results on the zoning of deformation and failure of rock under contained explosions [6].There are four zones: (1) cavity; (2) crushed zone; (3) radial cracking zone; and (4) elastic zone.Besides,  1 ,  2 , and  3 are the cavity, crushed zone, and radial crack zone radii, respectively.However, all of them are posttest results and are measured after any explosion which cannot describe the dynamic processes inherent to explosive cavity expansion.
Regarding a thin layer near the cavity as an impact body [28], the problem of an explosion in rock can be transformed to a problem of impact between the impact body and the surrounding medium.Detonation products transmit effective initial mechanical energy to the surrounding medium through the impact body, which leads to irreversible deformation and movement of the medium.The impact body is a shell whose mass is constant and whose position moves outwards with the cavity-expansion velocity.As supposed by Forrestal and Tzou [29], when the cavity-expansion velocity  is large enough, the radially cracked region diminishes and the response can be depicted by an elastic-plastic model (Figure 5).
Radius r, as a varying displacement of the impact body, is the position of the contact surface between the impact body and the surrounding medium.Therefore, the initial conditions of the impact body can be written as According to the principle of the conservation of energy, the kinematic energy of the impact body is where  denotes a relative effective energy,  V is the specific energy (kcal/kg),  is the mass of the charge (kg), and  = 4187 (J/kcal) denotes the mechanical equivalent of heat.
The equation of motion of the impact body is where M is the mass of the impact body, r is the displacement of the impact body, and P is the deformation resistance of the surrounding medium (material resistance).The material resistance force can be acquired by seeking the limiting load.When the load reaches the limit for an ideal plastic-rigid body (with a yield plateau), free plastic flow occurs.According to extremum principles and energy equilibrium conditions, the rate of work can be written as follows [27,30,31]: where  and , respectively, are the area and the volume of the deformed body,    are components of surface forces, ]  are components of the possible kinematic velocity fields,    = (1/2)(]   /  + ]   /  ) are components of strain rates,    are components of the actual stress,   are surfaces of velocity discontinuity, and [V   ] denotes the jump V  + − V  − on   .By virtue of the Saint Venant-von Mises equation [31], it is found that the vectors    and    are parallel and the quantity       , being the scalar product of parallel vectors, will be equal to the product of their moduli; that is,       =     , where   is the yield limit for shear and   is the intensity of the shear strain-rate.Then (4) can be written as According to (4) or ( 5), the upper bound for the limit load of    can be calculated from the possible kinematic velocity field.Therefore, the problem becomes one of how to accurately establish the possible kinematic velocity field which is based on experimental data.

Material Resistance to an Explosion in an Infinite Medium.
In a spherical coordinate system (Figure 6), a microelement measuring  0 in both  and  directions is cut from the model of a contained explosion.The sketch of a microelement near the cavity depicts the velocity components in the infinite medium velocity field, where V  , V  , and V  are components of the velocity in , , and  directions, respectively;  2 is the initial outer radius of the impact body, and r is the new outer radius of the impact body which forms the contact surface at some point.From  2 to r, the motion of the impact body is resisted by the surrounding medium which is in a plastic state.Then it is essential to find the material resistance.The limiting material resistance can be calculated as follows: (1) Considering the symmetry of the velocity field, V  = V  .According to the law of conservation of mass, that is, , the radial velocity of the medium in the plastic zone is V  =  2  2 / 2 V 0 .According to the incompressibility condition, V  / + 2(V  /) + V  / + V  / = 0, and then V  = V  = 0.
(2) Since the velocity field is known, it is easy to obtain the strain-rate components: (3) The intensity of the shear strain-rate is (4) Substitute the physical and geometric quantities into (5), manipulate algebraically, and then the limiting material resistance to an explosion in an infinite medium is Using (8), material resistance that resists the motion of the impact body generated by the explosion in an infinite medium is calculated.It is a function of the yield limit for shear, the initial position of the contact surface, and the new position of the contact surface at some point.

Material Resistance to Explosions in a Semi-Infinite
Medium.For explosions in a semi-infinite medium, it is essential to consider the effect of the free surface.There are two possible kinematic velocity fields: the first velocity field is a single-variable field where the radial velocity is related to the radial coordinate.It describes the motion when the free surface has been detected by the shock wave before the formation of a perforation body.However, the second velocity field is a constant field where the radial velocity remains constant over the whole deformation zone.It describes the motion after the formation of the perforation body which moves as one mass.
A funnel-shaped plastic region between the impact body and the free surface will be formed when the explosive is not buried deeply.In the funnel range, large displacements will occur as the medium is in an ideal plastic state; but, outside the funnel range, the deformation is restricted since the medium is seen as a rigid body.As shown in Figure 7, ℎ is the depth of burial of the explosive;  1 ,  2 are, respectively, the initial inner radius and outer radius of the impact body; and  0 is the interfacial angle between the plastic zone and the rigid zone.
In the first velocity field, the plastic zone is the cone-shape region shown in Figure 7 whose cone angle is  0 and whose height ranges from  2 to ℎ.The initial zone of the impact body ranges from the internal surface Σ 1 to the external surface Σ 2 .The area of surface Σ 1 cut by the cone is  1 = 2 2  1 (1−cos  0 ); the volume of the cone inside surface Σ 1 is  1 = 2/3 3  1 (1 − cos  0 ); likewise, the area of surface Σ 2 cut by the cone, and the volume of the cone inside surface Σ 2 , can be found.
Assuming that V  , V  , and V  are velocity components in the , , and  directions (Figure 6), the limiting material resistance is calculated as follows: (1) The initial velocity of the impact body is such that its average velocity is According to the law of conservation of mass, the radial velocity is Considering the symmetry of the velocity field and the incompressibility condition, the compatible equations of this velocity field are V  /+ 2(V  /)+V  /+V  / = 0 and V  = 0. Combined with the condition The strain-rate components are (3) The intensity of the shear strain-rate is (4) According to (5), the limiting material resistance to an explosion in a semi-infinite medium is Equation ( 11) is solved at a certain time, at which the material resistance corresponds to the initial position of the contact surface.However, during the whole explosion process, the impact body undergoes outward movement and the position of the contact surface develops accordingly.It can be assumed that the position of the contact surface is r at some point in the dynamic process, which is also the new outer radius of the impact body.Therefore, at some point, the limiting material resistance to explosions in a semi-infinite medium is calculated as However, in the second velocity field (Figure 8), the cone-shape region undergoes rigid body motion.The velocity within the whole deformation zone remains unchanged, which indicates that the rigid block is pressed out from the medium.
Initially, the velocities of the rigid block and the impact body are different; therefore there remains a resistance to the effects of the explosion.As the motion develops, the rigid block and the impact share the same velocity; the resistance to explosion therefore becomes zero.This physical process characterises the block motion after the formation of a perforation.
After mathematical manipulation, the limiting material resistance to an explosion can be written in the similar form to that of (11): Likewise, at some point, the limiting material resistance to explosions in a semi-infinite medium is

Calculation Method: Critical Depth of Burial
4.1.The Nondimensional Material Resistance.According to the literatures [1,4,[32][33][34], the cone angle of a standard blasting crater is /4.Therefore, it is assumed that the interfacial angle  0 = /4.If the position of the contact surface r is known, the material resistance in three velocity fields can be, respectively, acquired through ( 8), (12), and ( 14).Consequently, the relationships between material resistance and the contact surface position can be acquired.These are known as material resistance curves.Equations ( 8), (12), and ( 14) can be transformed to nondimensional forms.Assuming  0 is the radius of the explosive charge, for coupled charges, it is also the radius of the blast hole.Then the depth of burial ℎ = h ⋅  0 , the position of contact surface  = r ⋅  0 , and the initial external radius of the impact body is  2 = ã2 ⋅  0 .
From (8), the dimensionless material resistance in the velocity field in an infinite medium is From ( 12), the dimensionless material resistance in the first velocity field of a semi-infinite medium is written as From ( 14), the dimensionless material resistance in the second velocity field of a semi-infinite medium is Under the assumption that ℎ = 10 0 , there are three function curves   /  ∼ r (Figure 9): curve 1 is the material resistance curve in the first velocity field of a semi-infinite medium, curve 2 is the material resistance curve in the second velocity field of a semi-infinite medium, and curve 3 demonstrates the material resistance in the velocity field of an infinite medium.In Figure 9, curves 1 and 3 intersect at A, at which point the impact body begins to perceive the presence of the free surface when the position of the contact surface  reaches   .Besides, curves 1 and 2 intersect at B, which indicates that the impact body pushes the mass block out of the surface when the position of the contact surface reaches   .In the whole physical process, the material resistance is the minimum value among those on the three material resistance curves.To the left-hand side of A, the material resistance is given by curve 3; to the right-hand side of A, the material resistance is given by curve 1.
Since   is the force per unit area, the total force on the whole contact surface is where   is the area of the contact surface, in an infinite medium,   = 4 2 , and, in a semi-infinite medium,   = 2 2 (1 − cos  0 ).

The First Critical Depth of Burial.
According to ( 15) and ( 16), the total material resistance in an infinite medium  =   ⋅ 4 2 0 r2 .Therefore, based on (3), the equation of motion in an infinite medium is The initial conditions are calculated as Denoting 8 √ 3   0 / by  1 , the equation of motion may be transformed to As the differential equation cannot be solved directly, its transformed form is After separation of variables and integrating, After simplifying, Denoting (4.5/ 1 )(V 0 / 0 ) 2 by , a dimensionless impact factor, by solving (24), we obtain where () represents the Lambert W-function, and when the independent variable  > − −1 , the function increases progressively, and the function result is a real single value.Solving simultaneous equations ( 15) and ( 16) gives the medium thickness corresponding to point A in Figure 9.This is also the depth of burial of the explosive and can be written as By substituting ( 25) into (26), the minimum depth of burial of explosive required to prevent scabbing becomes From ( 2), the form of the dimensionless impact factor becomes The dimensionless impact factor can also be expressed as follows: It is a dimensionless combination of multiparameters, which characterises the relationships between the explosive density ( 0 ), specific energy ( V ), the effective energy ratio (), medium density (), yield strength (  ), peak pressure (  ), plastic wave velocity (  ), and maximum particle speed (V  ).
The initial impact body is a thin shell around the cavity, and its outer radius  2 =  0 (1 + ), where  ≪ 1 denotes the thickness of the shell.Therefore, ã2 , the dimensionless outer radius of the initial impact body, is 1 + , and it is approximately equal to 1. Consequently, (27) becomes When using a coupled charge,  0 is the radius of the explosive charge, and it is also the radius of the blast hole.Since the explosive density is assumed to be  0 , the mass of charge  = (4/3) 0  3 0 .Then the relationship between the depth of burial and the mass of the charge is Therefore Equation ( 32) denotes the scaled critical depth of burial for bulging, where the impact factor  and the dimensionless critical depth of burial h are determined by ( 28) and (30).Equations ( 30) and (32) give the mass of the charge for bulging or scabbing as  = (4 0 /3)/(()) 3 ⋅ ℎ 3 .Therefore, (4 0 /3)/(()) 3 is the functional expression of empirical coefficient   , which is the critical coefficient for scabbing.

The Second Critical Depth of Burial.
To calculate the second critical depth of burial, this research uses the initial conditions consisting of the displacement and velocity corresponding to those at point A in Figure 9.At point A, the impact body is affected by the existence of the free surface and the perforation begins to form.Then the material resistance (in motion) can be calculated from curve 1 between A and B, which corresponds to the first velocity field in a semi-infinite medium.When the position of the contact surface r reaches point B, the perforation ends and the cratering is complete.Therefore, the medium thickness h corresponding to point B is the critical depth of burial for cratering.
According to (11) and ( 13), the material resistance in the first velocity field of a semi-infinite medium is  =   ⋅ 2 2 0 r2 (1 − cos  0 ), and the mass of the impact body is  = (2/3)(1 − cos  0 )(r 3  − 1) 3 0 ; consequently, the equation of motion of the perforation is transformed to Combining the initial conditions r| =0 = r and ṙ | =0 = V  / 0 and denoting 17.64  /(r 3  − 1) 2 0 by  2 , the equations of motion are transformed to This differential equation cannot be solved directly, and its transformed form is After separation of the variables and integration, The initial conditions of perforation were given by ( 21), so According to the initial conditions on the impact body, the aforementioned equation may be transformed to Solving ( 36) and (38) simultaneously gives, after algebraic manipulation.
where  1 = 8 √ 3   0 / = 8 √ 3  /(4/3 ⋅ (r 3 − 1) 2 0 ) and  2 = 17.64  /(r 3  − 1) 2 0 .Since r = r as an initial condition, From point A (Figure 9), the relationship between the position of the contact surface and the medium thickness can be written as: From point B (Figure 9), the relationship between the position of the contact surface and the medium thickness can be written as Combining (39) with (42), which include equations of motion, boundary conditions, and initial conditions, an equation for the critical depth of burial for cratering is obtained: where h = () denotes the critical depth of burial for cratering and  denotes the dimensionless impact factor.However, (43) cannot be solved by using an analytic method, and there is no choice but to pursue an approximate solution by numerical method.Then, substituting the solution into (31) and manipulating algebraically, the scaled critical depth of burial for cratering can be written in the same form as (32).Similarly, from (32) and (43), the mass of the charge required for cratering or perforation is given by  = (4 0 /3)/(()) 3 ⋅ ℎ 3 .Therefore, (4 0 /3)/(()) 3 is the functional expression of empirical coefficient, and   is the critical coefficient for perforation.

Calculation of Critical Depths of Burial
Assuming that the key variables [35] are known, they consist of explosive density  0 = 1630 kg/m 3 , the effective energy ratio  ≈ 0.8, specific energy  V = 1010 kcal/kg, mechanical equivalent of heat  = 4187 J/kcal, medium density  = 2600 kg/m 3 , and material yield strength which ranges from 6 MPa to 2 MPa.
In Figure 10, as the material yield strength decreases, the dimensionless impact factor  increases, and the critical depth of burial ℎ/ 1/3 rises.It indicates that materials with higher yield strengths confer a greater blast-resistance and, therefore, have a smaller explosive critical depth of burial requirement to prevent bulging or cratering.These results are consistent with empirical conclusions.Some cases involving each critical depth of burial are summarised in Table 1.
From Table 1, the FCDOB ranges from 0.8 m/kg 1/3 to 1.4 m/kg 1/3 , and the SCDOB ranges from 0.5 m/kg 1/3 to 0.7 m/kg 1/3 .Actually, the scope of underground explosions can be empirically divided into three types, whose critical depths of burial are listed as follows [36,37]: those whose FCDOB ranges from 1 m/kg 1/3 to 2 m/kg 1/3 and whose SCDOB ranges from 0.4 m/kg 1/3 to 0.6 m/kg 1/3 .The critical depths of burial calculated by using the proposed theoretical method are only approximately equal to the empirical values found.In addition, two types of critical depth of burial have divided the scope of underground explosions into three regions-a, b, and c-which correspond to contained explosions, bulging explosions, and cratering explosions, respectively.

Conclusion and Recommendations for Future Research
This research transformed the problem of an explosion to one of an impact between a moving body and the surrounding medium.Then it derived the material resistance to explosive loading using limit load theorems and solved the equation of motion by combining the initial conditions and the boundary conditions.Finally, the research found two types of critical depth of burial for bulging and cratering and proved the correctness of the method through use of a case study.On the basis of the above calculations and analysis, the following was concluded: (1) Two types of critical depth of burial are both positively related to the dimensionless impact factor.The dimensionless impact factor  reflects the properties of the explosive performance and material strength, and it also presents the blasting ability of an explosive on its surrounding materials.The higher the explosive parameters, the more effective the explosive and the bigger the dimensionless impact factor; therefore, the critical depth of burial is greater.On the contrary, the higher the material yield strength, the stronger the material resistance to explosion and the smaller the dimensionless impact factor; therefore, the critical depth of burial is smaller.
(2) In the ℎ/ 1/3 ∼  plane, two critical depth of burial curves divide the plane into three regions corresponding to contained explosions, bulging explosions, and cratering explosions, respectively.Besides, the critical depths of burial calculated by the proposed theoretical method are approximately equal to the empirical results found.
(3) Based on the comparison with empirical results, the theoretical calculation method is verified to be reasonable.Besides, the theoretical method explains the physical essence of the empirical critical coefficients.
Since the critical depth of burial is a theoretical value, it is too difficult to implement an explosion which is located exactly on the critical depth of burial curve locus.At the same time, although there is much information about underground explosions, some key parameters, such as the effective energy ratio and material yield strength, are usually not given in most references.Therefore, it is difficult to verify the precision of the theoretical method.However, some verification tests,

Figure 4 :
Figure 4: The zoning of deformation and failure of rock under a contained explosion.

Figure 5 :
Figure 5: Outward motion of the impact body.

Figure 6 :
Figure 6: Velocity components in the velocity field of an infinite medium.

Figure 7 :
Figure 7: The first velocity field of explosions in a semi-infinite medium.

2 Figure 8 :
Figure 8: The second velocity field of an explosion in a finitethickness-medium.
The first velocity field in semi-infinite medium The second velocity field in semi-infinite medium The velocity field in infinite medium

Table 1 :
Calculated results: two types of critical depth of burial.