Finite Element Simulation of Medium-Range Blast Loading Using LS-DYNA

This study investigated the Finite Element simulation of blast loading using LS-DYNA. The objective is to identify approaches to reduce the requirement of computation effort while maintaining reasonable accuracy, focusing on blast loading scheme, element size, and its relationship with scale of explosion. The study made use of the recently developed blast loading scheme in LS-DYNA, which removes the necessity to model the explosive in the numerical models but still maintains the advantages of nonlinear fluid-structure interaction. It was found that the blast loading technique could significantly reduce the computation effort. It was also found that the initial density of air in the numerical model could be purposely increased to partially compensate the error induced by the use of relatively large air elements. Using the numerical approach, free air blast above a scaled distance of 0.4m/kg was properly simulated, and the fluid-structure interaction at the same location could be properly duplicated using proper Arbitrary Lagrangian Eulerian (ALE) coupling scheme. The study also showed that centrifuge technique, which has been successfully employed in model tests to investigate the blast effects, may be used when simulating the effect of mediumto largescale explosion at small scaled distance.


Introduction
The effect of blast loading has become a critical issue in the design, protection, and rehabilitation of engineering structures following the terrorist attacks on civil infrastructures in recent years.But analyzing the effects of blast loading on civil infrastructures is a difficult task, as it involves highly nonlinear fluid dynamics, structural dynamics, and fluid-structure interaction.At present, most blast resistant analyses make use of simplistic blast loading and structure models [1], but it is difficult to assume rational blast loading for a complicated structure that includes sophisticated surfaces and corners.Advanced numerical simulations thus become the economic and reliable approaches for this purpose.In this type of numerical simulations, generally all the critical components in the problem are properly modeled, including explosive, air, fluid-structure interaction, and structures.Many such numerical simulations can be found in the literature [2][3][4][5][6][7][8][9], and several computer codes have been designed to carry out such analyses [10,11].
What remains a critical issue is the computation effort required to accomplish such an analysis.For example, simulating the effect of blast loading on a full-scale bridge could require more than 10 million Finite Elements.The bottleneck in this type of computer simulation mainly comes from the requirement of element size for the explosive and air domains.Some studies showed that, in order to capture the accurate release of blast energy from explosives, element size as small as 0.5-1.0mm should be used [7,12].As for the blast wave propagation in the air domain, generally element size in the order of few mm was used to obtain reasonable results when the scaled distance is smaller than 1 m/kg 1/3 [2,13,14].
This study investigated the Finite Element simulation of blast loading using LS-DYNA, focusing on the blast effect on a structure surface at a medium scaled distance ( = 0.4 m/kg 1/3 to 1.0 m/kg 1/3 ).LS-DYNA has been commonly 2 Shock and Vibration used to analyze the effects of blast loading on structures [5,6,8,16], and this study hopes to shed some light on the application.
The presentations of this study are organized as follows: after the discussion of the theoretical background, the approaches for blast loading simulation in LS-DYNA are introduced.The findings from the simulation of free air blast are then presented, and they are employed to simulate the reflected blast pressure on a rigid surface.Finally, the findings from this study are validated against a model test.

Theoretical Background
2.1.Blast Wave Propagation in the Air.Blast wave propagation in the air is governed by three fundamental equations of fluid dynamics: mass conservation, momentum conservation, and energy conservation equations.These three equations can be combined to yield the Navier-Stokes equation that governs the blast wave propagation in the air [17].In 1958, von Neumann [18] gave the solution of air blast from a point source explosion in an infinite three-dimensional space.He assumed that the air was ideal gas, the equation of state (EOS) of which was  = , and the total energy for the air was  = [/(−1)], with  being the specific gas constant defined as  =   −  V and  being the heat capacity ratio defined as  =   / V .Here   is the specific heat under constant pressure and  V is the one under constant volume.He showed that the free air blast pressure could be expressed as 0  (, , ) . ( Here  0 is the initial density of the air before blast,  is time, and  is a parameter representing the location of the point of interest in the free air domain.(, , ) is a complicated function of , , and , the detail of which can be found in von Neumann [18].This solution is valid for 1 <  < 2. Equation (1) showed that, by increasing the initial air density  0 , the blast pressure in the air due to explosion would increase.This theoretical fact can be employed to increase the blast pressure in numerical simulation and to partially compensate the loss of energy, or numerical damping, due to large element size.

Blast Wave Propagation and Finite Element Size.
Extensive experiences of numerical simulations have shown that large element size may lead to smaller free air blast pressure than the theoretical and experimental values.Apparently, in blast simulation, significant "numerical damping" exists due to relatively large element size.There exist a few studies investigating the required element size in the simulations of blast waves and their effects on structures.Chapman et al. [2] found that a grid size of less than 3 mm was necessary to simulate the blast pressure and specific impulse at a scaled distance of 0.95 m/kg 1/3 .The amount of explosive simulated was 0.075 kg TNT.Luccioni et al. [19] compared their simulated peak pressures with those presented by Kinney and Graham [20] and concluded that a grid size of 100 mm was adequate to simulate the peak pressure due to the explosion of 100 kg TNT at a scaled distance of larger than 1.2 m/kg 1/3 .The parametric study by Shi et al. [13] showed that a grid size of 100 mm was able to capture the peak pressure and specific impulse at a scaled distance of 2.0 m/kg 1/3 due to an explosion of 1000 kg TNT.Recently, Wu et al. [12] used an element size of 10 mm to model the air blast effect at a scaled distance of 2 m/kg 1/3 from an explosion of 0.2 kg TNT.These analyses were all carried out using the Finite Element code AUTODYN.It can be seen that large air element size may be employed if the targeted structure is at large scaled distance.
The relationship between element size and amount of explosive is related to the blast wavelength   , which is defined as the length of positive pulse in the air due to an explosion according to UFC 3-340-02 [1].There exists a unique relationship between scaled wavelength (unit:   / 1/3 ) and scaled distance of incident blast due to a free air explosion [1].This means that, for the same targeted scaled distance, larger amount of explosive resulted in larger wavelength, which may in turn relax the requirement on element size [21].This can also be seen from the results of available studies [2,13,19].However, it is not known whether the same scale ratio can be used to determine the element size for large-scale explosion simulation based on that of a small-scale explosion simulation.For example, if a grid size of 3 mm is adequate to simulate the blast effect on a structure due to 0.1 kg TNT, it is not known whether, for the same targeted scaled distance, a grid size of 30 mm is enough for an explosion of 100 kg TNT.

Simulation of Blast Loading in LS-DYNA
Conventionally, blast loading may be simulated by employing the Multimaterial Arbitrary Lagrangian Eulerian (MM ALE) solver in LS-DYNA.In this approach, Lagrangian formulation is used to model solid continua and structures, and Eulerian formulation is used for large distortions of fluids, gases, and explosives.The two domains are coupled through appropriate schemes.The explosive must be properly modeled using very small Finite Element size [7] and a large number of fluid elements must be simulated between the explosive and structures to transfer the blast wave.It is therefore extremely demanding to complete a full-scale simulation involving large structures.
There also exists a simplified approach, which imposes pressure time-history on structure surfaces based on the CONWEP reflected pressure on a rigid surface [1].A disadvantage of this approach is that wave reflection and superposition cannot be accounted for as would occur at the corners of a structure.The reflected pressure time-history on a deformable surface is also not the same as that on a rigid surface.Børvik et al. [22] demonstrated that this approach resulted in large errors in the simulation of a simple container structure subjected to external blast loading.
The combination of the CONWEP incident pressure and the MM-ALE solver can be used to obtain acceptable results with less computational efforts [14].This coupling method is explained in Figure 1.Only the air immediately surrounding the target structure is modeled with an ALE domain.Incident blast pressure time-history based on CONWEP is applied to a layer of ALE elements (the ambient layer shown in Figure 1), which faces the explosive charge and acts as a source for the adjoining air elements.The blast wave then propagates through the air domain and eventually interacts with the structure.The loading on the structure due to wave reflection and superposition can then be captured by the ALE domain.This capacity can be realized in LS-DYNA by employing the key word LOAD BLAST ENHANCED (LBE) and defining a layer of special elements (the ambient layer).Figure 2 compares the blast pressure inside the ambient layer with that of CONWEP in a blast simulation [14].The sizes of the air elements and the ambient layer elements were all 3 mm.Almost identical results were obtained.Different thicknesses of the ambient layer were also tested, which resulted in identical air pressure in the layer.It is noted that similar approaches have been employed for blast simulations in previous studies [22].Similar to the LOAD BLAST ENHANCED in LS-DYNA, the "inflow" of blast pressure was realized by defining a special layer of Eulerian Finite Elements.

Simulation of Free Air Blast
In this section, explosions of different amounts of spherical TNT in an infinite space filled with ideal gas were simulated using the new blast loading scheme in LS-DYNA.The results were compared with those provided by UFC 3-340-02 [1].In all of the analysis that employed the LBE approach, the ambient layer is 4 mm thick and stands between the block of air elements and the point source of explosion.
The air in all the numerical simulations was assumed to be ideal gas and modeled using the linear solid elements in LS-DYNA.The elements are 8-noded linear hexahedron elements integrated using 8-point Gaussian method.The equation of state for the air is given as Here  is the ratio of specific heat and was assumed to be 1.4 according to previous numerical experiences [14]. 0 is the initial density of air.The initial value of  0 was given as 0.25 MPa according to Schwer [14] and assuming an initial air temperature of 20 ∘ C. It can be seen that the EOS in ( 2) is basically identical to the one employed in von Neumann [18].It has also been employed by others to model air in their blast simulations [5,6,8,16].
Null material has no shear stiffness and hourglass control must be used with great care.In some applications, the default hourglass coefficient in LS-DYNA might lead to significant energy losses.The hourglass control scheme that was used by Schwer [14] to successfully model ideal gas was employed in this study.

Effects of Initial Air Density.
The theoretical analysis by von Neumann [18] shows that when the air is assumed to be ideal gas, the blast pressure increased with an increase in the initial air density.In other words, if air density larger than the actual value is used in a blast analysis, it would overpredict the air blast pressure if everything else is accurate.This overprediction may be employed to compensate the numerical damping resulting from relatively large element size.In this set of parametric studies, the length, width, and height of the hexahedron elements were all 4 mm, except in one case where the grid sizes were 3 mm.The air domain was simulated as a rectangular block of air elements, the dimension of which was 14.6 cm × 30 cm × 55 cm, as shown in Figure 3(a).A spherical TNT of 8 kg was assumed to be located 30.2 cm (scaled distance = 0.151 m/kg 1/3 ) away from the ambient layer.Three initial air densities were tested: the original air density  0 = 0.00129 g/cm 3 , 1.2 of the original value, and 1.4 of the original value.
Figures 3(b) and 3(c) show the relationships between scaled distance and incident peak blast pressures and positive specific impulses.It can be seen that, with an increase in the initial air density, both the simulated peak blast pressure and specific impulse increased.Compared with the values in UFC 3-340-02 (2008), an air density factor of 1.2 yielded the best results with a grid size of 4 mm.A proper increase in the initial air density had a similar effect to a reduction in element size.The numerical simulation was not as good at small scaled distance (<0.2 m/kg 1/3 ).Similar results were also reported by other investigators in their simulation of free air blast [13,19].
It is noted that, with 8 kg of TNT, the numerical simulations with a density factor of 1.2 and a grid size of 4 mm were adequately accurate for the range of scaled distance around 0.4 m/kg 1/3 .This is very relevant for the simulation of blast effects on structures due to terrorist attacks.
Other parameters in the EOS for air in (2) were also varied to investigate their effects on the blast pressure and specific impulse, but the results are not presented here due to space limitation.It was shown that although both the specific ratio of heat  and the initial air temperature  may be varied to obtain reasonable peak blast pressure with relatively coarse mesh, their effects on the specific impulse were not satisfactory.

Amount of Explosive and Element Size.
As previously discussed, the required element size to capture blast wave propagation is related to the amount of explosive, as a unique relationship between scaled distance and scaled wavelength exists.This fact is demonstrated in Figure 4, which shows the numerical simulations of three free air blasts.The amounts of TNT employed were 0.025 g, 1.6 kg, and 8.0 kg, respectively; the element size was still 4 mm; and the air density factor was 1.2.Apparently, the case with 0.025 kg TNT was the worst, while the case with 8.0 kg of TNT yielded the best results in terms of incident peak pressure and specific impulse.The accuracy with 1.6 kg of explosive was also acceptable.The results show that smaller element size is needed while simulating small-scale explosion.
However, proportional scaling of element size based on the wavelength did not seem to work in the numerical analyses.Figure 5 compared the numerical simulations of three cases: 1.2 kg TNT with 4-mm element size, 32.4 kg TNT with 12 mm element size, and 200 kg TNT with 22 mm element size.The scaled element size by the weight of explosive was basically the same in these three cases.However, the 200 kg case resulted in much different incident peak pressures and specific impulses compared to the UFC 3-340-02 values, even with increased initial air density.The results for the case 32.4 kg TNT were also not satisfactory.This is because the mechanisms of numerical damping and overprediction due to large initial air density are not the same.For a given element size, accurate air blast pressure may not be obtained by increasing the initial air density.The element size has to be reduced with a proper increase of the initial air density to result in accurate blast pressure.
The results in Figure 5 showed that, for relatively small standoff distance that corresponds to a scaled distance smaller than 1.0 m/kg 1/3 , the numerical simulation would be very demanding if relatively large amount of explosive is used.For example, for an explosion of 200 kg TNT and a standoff distance of 3 m, in order to properly simulate the specific impulse, a very large number of air elements have to be used.The computation effort would be very large even if the LBE blast loading scheme is used.

Reflected Blast Pressure on a Rigid Surface
The geometrical model of this series of simulations was the same as the one in Section 4.1.However, a rigid boundary was placed inside the air domain.The rigid boundary was modeled with one layer of extremely stiff solid element which overlapped with the air domain and was fixed on one side.Five different amounts of TNT, 0.088 kg, 0.176 kg, 0.264 kg, 0.353 kg, and 8 kg, were tested, the scaled distances of which to the rigid boundary ranged from 0.35 m/kg 1/3 to 0.68 m/kg 1/3 .The air density scaled factor was also 1.2.The simulations of free air blast with the same amounts of TNT at the same scaled distances were also carried out, which resulted in very close peak pressures and specific impulses compared to the UFC 3-340-02 values.These results are not presented herein to save space.
Figure 6 shows the results and their comparisons with the corresponding values in UFC 3-340-02 [1].The incident blast wave was perpendicular to the rigid surface.Except for the reflected specific impulse at 0.35 m/kg 1/3 scaled distance, the other simulation results were very close to those in UFC 3-340-02 [1].The reflected specific impulse at 0.35 m/kg 1/3 scaled distance was about 20% smaller.The smaller reflected specific impulse at this location could be due to the small scaled distance [13,19].

Simulation of Blast Effect on a Steel Plate
Boyd [15] carried out two identical experiments (E14 and E15) to investigate the effect of air blast on a steel plate.Figures 7(a) and 7(b) show the experiment setup.A 5 mm thick steel plate was supported by four concrete blocks.250 g of Pentolite explosive was denoted 50 cm above the steel plate.The steel plate was attached to a steel ring and placed directly on the concrete blocks.The steel plate was instrumented for vertical displacement at the center, accelerations at two locations, and overpressures at another two locations, as shown in Figure 7(b).
Due to symmetry, only half of the system was modeled in LS-DYNA as shown in Figure 7(c).The air domain above the steel plate was modeled by two approaches: one with 10 mm element size and an air density factor of 1.0 and the other with 15 mm element size and an air density factor of 1.2.Larger elements could be used herein as the scaled distance of interest was larger than 0.8 m/kg 1/3 , which were selected by comparing the incident air pressure from numerical simulation to that by UFC 3-340-02 [1].The thickness of air domain above the steel plate was 255 mm, and on its top laid Scaled distance (m/kg 1/3 ) Peak pressure (MPa) Scaled distance (m/kg 1/3 ) Incident impulse (MPa * ms) the ambient layer.The steel plate was meshed with Lagrangian solid elements and coupled with the Eulerian elements for the air.In order to capture the bending response, 5 layers of solid elements were used along the thickness of the steel plate.The steel plate was modeled by the isotropic elastoplastic model, with Young's modulus of 203 GPa, Poisson's ratio of 0.3, and a yield strength of 270 MPa according to Boyd [15].The steel ring was simulated by elastic material, with Young's modulus 10 times larger than that of the steel plate.Frictional surfaceto-surface contact was applied between concrete blocks and the steel ring.250 g Pentolite explosive was converted into 260 g equivalent TNT explosive and used in the analysis [15].
The results of both simulations matched well with the experimental data, as shown in Table 1.The peak reflected pressures simulated were smaller than the measured ones.However, Boyd [15] reported that the traces from the pressure   transducers exhibited a large amount of noise at the peak.
The measured overpressure jumped and reduced immediately, causing very small difference in the specific impulse, which indicated that it could be due to some noise in the measurement.

Discussions
A parametric study was carried out in this study to investigate a numerical approach for simulating air blast and fluid-solid interaction using LS-DYNA.It was found in this study that, with the ideal gas EOS in (2), some increase in the initial air density may increase both the peak blast pressure and the specific impulse, so that relatively large air elements may be used to capture the air blast pressure at medium scaled distance ( = 0.4 m/kg 1/3 to 1.0 m/kg 1/3 ).It is not known at this stage whether the same conclusion still applies if more sophisticated EOS for air is used and more study is necessary to obtain the conclusion.As can be seen in the results presented above, above a scaled distance of 0.4 m/kg 1/3 , acceptable accuracy of blast pressure on a structure surface was obtained as long as the free air blast pressure at the same scaled distance was accurate.Therefore, in a numerical modeling employing the approach discussed in this study, the parameters and grid size for air may be determined in advance by free air blast simulation.This study also shows that, for medium-range explosion at a scaled distance smaller than 1.0 m/kg 1/3 , the numerical simulation would be very demanding if the effect of relatively large amount of explosive on a large structure is simulated.One alternative is to use the centrifuge technique that has been successfully employed in experiments to simulate the effects of blast on geotechnical structures [23,24].Using this technique, the gravity acceleration would be increased to  times of the original value, while the model size would be reduced by  times.According to the scaling law of blast effect, the amount of explosive in the model would then be 1/ 3 of the prototype amount.For example, if  = 5, 1.6 kg of TNT could be used to simulate the effect of 200 kg of TNT in the prototype scale, while the model size could be reduced by 5 times.The element size can be as large as 4 mm using an air density factor of 1.2.According to the scaling law of centrifuge technique, the stress and strain in a centrifuge model are the same as those in the prototype structure.Therefore, as long as the structural materials are properly modeled in the numerical simulations, centrifuge technique could result in accurate estimation of structural responses under blast loading.
The study also found that, at small scaled distance (<0.4 m/kg 1/3 ), the new blast loading scheme also resulted in smaller peak pressure and specific impulse, even if they were accurate in the boundary layer of inflow pressure (the ambient layer).This is consistent with previous simulations that modeled the explosive and the surrounding air [13,19].However, compared to the conventional approach, the new blast loading scheme has the potential to duplicate the peak blast pressure and specific impulse above a scaled distance around 0.4 m/kg 1/3 .This is very relevant for practical applications.
Finally, it should be pointed out that the new blast loading scheme cannot simulate the effects of fragments from some sources of explosions.

Conclusions
In this study, a recently developed scheme for blast loading in LS-DYNA was used to investigate the effect of mediumrange blast loading on a structure surface.In this scheme, a special layer of Eulerian elements is employed to receive the inflow of blast pressure from a point source of explosion.The blast pressure time-history at the ambient layer, namely, the pressure inflow, is determined using the CONWEP prediction.The blast pressure then propagates through the adjacent air domain to the targeted structure.
The objective of this study is to validate this blast loading scheme, to investigate the combined effect of element size and initial air density, and to shed light on the simulation of medium-range blast loading on structures using LS-DYNA.The following conclusions were obtained from this study: (1) The new blast loading scheme could be used to properly simulate air blast, provided that proper element size corresponding to a specific amount of explosive is used.This blast loading scheme can significantly reduce the number of air elements and totally eliminate the explosive elements, thus lowering the computational effort.
(2) The initial density of air in the numerical model could be purposely increased to some extent when the equation of state for ideal gas is used to model the air domain.The increase of air density could partially compensate the error induced by the use of relatively large air elements.The exact ratio of density increase can be obtained by comparing the free air blast pressures with the values in UFC 3-340-02 [1].
(3) Although relatively larger element size could be used in the modeling of large air blast, the required element size was not proportional to the cubic-root of the weight of explosive , although there exists unique wavelength   / 1/3 as a function of scaled distance according to UFC 3-340-02 [1].Centrifuge technique, which has been successfully employed in model tests to investigate the blast effects, may be used when simulating medium-to large-scale explosion.
(4) Above a scaled distance around 0.4 m/kg 1/3 , free air blast can be properly simulated by the proposed numerical approach, and the fluid-structure interaction at the same location could be properly duplicated using proper Arbitrary Lagrangian Eulerian (ALE) coupling scheme.

Figure 1 :Figure 2 :
Figure 1: Illustration of the LOAD BLAST ENHANCED approach in LS-DYNA.

6 UFC 3 - 3 )Figure 4 :
Figure 4: Requirement of amount of explosive on element size: (a) peak incident pressure and (b) specific impulse.

Figure 5 :
Figure 5: Relationship between amount of explosive, element size, and initial air density: (a) peak pressure and (b) specific impulse.