Considering explosive charge shape and embedded depth in the design of concrete shelter thickness

Yi Fn , Li Chen , Jin Hong , Runqing Yu , Hengbo Xing , Qin Fng

a State Key Laboratory of Disaster Prevention & Mitigation of Explosion & Impact, Army Engineering University of PLA, Nanjing, Jiangsu, 210007, China

b Engineering Research Center of Safety and Protection of Explosion&Impact of Ministry of Education,Southeast University,Nanjing,Jiangsu,211189,China

c The People's Liberation Army 95338 Troops, Guangzhou, Guangdong, 510403, China

Keywords:Blast resistance Charge shape Embedded depth Structural design

ABSTRACT Cost and safety are important considerations when designing the thickness of a protective reinforced concrete shelter.The blast perforation limit(BPL) is the minimum concrete shelter thickness that resists perforation under blast loading.To investigate the influence of the depth of embedment (DOE) and length-to-diameter ratio(L/D) of an explosive charge on the BPL, the results of an explosion test using a slender explosive partially embedded in a reinforced concrete slab were used to validate a refined finite element model.This model was then applied to conduct more than 300 simulations with strictly controlled variables,obtaining the BPLs for various concrete slabs subjected to charge DOEs ranging from 0 to ∞and L/D values ranging from 0.89 to 6.87.The numerical results were compared with the experimental results from published literature, further verifying the reliability of the simulation.The findings indicate that for the same explosive charge mass and L/D, the greater the DOE, the larger the critical residual thickness (Rc, defined as the difference between the BPL and DOE) up to a certain constant value;for the same explosive charge mass and DOE,the greater the L/D,the smaller the Rc.Thus,corresponding DOE and shape coefficients were introduced to derive a new equation for the BPL,providing a theoretical approach to the design and safety assessment of protective structures.

1.Introduction

When struck by a precision-guided munition (PGM), the structure of a protective concrete shelter is subjected to the combined effects of damage from penetration and explosion, as shown in Fig.1.If the shelter is too thin, the PGM will penetrate through it before exploding, causing severe damage to the personnel and equipment within [1].If the shelter envelope is too thick, it will withstand any damage effects from the PGM, but at significantly increased construction difficulty and cost.Therefore,it is important that the correct thickness be chosen for the envelope of the protective structure during design.

The combined effects of penetration and explosion damage should be considered in the design of the shelter envelope thickness.Many factors influence shelter design, including the content and characteristics of the concrete materials and any reinforcement(e.g., rebars and fibers), the characteristics of the projectile body(i.e., material, shape, and weight), the parameters related to the projectile strike (e.g., speed and angle), and the characteristics of the explosive(e.g.,material,weight,shape,and initiation point).In recent years, advances in science and technology have resulted in an increase in the length-to-diameter ratios (L/Ds) of PGM warheads to values in excess of 10.In addition, the continuous development of earth-penetrating weapons has progressively increased their potential penetration depth.For example, the tandem warhead contains a forward shaped charge and a follow through warhead to induce deeper damage in concrete structures [2].As a result of their excellent penetration performance, tandem warheads are widely used for attacking aircraft shelters, underground fortifications, and other important targets.Owing to the variety of factors influencing the penetration and explosion processes, they are often studied separately.As a result, when studying the explosive damage to a concrete slab, a blast hole can be prefabricated to represent the penetration depth of the projectile.

Fig.1. Structure for the protection from PGM.

The damage to a concrete shelter when subjected to blast loading provides an important basis for the design of the shelter envelope thickness.A particularly important index during this design stage is the blast perforation limit (BPL).Researchers have proposed different calculation methods to determine the BPL for concrete slabs based on data obtained from contact explosion tests in which the depth of embedment (DOE) of the explosive charge was 0 [3-9].These previously proposed equations all account for the thickness of the concrete slab and the mass of the explosive charge, but other influencing factors are rarely considered.For example, the shape of the explosive was considered only by the United Facilities Criteria (UFC) [5]and Remennikov et al.[6,7],though it has been found to exert an important influence on the destruction of concrete structures.Indeed, in close-range blast tests,Chen et al.[10]and Hu et al.[11]found that the reflected peak overpressure and impulse generated by the detonation of a cylindrical explosive charge at both ends were 3.5 and 2 times higher,respectively,than those generated by the detonation of a spherical explosive charge of the same mass, and the damage to a concrete column was significantly more severe.Shi et al.[12]conducted close-range blast tests, observing that for detonations of the same cylindrical explosive charge mass, a largerL/Dcorresponded to more severe spalling damage to the concrete slab.Using theoretical analyses and experiments, Orlenko [13]found that when a constant-diameter cylindrical explosive charge was detonated with its axis perpendicular to the surface of a metal plate and the detonation initiation point at the far end, the resulting impulse acting on the plate increased with the length of the explosive up to anL/Dof 2.25.Therefore, the prediction of BPL must consider the influence of the explosive charge shape.

In addition, it should be noted that the axis of the slender cylindrical explosive charge scenario addressed in UFC [5]is parallel to the concrete slab.Therefore,despite the consideration of charge shape in the UFC method,it is unsuitable for scenarios in which the axis of a slender cylindrical explosive charge is perpendicular to the concrete slab, such as the PGM penetration scenario described above.Remennikov et al.[6,7]found that for a constant explosive charge diameter, the impulse acting on the concrete slab did not increase onceL/Dwas greater than 0.5.This ratio is much smaller than theL/Dof 2.25 identified by Orlenko[13]because Remennikov et al.[6,7]ignored the influence of the detonation initiation point[14], assuming instead an instantaneous detonation.In other words, they did not consider the propagation of the blast wave through the explosive charge,which results in a different scenario from the single or multi-point detonation commonly applied in PGM warheads.

Finally, present methods for the calculation of BPL based on contact explosions have failed to consider the influence of the explosive charge DOE.According to the TM5-855-1 manual (TM5)[15]and Mu et al.[16], the greater the DOE in a semi-infinite concrete medium, the greater the ground shock at a given distance from the center of the explosive charge.Furthermore,Xu et al.[17]and Lai et al.[18]conducted explosion tests on concrete slabs with finite thicknesses,finding that the greater the DOE,the greater the damage to the slabs.The US Army Explosives and Demolitions manual(EDM)[19]provides a method for calculating the explosive mass required to penetrate a finite-thickness concrete slab.It suggests that the destruction caused by an explosive with a DOE at mid-thickness of a concrete slab can be 1.8 times stronger than that caused by an equivalent contact explosion.However, this method fails to consider the influence of the shape of the explosive charge.Duan et al.[20]conducted tests using fully embedded cylindrical explosive charges with a maximumL/Dof 10 in plain concrete thick slabs.After analyzing the resulting data, they obtained a relationship between the BPL andL/D.Nevertheless, the method they proposed is not applicable to cases in which the DOE is small,as the DOE values in their tests were all greater than half the thickness of the slab.

In previous research, the influences of the explosive charge shape and DOE on the blast damage to concrete were typically studied separately.Indeed, no previously proposed calculation method has accounted for both the influence of the DOE and the shape of the explosive on the BPL,mainly because of the difficulties in testing multiple variables.Blast testing is costly,dangerous,and time-consuming.Furthermore, for a given explosive charge shape and DOE, a concrete slab of a certain thickness may or may not be damaged by perforation in a given test.To obtain the BPL, it is therefore necessary to pre-design and manufacture a large number of slabs with various similar thicknesses.Because of the difficulties associated with conducting such a large number of blast tests,previous researchers have only been able to obtain limited data describing the BPL.In addition, the results from different studies cannot be directly compared considering differences in the properties of the concrete materials and reinforcement used.Indeed,even when using concrete slabs from the same batch, the results can vary slightly owing to the heterogeneity of the concrete material, which can be problematic for BPL research.

Nevertheless, with the advent of high-speed, low-cost computational techniques, finite element software has been increasingly used to study the damage to concrete structures subjected to blast effects and provide much richer information than typically available from experimental results.This study accordingly constructed a refined finite element model (FEM) using a modified Karagozian& Case (K&C) concrete material model that was validated using results obtained from a physical testing program.Over 300 simulated blast tests with strictly controlled variables were then conducted to obtain a large quantity of data describing the BPL of concrete slabs according to DOE values from 0 to ∞andL/Dvalues from 0.89 to 6.87.A BPL calculation method was then proposed based on the results considering the influence of both the DOE and charge shape.

The remainder of this paper is structured as follows:the FEM is developed in Section 2 and is then validated using blast testing results in Section 3; in Section 4, the numerical simulation results are compared with experimental results in published literature for further verification of the FEM;the verified FEM is then applied to conduct a systematic analysis of the influence of DOE andL/Don the BPL in Section 5 and thereby obtain an empirical equation for the latter;finally, concluding remarks are presented in Section 6.

2.Numerical simulation test program

2.1.FEM of the concrete shelter

As stated in the introduction, the effects of an explosion on a concrete slab are influenced by various factors.Different explosive charge shapes and DOEs will directly affect the intensity and distribution of energy transmitted into the concrete upon detonation,thus affecting the BPL.The purpose of this study was therefore to analyze the influence of the explosive charge shape and DOE on the BPL of a finite thickness concrete slab.As a result,concrete material characteristics were not considered in this study;that is,rebars and fibers were not considered in the numerical simulations, and the parameters of the concrete material model remained unchanged.

A one-quarter FEM of the concrete slab was constructed in LSDYNA with the geometric dimensions and boundary conditions set as shown in Fig.2, whereHdenotes the total thickness of the concrete slab, DOE refers to the depth of embedment of the explosive, and the residual thickness of the concrete slab is calculated as:R=H-DOE.The TNTcharge was modeled as a rectangular prism with a length L and a square base side length a,and its point of detonation was at its geometric center.In this model, the Arbitrary Lagrangian Eulerian (ALE) element (modeling the TNT and Air) overlapped with the Lagrange element (modeling the concrete); the coupling between the ALE element and the Lagrange element was defined using the keyword*CONSTRAINED_LAGRANGE_IN_SOLID to simulate the damage effects of the blast wave on the concrete.

Using LS-DYNA,Yu et al.[1]previously investigated the damage to concrete slabs owing to the explosion of an embedded cylindrical explosive charge.Through an analysis of element convergence,they set the sizes of the concrete, air, and explosive elements to 5 mm,and found that the numerical simulation results were in good agreement with the experimental results.Based on a numerical analysis of contact explosions on concrete, Li and Hao [21]and Hong et al.[22]also found 5 mm to be a suitable element size for describing the local failure of concrete.Therefore, all elements in this study were 5 mm in size.

It should be noted that when setting the transmission boundary on the side of the concrete slab, the elements of the outermost surface of the concrete on that side were placed into a separate PART.An elastic material model was assigned to this PART by defining the keyword *MAT_ELASTIC with the same density, Poisson’s ratio, and Young’s modulus as the concrete material.A transmission boundary was then applied to the outer surface of this PART using the *BOUNDARY_NON_REFLECTING setting.Thus, the elastic PART was co-noded with the concrete material elements along the edge of the slab.This method of using an elastic PART to simulate an infinite domain has been successful in practice,as it can suppress reflected waves better than directly applying a transmission boundary to the surfaces of the concrete edge elements.

2.2.Modified K&C model

It is necessary to select suitable concrete, TNT, and air material models for the numerical simulation of concrete slab failure modes.The most commonly used concrete material models are the Holmquist-Johnson-Cook (HJC), Riedel-Hiermaier-Thoma(RHT),and K&C models[23].As the spalling behavior of a concrete slab is mainly caused by the reflection of tensile waves off the free surface [24], the accuracy of the tensile properties in the concrete material model is critical.However, the HJC and RHT models are inadequate when describing the tensile properties of concrete[25-27], making it difficult to effectively describe the spalling behavior of concrete slabs.

To date,many researchers have simulated the spalling behavior of concrete slabs under blast loading using the K&C model and obtained reasonable prediction results [21,28].However, the current K&C model still has three main shortcomings:1)the automatic generation of the dynamic increase factor in tension(DIFt)is low;2)the tensile fracture strain increases with the strain rate though test results suggest that it should be constant at different strain rates;and 3)there is no softening stage in the stress-strain curve under triaxial tension.

In view of the above shortcomings,Hong et al.[22]improved the K&C model in three aspects: “the tensile damage accumulation,relationship between the yield scale factor and damage function,and the dynamic increase factor for tension.”Hong et al.[22]then simulated a contact explosion using the improved model.The simulation results for surface compressive damage and the spalling on the back of the concrete slab were in good agreement with equivalent test results, verifying the validity and accuracy of the improved model.

This study therefore adopted the improved K&C model proposed by Hong et al.[22]using the main material parameters listed in Table 1.The physical fracture or crushing of the concrete wassimulated using an erosion algorithm [29]with the erosion criterion governed by the tensile damage.According to the modified definition of tensile damage, an element was deleted when the tensile damage reached 85 times λm;i.e.,when the fracture tensile strain reached 0.01.For more details of the material model, please refer to Hong et al.[22].

Table 1Main material parameters for the concrete model.

Fig.2. Numerical model of an embedded explosive detonation in a concrete slab (mm).

Finally, *MAT_HIGH_EXPLOSIVE_BURN was used for the material model of the TNT,and*MAT_NULL was used for the air material model.The parameters of these two models were set as the same as in Hong et al.[22].

2.3.Numerical simulation test scheme

The variableRcis defined as the critical residual thickness of the concrete slab, or the minimum residual thickness required to prevent the slab from being perforated by the explosion.It is calculated asRc= BPL - DOE.As can be seen in Fig.3(a), when the residual thickness R ＜Rc, the ejection crater formed by the explosion connects with the spalling crater, thereby completely perforating the concrete shelter.In Fig.3(b), R = Rc; in this criticality case, the concrete shelter is exactly not perforated, and the corresponding slab thickness H=BPL.When R>Rc,as shown in Fig.3(c),there is no perforation of the concrete shelter at all, as the ejection crater does not connect with the spalling crater.

Five explosive charge sizes were evaluated in the simulations,as shown in Table 2.The analysis was conducted using two different masses of TNT (C),100 g and 200 g, andL/Dvalues ranging from 0.89 to 6.87.To calculateL/D, the square base of the explosive was converted to a circle of equal area, and the diameter of this circle was used asD.

The typical numerical simulation results shown in Fig.3 correspond to explosive charge no.1 in Table 2 with a DOE=35 mm;the slab thicknesses in Fig.3(a)-3(c) are 165 mm, 175 mm, and 195 mm, respectively.For a given explosive charge size with aparticular DOE, the value of R was increased from perforation, as shown in Fig.3(a), by increasing the slab thickness until criticality was reached, as shown in Fig.3(b), at which time R = Rcand H =BPL, where BPL is the minimum concrete shelter thickness that resists perforation under blast loading.The DOE was then gradually increased from 0 until Rcno longer increased for a given explosive charge size,at which time the calculation was stopped.These steps were repeated for the different explosive charge parameters shown in Table 2.In the end, more than 300 simulations were performed using this procedure.

Table 2Explosive charge parameters used in the simulations.

3.Field test and FEM validation

3.1.Field test

Fig.3. Schematic and typical simulated cross-sections of concrete slabs exhibiting: (a) Perforation (R ＜Rc), (b) Criticality (R = Rc), (c) Spalling only (R > Rc).

Fig.4. Field test: (a) Setup; (b) Photograph (mm).

In order to verify the FEM introduced in Section 2,blast testing was performed as shown in Fig.4.The explosive charge used in the field test was a 50×50×500 mm3rectangular prism-shaped TNT explosive charge with anL/Dof 8.86,mass of 2000 g,and a DOE of 230 mm.The detonation point was located in the geometric center of the charge.To meet the requirements of the research project and accommodate the limitations of the test equipment, the target RC slab was square in shape.To minimize boundary effects,the desired ratio of slab width to explosive charge diameter was determined to be 31.9 with reference to relevant test literature[17,18,20];thus,the length and width of the target RC slab were both 1800 mm.The 650 mm thick concrete slab was placed horizontally and supported 700 mm above the ground by steel supports on two opposite sides with a 100 mm overlap,as shown Fig.4(b);the remining two sides were unsupported.A 75 mm diameter cylindrical hole, 230 mm deep was prefabricated at the center of the top surface of the slab by installing a sealed tube blockout prior to casting.The concrete material had an average compressive strength of 47.7 MPa.Each cubic meter contained 417 kg of normal Portland cement (P·O 42.5), 582 kg of sand, 175 kg of water, and 1237 kg of coarse aggregate gravel sized 5-30 mm for a water-cement ratio of 0.42 and sand percentage of 32% by mass.The concrete slab was reinforced with 8 mm diameter threaded rebar (HRB400) with a yield strength of 400 MPa.These bars were spaced at 100 mm in both horizontal directions and placed in three layers equally spaced vertically with a cover layer of 30 mm.In addition,vertical hooked rebars were placed at every other intersection of the horizontal rebar layers, corresponding to a spacing of 200 mm.

3.2.Validating an FEM of the field test

The FEM shown in Fig.5 was constructed using the same dimensions as the test slab introduced in Section 3.1, while the remaining settings (e.g., material parameters, etc.) were the same as the model developed in Section 2 and shown in Fig.2.This FEM was used to verify the proposed modeling approach using the results of the field test.As rebar were present in the test slab,it was necessary to add them to the FEM using beam elements with an elastic-plastic material model (*MAT_PLASTIC_KINEMATIC).The rebar had a density of 7830 kg/m3, a modulus of elasticity of 200 GPa,and a Poisson’s ratio of 0.29.The bond between the rebar and the concrete was simulated using the co-node method.

3.3.Comparison of test results and verification of FEM

The field test results are shown in Fig.6.Fig.6(a) shows the crater on the top surface of the slab, which resembled a 300 mm deep funnel.The sides of the funnel were shaped like stair steps,with different steps corresponding to different diameters.The damaged surface was not smooth, so a range of diameters was recorded for each step.The same damage pattern was also observed in the numerical simulation, as shown in Fig.7(a).The simulation results indicated that the first step from the bottom of the crater had the smallest diameter of 165 mm; this damage was caused by the compressive force of the explosion acting on the concrete.The diameters of the second and third steps were 510 mm and 1020 mm, respectively; these were caused by crack penetration along the horizontal and vertical directions.Comparing Fig.6(a)and Fig.7(a), it can be observed that the shapes and sizes of the craters in the field test and numerical simulations were basically the same.

It can be seen from Fig.6(b)that spalling occurred on the bottom surface of the slab, causing some fragments to fall off.However,many fragments on the sides of the spalling crater in the field test did not detach completely from the concrete slab owing to the mechanical interlock between the fragments and the rebar [30].Fragments that were apparently loose but remained attached to the bottom surface were mostly distributed in a circular area with a diameter of 300-400 mm, which is very similar to the numerical simulation results in Fig.7(a).

Fig.5. FEM used for verification of the field test (mm).

Fig.6. Field test results: (a) Top surface; (b) Bottom surface; (c) Side faces of the slab (mm).

Fig.7. Numerical simulation results obtained using field test FEM: (a) Section plane, (b) Side face (mm).

Obvious vertical cracks with a maximum width of 10 mm can be observed on the sides of the concrete slab in Fig.6(c).There were two main cracks on each side of the slab, basically symmetrical about its centerline and separated by 117-354 mm.These cracks were caused by the complex reflection and superposition of stress waves in the concrete slab[31].The morphologies of the cracks on the four sides are different; this could be caused by the different support conditions on the bottom surface or by the uneven distribution of the concrete material.Comparing Fig.6(c)with Fig.7(b),the simulation results are basically consistent with the field results in terms of the sizes and the distribution of the cracks.

In summary, a comparison of Fig.6 and Fig.7 shows that the numerical simulation results are in good agreement with field test results, suggesting that the FEM constructed and applied in this study is reliable.

3.4.Influence of boundary condition

To verify and demonstrate the validity of the transmission boundary applied on the side of the concrete slab in the FEM,models with and without this boundary were evaluated based on the model presented in Section 3.2.The calculation results with the transmission boundary are shown on the right-hand side of each image in Fig.8,and the calculation results with the unconstrained boundary condition are shown on left-hand side of each image.It can be seen in Fig.8 that the overall damage to the concrete slab was significantly reduced when applying the transmission boundary,indicating that the tensile waves reflected off the lateral free boundary exerted considerable influence on the damage effect.According to the section plane in Fig.8(a), when the transmission boundary was applied, the extension of the horizontal cracks remained basically unchanged, but the number of vertical cracks decreased such that they nearly disappeared.As a result, the diameters of the ejection and spalling craters were significantly reduced,indicating that the horizontal tensile waves reflected off of the unconstrained boundary on the side of the slab increased the diameters of the ejection and spalling craters.The side view in Fig.8(b) shows that the vertical cracks disappeared when the transmission boundary was applied, also suggesting that the vertical cracks observed on the side of the slab in the test were caused by the reflection of horizontal tensile waves off of the unconstrained boundary.

Boundary conditions indeed have an important influence on the explosion damage to concrete slabs.Therefore,when designing an explosion or blast test, the width or diameter of the concrete slab should be sufficiently large or the effect of the tensile wave reflection off the unconstrained boundary should be minimized by other means, especially when testing with embedded explosives.

4.Comparison of numerical simulation results with experimental results from published literature

4.1.L/D ≤1

Fig.9 shows the slab damage induced by the detonation of an explosive charge with parameter set no.2,as defined in Table 2,at criticality(R=Rc)using three different DOEs.It can be observed in the figure that as the DOE increased, the value ofRcincreased, as did the diameters of the ejection and spalling craters.This occurred because more energy acted on the concrete after the explosion as the DOE increased [15].

In order to further evaluate their reliability, the numerical simulation results obtained in this study were compared with test data and the results of prediction methods provided in previous literature.Note that only cases in which DOE=0 andL/D≤1 were considered in this comparison.

Fig.8. Influence of boundary condition(the right side of each image shows the transmission boundary condition and the left side shows the unconstrained boundary condition)on the: (a) Section plane; (b) Side face.

Fig.9. Numerical simulation results with m=200 g,L/D=0.89,and(a)DOE=0 mm,Rc = 190 mm; (b) DOE = 250 mm, Rc = 275 mm; (c) DOE = 475 mm, Rc = 305 mm.

Fig.10. Comparison of contact explosion numerical simulation results with the results obtained using different prediction methods (the purple region indicates that the higher the compressive strength,the lower the BPL for the same explosive charge mass as per UFC [5]).

The previously published test data used in this paper are displayed as hollow and solid points in Fig.10.Ichino et al.[32]and Ohkubo et al.[33](hollow points) used plain concrete without rebar or fibers, as in the FEM presented in this paper.The other studies[4,8,34-37](solid points)all used reinforced concrete.As it can be quite challenging to obtain an accurate BPL from a physical test, it was necessary to process these test data.For the same test conditions(e.g.,reinforcement method,concrete strength,etc.)and explosive charge mass, if one slab was perforated while another slab with the same thickness was not,this thickness was considered to be the BPL.For the same explosive charge mass and similar thicknesses,if the thinner slab was perforated and the thicker slab was not, the average of the two thicknesses was taken as the BPL.For two slabs of the same thickness, if perforation occurred when the explosive charge mass was large and not when it was small,this thickness was taken as the BPL, and the average explosive charge mass was taken as the corresponding mass.It should be noted that in order to ensure reasonable data,the difference between adjacent thicknesses or masses was limited to less than 100% when processing test data.In addition, the mass of each explosive charge type used was converted into equivalent TNT mass using the ratio of its heat of detonation to that of TNT as the equivalence coefficient[5].Two of the previous studies[35,36]shown in Fig.10 used cubeshaped explosives, while the remaining studies used cylindrical explosives(L/D≤1).Note that whenL/D≤1,the influence ofL/Don the perforation of the concrete slab was ignored[5].

The three most commonly used methods for calculating the BPL under contact explosion conditions(denoted as BPL0)were initially presented by McVay[3], Morishita et al.[8], and UFC [5].They are all suitable for reinforced concrete slabs.According to McVay [3],BPL0(in m) can be calculated by

whereRis the distance between the center of the explosive and the surface of the concrete slab(in m)andCis the equivalent TNT mass(in kg).

The equation for BPL0(in m)presented by Morishita et al.[8]is

According to UFC [5], BPL0(in ft) can be calculated as follows:

The results obtained using the three prediction methods given by Eq.(1)-Eq.(3)are plotted in Fig.10.The UFC[5]approach given by Eq.(3) considered the effect of the compressive strength of the concrete and therefore serves as the purple region, in which the higher the compressive strength, the lower the BPL for the same explosive charge mass.It can also be seen from the figure that the curves for the McVay [3]and Morishita et al.[8]methods-respectively given by Eq.(1)and Eq.(2)-overlap exactly and are both located in the purple region.

As can be seen from the comparison in Fig.10, the test data points are all basically in the purple region.For the same explosive charge mass,the BPL of a concrete slab with rebar was significantly smaller than that of a plain concrete slab because the restraint provided by the rebar and the mechanical interlock of the concrete fragments can enhance perforation resistance [30], as shown in Fig.11.It can also be seen in Fig.10 that the simulated data points are located above the upper boundary of the purple region,whereas the plain concrete test data points are located inside the purple region because the UFC equation is only applicable to reinforced concrete slabs.In addition,the BPL of the plain concrete slab determined by the experimental results is slightly smaller than that determined by the FEM; this may be because the full extent of actual damage indicated by the finite element simulation was not considered in the experimental measurement of the observed damage, which was smaller [30].This suggests that the BPL obtained using the FEM is on the conservative side and thus demonstrates that the calculation results are reliable.

4.2.L/D > 1

Fig.12 shows the slab damage induced by the detonation of an explosive charge with parameter set no.5,as defined in Table 2,at criticality(R=Rc)using three different DOEs.It can be seen in the figure thatRcincreased with increasing DOE.

Fig.11. Section plane of the concrete slab at the BPL in Yamaguchi et al.[34](m = 200 g, L = 58 mm, D = 58 mm, DOE = 0 mm, Rc = 100 mm).

Fig.12. Numerical simulation results with m=100 g,L/D=6.87,and(a)DOE=0 mm,Rc = 70 mm; (b) DOE = 250 mm, Rc = 180 mm; (c) DOE = 450 mm, Rc = 200 mm.

Fig.13. Damage on the back of the concrete slab at Rc = 250 mm in Duan et al.[20].(m = 110 g, L = 93.5 mm, D = 30 mm, DOE = 450 mm, Rc = 250 mm).

Duan et al.[20]evaluated the damage owing to the detonation of a cylindrical TNT charge embedded inside a plain concrete slab at different DOEs to investigate the BPL of the concrete slab according to the corresponding explosive charge mass.In the ideal case exactly at the BPL, a very small hole would appear in the center of the slab, as shown in Fig.13.The experimental method employed by Duan et al.to find the mass of the explosive charge corresponding to the BPL adjusted the length of the explosiveLfor the sameH, DOE,R, andD.This differs from the numerical simulation method used in this study,which found theRccorresponding to the BPL by adjusting the residual thickness of the concrete slabRfor the sameH,DOE,L,andD.Briefly,the explosive charges used by Duan et al.hadDvalues of 25-40 mm, L values of 53.5-270.3 mm, and DOEs of 300-520 mm;they were detonated from the end nearest the top surface of the slab.The slab H values employed were 700-900 mm, and the slab width was 1000-1500 mm.Notably,adjacent concrete slabs for subsequent tests were spaced apart using 20 mm thick wood plates along their sides.

To evaluate the reliability of the numerical simulations whenL/D>1, their results were compared with the test data published in Duan et al.[20](referred to as the“test results”in this section),as shown in Fig.14.It can be observed in the figure that the maximum L/D of the explosive charges was 10.4, suggesting that L/D exerts considerable influence on the BPL.When the explosive charge diameters D were the same, the larger the L/D, the larger the Rc; for the same L/D, the larger the D, the larger the Rc.

The hollow points in Fig.14 represent the numerical simulation results obtained in this study.When the DOE andL/Dremained the same, the larger the DOE, the larger theRc.As theRcgradually converged to a constant value, the data points with the largest vertical coordinates corresponding to the sameL/Doverlap.WhenD=34 mm,the numerical simulation results show that the larger theL/D,the larger theRc,which is the same pattern observed in the test results.

The simulations conducted in this study considered cases with DOE values ranging from 0 to ∞, while only the cases with a large DOE (DOE >L) were considered by Duan et al.Therefore, the data obtained when DOE ＜l should be excluded when comparing the simulated and test results.It can be observed in Fig.14 that both the simulated and test results show the same change trend:Rc/Dincreases with increasingL/D.However,for the sameL/Dand similarDvalues,theRc/Dvalues obtained in the tests were generally larger than those obtained using the numerical simulation.There are three main reasons for this result:1)the explosion was initiated at one end in the tests,so the strength of the stress wave transmitted from the bottom of the explosive to the concrete medium was greater [39,40]; 2) the test slabs were spaced apart using wood plates along their sides,so the stress waves reaching the sides of the test slabs were reflected as tensile waves,increasing the damage to the slabs; 3) the strength of the concrete used in the tests(fc= 35 MPa) was lower than that used in the numerical simulations.

So, I was very annoyed7 by those things. And I even told Ivy to tell Mamun to stop these foolishness. After my exams were over, I had a break. So I used to go to the roof and read books to spend my time. Mamun used to come to their roof also and both roofs where so close to each other that you can just jump from one to another.

Fig.14. Comparison of test results (from Duan et al.[20]) with numerical simulation results from this study.

5.Analysis of parameter influence and development of empirical model

Based on the analysis presented in Section 4,it can be concluded that the DOE and shape of the explosive charge exert considerable influence on the damage to a protective concrete shelter.This section therefore presents a systematic analysis of the influence of DOE andL/Dto provide guidance for the design of concrete shelter thicknesses.

5.1.Influence of DOE

According to Fig.15, as the DOE increases,Rcfirst increases rapidly before slowing and finally reaching a constant maximum value.Under the samemand DOE, the larger theL/D, the smaller theRcbecause a largerL/Ddecreases the impulse from the bottom of the explosive acting on the concrete medium[7].If DOE=0,its influence can obviously be ignored, allowing the influence of the shape of the explosive to be more clearly observed(see Fig.16).In this case, the larger theL/D, the smaller theRc.

To better describe and analyze the influence of the DOE, τ was defined in this study as the DOE coefficient (see Fig.17).For the simulation results obtained in this study, τ was calculated as the ratio of the Rcunder any given DOE to the Rcwhen DOE = 0, i.e.,τ = Rc/Rc0.For tests in which ground shock was measured, τ was calculated as the ratio of the ground shock magnitude under any given DOE to the ground shock magnitude when DOE=0.Note that the ground shock referred to here includes the peak pressure and impulse of the stress wave, the peak (particle) velocity, the peak(particle) acceleration, and the peak(particle) displacement [15].

Fig.18 compares the τ calculated using the numerical simulation results obtained in this study with the τ calculated using experimental results presented in the literature.Notably,compared with the use of DOE/,using DOE/L or DOE/D as the abscissa in Fig.18 not only describes the influence of DOE on τ,but also the influences of different explosive shapes on τ.Moreover, the ability to refer to DOE/L = 1 helps to distinguish partially embedded and fully embedded cases.As a result,DOE/L was selected as the abscissa for this comparison.

Fig.15. Rc values for different DOEs from the numerical simulation results.

First,the simulation results were analyzed.According to Fig.18,τ increased with increasing DOE/L,such that the larger the L/D,the larger the τ.When L/D = 6.87,τ reached a maximum value of 2.9.When DOE/L increased from 0 to 1,τ increased rapidly because the length of explosive charge exposed to the air decreased, causing more and more explosive energy to be transmitted into the concrete.For DOE/L > 1, the increase in τ slowed, and once DOE/L reached a certain value, τ reached its maximum and remained constant thereafter.It can also be seen in Fig.18 that the larger the L/D,the smaller the DOE/L corresponding to the maximum τ.When L/D = 0.89, i.e., when the shape of the explosive was cubic, the numerical simulation results for two models with different masses were basically the same, with τ reaching a maximum of ~1.6.

Next,the test results were analyzed.Mu et al.[16],TM5[15],and Haas and Rinehart [41]all studied the influence of the explosive charge DOE on the ground shock magnitude using blast testing.The ground shock magnitude is directly related to the blast damage to the concrete slab,so it has a considerable impact on the BPL.Thus,ground shock magnitude offers a meaningful basis for the comparison of previous test results with the numerical simulation results.Both Mu et al.[16]and TM5 [15]used spherical explosives embedded in concrete materials for blast testing, with the detonation initiation in the center of the charge.When processing the data, the diameter of the spherical explosive was taken as theLof an equivalent cylindrical explosive.Haas and Rinehart [41]used a cylindrical explosive with anL/Dof 3 and an initiation point at far end from the surface of the test slab,which was a marble specimen.

As can be seen in Fig.18, the results from Mu et al.[16]were obviously larger than those from TM5[15].When DOE/L=12,the τ obtained by Mu et al.[16]was about to reach its maximum value of 4.7,while that from TM5[15]was only 1.4.Mu et al.explained that there were two main reasons for this substantial difference between their test results and those of TM5.The first reason is that there can be different definitions of a closed explosion.The second reason is that the upper part of the concrete in the TM5 tests had a high-strength concrete bullet-shielding layer.However, τ is based on the magnitude of the ground shock at the time of contact explosion, which is independent of the “definition of a closed explosion”; moreover, the use of a “high-strength concrete bulletshielding layer” should actually cause an increase in the value of τ.Therefore, Mu et al.’s reasoning cannot be used to explain the observed difference between the τ values obtained from the studies compared in this paper.In addition, Mu et al.[16]measured the peak stress wave pressure and the peak particle acceleration at different scaled distances from the center of the explosive ranging from 0.255 to 2.504 m/kg1/3.The τ calculated by directly using the raw data measured from their tests was highly scattered,and so the arithmetic mean of τ was used in Fig.18.In fact,the τ values ranged from 1.34 to 8.48, with the lowest overall τ for the peak stress occurring at a distance of 0.515 m/kg1/3.The τ values for DOE/L=0.5 and DOE/L=12 were 1.34 and 1.99,respectively.However,when considering the influence of DOE, TM5 accounted for the peak pressure and impulse of the stress wave, the peak (particle)velocity, the peak (particle) acceleration, and the peak (particle)displacement.Therefore,the conclusions of TM5[15]regarding the DOE need to be analyzed in further detail using specific experimental data and compared with the conclusions provided by Mu et al.[16].Furthermore,as Haas and Rinehart[41]only considered the difference in the impulse of the stress waves, their experimental results were similar to those of Mu et al.and differed from those of TM5 as well.

Fig.16. BPL failure modes of a concrete slab (DOE = 0 m; m = 100 g) with (a) L/D = 0.89, (b) L/D = 2.07, (c) L/D = 6.87.

Fig.17. Schematic diagram defining the DOE coefficient τ.

Fig.18. Influence of DOE on τ.

Finally, the numerical simulation results were compared with the experimental results.Since spherical and cubic bodies are similar in shape, it can be seen in Fig.18 that the numerical simulation results using a cubic explosive were very similar to the results provided by TM5.However, the results in TM5 can only be used as a reference because they describe the influence of DOE on ground shock, while the numerical simulation results directly quantify the influence of DOE on BPL.In addition, the numerical simulation results show that the DOE had an enhancing effect on the BPL after DOE/L> 5.In TM5 [15], when DOE/L= 1.5, the DOE ceased to exert any influence after reaching its maximum value;this was a result of the difference in the “definition of a closed explosion.” Indeed, Mu et al.suggested that a closed explosion is achieved when DOE/L = 12, which is similar to the numerical simulation results obtained in this study.

To quantify the influence of DOE on BPL, the numerical simulation results forL/D=0.89 were selected as the reference standard,and an empirical equation set was proposed using linear curve fitting segments as follows:

As shown in Fig.18, the trend of change in DOE/Lcalculated using Eq.(4)was basically the same as in the numerical simulation results.In addition, the results calculated using the empirical equation were greater than or equal to the numerical simulation results.Therefore, Eq.(4) can be used to quickly calculate the BPL under different DOEs when the influence of explosive charge shape is not considered.

5.2.Influence of L/D

According to the analysis in Section 5.1, both the DOE andL/Dexert considerable influence on the BPL, making it difficult to separately quantify the influence of either.Thus,to better describe the influence ofL/D, μ is defined in this study as the shape coefficient(see Fig.19).For the same DOE and explosive charge mass,μ is the ratio of the Rcvalue corresponding to any given L/D to the Rcvalue when L/D = 0.89.

Fig.20 describes the influence of explosive shape(L/D)on μ.The abscissa in Fig.20 is given as DOE/to eliminate the effects of the different explosive charge masses used in the different FEMs.This abscissa is different from that employed in Fig.18 because it does not need to reflect the influence of charge shape.It can be observed in Fig.20 that the value of μ was constant at 1 when L/D=0.89.For the same L/D(L/D>1),μ gradually increased with increasing DOE,then remained unchanged.For the same DOE, μ gradually decreased as L/D increased.To facilitate engineering applications,the calculated results for μ were fitted to obtain:

Fig.19. Schematic diagram defining the shape coefficient μ.

Fig.20. Influence of L/D on μ.

Fig.21. Comparison of numerical simulation and equation-fitted results for μ.

for 0.89 ≤L/D≤6.87.Fig.21 shows a comparison of the results obtained using Eq.(5) and the numerical simulation results,showing a good fit with differences of less than 5% in most cases.

5.3.Empirical equation for BPL

Based on the analysis presented in Sections 5.1 and 5.2, the relationship between BPL0and BPL can be obtained as follows:

Fig.22. Schematic of application of the empirical equation.

in which μ and τ＇describe the influence of the shape and DOE of the explosive, respectively, while BPL0reflects the influence of the concrete material characteristics.

According to Fig.22, if the properties of a particular concrete slab are known(e.g.reinforcement content,concrete strength,etc.),Eq.(6) can be used to predict the BPL for a specified DOE and explosive charge shape as long the BPL0of the same concrete slab when subjected to a contact explosion of a cubic charge is known.The value of BPL0can be obtained using experiments or estimated by various prediction methods, such as that presented by UFC [5].Alternatively, if the DOE and shape of the explosive charge are known,BPL0can also be extrapolated from Eq.(6),in turn allowing for the calculation of the BPL for the same explosive charge mass at other DOEs and in other shapes.To further extend the application of Eq.(6), if the damage effects (perforation or non-perforation) are known for a certain concrete slab, the damage induced by an explosive charge with the same mass at other DOEs and in other shapes can also be inferred.

When using Eq.(6),two aspects require particular attention.The first aspect is the influence of the slab boundary condition.It was determined in Section 3.4 that an unrestrained boundary can significantly enhance the damage induced by the detonation of an explosive charge embedded in a concrete slab.However,in practice,the influence of the boundary condition can be ignored if the width of the concrete slab is sufficiently large,such as found in large-scale projects.However,for other tests in which the size of the concrete slab is limited, the influence of the boundary condition cannot be ignored.The literature provides two common methods for reducing the influence of the boundary condition:one separates rectangular slabs with wood panels[20,42],and the other wraps the outside of a cylindrical slab with a steel plate[17,18].However,the influence of these methods on the damage induced in concrete slabs by the detonation of embedded explosive charges still requires further investigation.

The second aspect that requires attention is the influence of the explosive charge shape.Common shapes used for testing include rectangular prisms,cylinders,spheres,and hemispheres.In general,when the bottom surface of a rectangular prism is square and its area is the same as that of a cylinder, the strength of the stress waves transmitted into the concrete through the base of the prism can be considered to be equal to those transmitted through the base of a cylinder of the same height.However, if the bottom of the rectangular prism is not strictly square, its influence can become more complicated; this condition is beyond the scope of the present study.In addition, note that the range ofL/Dfor Eq.(6) is 0.89 ≤L/D≤6.87.A sphere is close to a cylinder withL/D=1,so Eq.(6) is also applicable to a spherical charge.While a hemisphere is close to the shape of a cylinder withL/D= 0.5, further analysis remains required to determine whether the proposed equation is suitable for hemispheric charges.

6.Conclusions

To investigate the influence of explosive charge shape and DOE on the BPL of a concrete shelter, more than 300 numerical simulations were conducted in this study using an FEM with strictly controlled variables.The following conclusions were drawn from the results.

(1) When blast testing using embedded explosives, the stress waves in the concrete arrived at the lateral free boundary of the slab and were reflected as tensile waves, significantly increasing the diameter of the ejection and spalling craters.Therefore, the width (or diameter) of the concrete slab and its boundary conditions are critical considerations.

(2) For the same explosive charge mass andL/D,the greater the DOE,the larger theRc(whereRc=BPL-DOE).However,Rceventually converged to a constant value.For the same explosive charge mass and DOE, the larger the L/D, the smaller the Rc.

(3) To separately quantify the influence of the explosive charge shape and DOE,this study defined two coefficients,the DOE coefficient τ＇and the shape coefficient μ.A method for calculating the BPL of a concrete shelter was then proposed using these coefficients, and its scope and constraints in application were discussed.It is expected that the results of this study will provide guidance for the improved design of protective concrete structures.

As stated in the introduction, the combined effects of charge penetration and explosion damage should be considered in the design of concrete shelter envelope thickness.Note that penetration will produce pre-damage that helps to dissipate explosion energy into the air, thereby reducing that transmitted into the concrete and resulting in less damage to the concrete slab owing to the explosion.However, when an explosive charge is directly placed in a prefabricated hole, the embedment effect will be quite significant as no penetration damage is inflicted on the surrounding concrete first,increasing the damage to the concrete slab owing to the explosion.Therefore, the use of prefabricated holes in the explosion test may result in more serious damage to the concrete than performing the penetration test first, then the explosion test; that is,the empirical model proposed in this paper should be considered conservative.Understanding the full effects of penetration on the BPL require further experimental and numerical analysis.

In future work, the authors intend to perform additional numerical analyses using axisymmetric models to analyze damage to plain concrete slabs so that different explosive shapes can be considered.Furthermore,3D models must be developed to analyze concrete slabs reinforced by steel bars or steel fibers.Finally, concrete slabs subjected to partially embedded explosions should be physically bisected to allow for the observation of richer damage details in the section and identify the scope of real damage,further validating the FEM employed in this study.

Funding

This work was supported by the National Natural Science Foundation of China (Grant No.51978166).

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgments

We would like to thank Editage [www.editage.cn]for English language editing.