Nonlinear dynamic response and deformation analysis of soil under the explosion shock loading

2.1. Experimental equipment and approach

This study used the shock wave pressure produced by and ground acceleration induced by an explosion to verify the accuracy of the numerical analysis model. The goal of the explosion experiment was to measure the free-field air explosion pressure that was in contact with the explosion and the ground acceleration to verify the numerical analysis model. In the onsite experiment, this study used 113.389 g (0.25 lb) of trinitrotoluene (TNT) in a cuboid form standing vertically and making contact with the ground. Soil was sampled from the explosion field. Its water content was 10.27 %, and the soil density was 1.78 g/cm³. According to the Unified Soil Classification System, the soil was classified as SC, whereas according to engineering classifications, it was sand with the characteristics of clay. The experimental analysis of direct shear revealed that the shear force parameter cohesion ($c$) and angle of friction ($\varphi$) were 0.95 kg/cm² and 36°, respectively.

Fig. 1 presents the overall configuration of the experimental explosion site and measurement device. To measure the near-ground free-field shock wave pressure during the experiment, blast meters were placed 200 and 300 cm from the explosion center.

Fig. 1Overall configuration of the blast pressure and ground acceleration experimental field and signal acquisition system

The height of the sensor was the same as that of the TNT, approximately 9 cm from the ground. The pressure sensor was a blast pressure pencil probe (model 137A23; PCB Piezotronics, NY, USA). This blast meter can measure the blast pressure in a free field. It has a pen-like metal stick shape, and the sensor is located approximately one-third of the distance from the front of the metal stick. The maximum measurement of blast pressure was 6,895 kPa. To measure the ground acceleration triggered by the explosion, a vertical axial accelerometer was placed 200 cm above the explosion center. The transmission of signals generated by the blast meter and accelerometer was performed by connecting these devices to an oscilloscope through a signal conditioner. Other equipment included power supply and data acquisition systems. The transmission of signals quantity examines blast pressure and ground acceleration via the oscilloscope through a signal conditioner at the same time.

2.2. Establishing the finite element model

According to the research question, analysis type, and algorithm, this study used the solid elements of an eight-node hexahedron together with the multimaterial arbitrary Lagrangian Eulerian (MMALE) algorithm to establish a 3D spatial solid structure model for describing fluid-solid interaction. Each element in the model is defined as having eight nodes. The range of values for the node locations is represented by ±1. Each node has nine degrees of freedom. Their volumes under compressive stress or large deformation are larger than 0, which is suitable for analyzing explicit dynamic temporal and spatial problems [20].

The Lagrangian algorithm is based on material coordination. By contrast, the Eulerian algorithm is based on spatial coordination. In the MMALE algorithm, the Lagrangian algorithm is adopted to analyze solid materials, whereas the Eulerian algorithm is employed for fluid materials. This algorithm overcomes the problem of computation suspension in numerical analysis due to large mesh deformation. Furthermore, it controls and tracks the movement behavior of the boundary structure and is thus suitable for the instant dynamic analysis of fluid–solid interaction. Moreover, it satisfactorily analyzes large displacement and deformation from explosions, blasts, and high-speed collisions. To analyze the dynamic reaction of fluid–solid interaction to shock waves in the soil, the Lagrangian algorithm was used to analyze the solid element, whereas the Eulerian algorithm was used to analyze the liquid element [21].

Fig. 2 displays the one-quarter initial symmetrical model used in numerical simulation analysis for ground contact explosions. The numerical analysis model was established with a unit system in g, cm, μ-second, MPa. For the instant dynamic analysis of fluid–solid interaction, the Eulerian algorithm was used to analyze the elements of TNT and air, and the Lagrangian algorithm was used to analyze soil. In this study, 3D SOLID164 elements were used to construct the numerical analysis model. The meshes of the fluid and the solid were independent, and the density of the solid mesh was twice that of the fluid. The interaction between the fluid mesh and solid mesh was constructed through overlapping [22]. After the CONSTRAINED_LAGRANGE_IN_SOLID program setting was used, the numerical analysis model of the MMALE algorithm was established.

Fig. 2One-quarter initial symmetrical numerical simulation model of the explosion on the ground

The explosion experimental condition was set to no reflex boundary. Because of the symmetricity of the model, only one-quarter of the model was used for analysis. The numerical model boundary condition was set to no reflex boundary for all surfaces apart from the symmetrical level. Air was defined as the ideal gas, and the air material model was 320 × 180 × 420 cm³. The rectangular TNT (3.28 × 3.28 × 9.3 cm³) weighed 113.389 g (0.25 lb) and had a density of 1.63 g/cm³. It was placed at the center of the model in contact with the ground. The point of explosion was at the center. The soil model measured 320 × 180 × 320 cm³. The application of finite element analysis requires consideration of the complexity of an analysis model and of the mesh density. A TNT free-field explosion simulation was analyzed for mesh density convergence under an explosive load. The analysis results revealed the optimal mesh density to be 0.5 times that of the side length of the TNT; the time-step control parameter was set to 0.3 [18]. Per the conversion analysis of the numerical model, the smallest width of the TNT was adopted as a reference for the mesh size, and the Eulerian mesh was set to 0.5 times the smallest width of the TNT. The minimum mesh scales of the rectangular TNT (3.28 × 3.28 × 9.3 cm³), air, and soil were 1.64, 1.64, and 3.28 cm, respectively.

2.3. Material parameters and equation of state

The material constitutive law is expressed through the stress-strain relationship. When a material undergoes large displacement or deformation, the corresponding equation of state is required for analysis to reflect the actual dynamic characteristics of the material. To analyze an explosion, as well as the basic parameters of a material, the corresponding equation of state is required. Table 1 presents the material parameters of air, TNT, and soil. The air model used the software package’s MAT NULL material mode corresponding to the EOS LINEAR POLYNOMIAL equation of state to analyze the free-field shock wave characteristics. The equation is shown in Eq. (1). The parameters of the air model include mass density, pressure cutoff, dynamic viscosity coefficient, initial internal energy, Young’s modulus, and Poisson’s ratio [3, 20]:

with $\mu =\frac{1}{V}-1$, where $E_0$ is the initial energy of a unit volume; $\mu$ is the dynamic viscosity parameter; $V$ is the relative volume; and $C_1$-$C_6$ are constants.

For the TNT model, the Jones-Wilkins-Lee equation of state [Eq. (2)] was adopted from the dynamite manual of the Lawrence Livermore National Laboratory in the United States; this equation corresponds to the MAT HIGH EXPLOSIVE BURN material model in the LS-DYNA equation. It is used to study high speeds, temperatures, and pressures as well as the rapid release of energy during large explosions. For the analysis parameters, this study referred to the LS-DYNA manual and the Lawrence Livermore National Laboratory dynamite manual. The TNT material parameters were mass density, detonation velocity, Chapman-Jouguet pressure, beta burn flag, bulk modulus, shear modulus, and yield stress [17, 23]:

where $p$ is the pressure; $V$ is the relative volume; $E_0$ is the initial internal energy per unit reference specific volume; and $A$, $B$, $R_1$, $R_2$, and $\omega$ are constants representing characteristics of the explosive.

The soil composition model, stress and strain analysis model, and material yielding principle were employed to analyze the characteristics of porous materials and dynamic reactions after a blast. The Krieg yielding principle is based on isotropic plastic theory, which can be used to analyze the crushing or compacting behaviors of porous materials under pressure. The stress-strain relationship in a material before it yields can be analyzed on the basis of linear plasticity. After plastic deformation, the relationship follows the plastic mechanics principle. Damage to a material is controlled by shear force. To analyze the dynamic phenomenon of porous materials, this study divided the behavior of such materials after they yielded into hydrostatic pressure and shear force [24].

Table 1Material parameters and equation of state

Element	Material and equation of state parameters (unit: g, cm, μ-second, MPa)
Air	MAT_NULL
	Density, $\rho$	Pressure cutoff, $P_C$	Young’s modulus, $Y_M$		Poisson’s ratio, $P_R$	Initial internal energy, $E_0$
	0.00129	0.0	0.0		0.0
	EOS_LINEAR_POLYNOMIAL
	Constants
	$C_0$	$C_1$, $C_2$, $C_3$, $C_6$		$C_4$	$C_5$	Initial relative volume, $V_0$
		0.0		0.4	0.4	1.0
TNT	MAT_HIGH_EXPLOSIVE_BURN
	Density, $\rho$	Detonation velocity, $D$		Chapman-Jouget pressure, $P_{CJ}$		Initial internal energy, $E_0$
	1.63	0.693		0.21		0.07
	EOS_JWL
	Constants representing characteristics of the explosive				OMEGA, $\omega$	Initial relative volume, $V_0$
	$A$	$B$	$R_1$	$R_2$
	3.712	0.03231	4.15	0.95	0.3	1.0
Soil	MAT_SOIL_AND_FOAM
	Density, $\rho$	Elastic modulus, $E$	Shear modulus, $G$		Bulk modulus, $K_v$
	1.78	$E=$ 43.75 MPa	$G=$ 14.68 MPa		$K_v=$ 729 MPa
	Shear-yield surface parameters
	$a_0=$ 9.025×10-12, $a_1=$ 1.38×10-11, $a_2=$ 5.28×10-12

In LS-DYNA code, the shear failure principle of the soil and rock model reflects the relationship between the average stress and failure strength. The yielding function is shown in Eq. (3), the second invariant of stress deviation is shown in Eq. (4), and the axial yield stress is shown in Eq. (5). The shear yield surface parameters are $a_0=c^2$, $a_1=2c\mathrm{t}\mathrm{a}\mathrm{n}\varphi$, and $a_2={\tan }^2\varphi$. A soil composition model must simulate the consolidation or compaction characteristics of porous materials. Material characteristics were used to analyze the dynamic reaction of soil after a blast force. A material model using MAT_SOIL_AND_FOAM was used to calculate the shear failure of the porous material by using the average stress and failure strength [20].

The basis of material failure was defined with the element erosion criterion MAT_ADD_EROSION. Soil is a low-resistance material. According to the relationship between material stress and strain, this study defined the maximum shear strain (${\gamma }_1$) and failure shear strain (${\gamma }_{\mathrm{m}\mathrm{a}\mathrm{x}}$) as elemental erosion failure conditions. As the conditions of the erosion criterion were established, and the material failure criterion ${\gamma }_1\ge {\gamma }_{\mathrm{m}\mathrm{a}\mathrm{x}}$ was satisfied, the element was deleted from the calculation model [20]:

where $f$ is the yield function; $J_2$ is the second invariant of deviatoric stress; $p$ is the pressure; $a_0$-$a_2$ are the shear yield surface parameters; and ${\sigma }_y$ is the yield stress.

3.1. Numerical analysis model for verifying the explosion experiments

The main objective of this study was to analyze the dynamic reaction and deformation of soil subject to blast force. An explosion generates shock waves. We therefore employed shock wave pressure and ground acceleration to analyze the characteristics of shock waves. In the ground explosion experiment, the free-field near-ground blast pressure and ground acceleration were measured to verify the numerical analysis model. The relative error (%) of the numerical analysis and experimental results was equal to (numerical analysis value − experimental measurement value) / experimental measurement value × 100 %.

To verify the blast effect, the maximum ground acceleration in the vertical direction 200 cm from the explosion center was used. Fig. 3 shows the vertical acceleration curves from the explosion experiment measured at 200 cm from the explosion center over time and those from the numerical analysis. The maximum values from the experiment and numerical analysis were 6.806 gal (m/s²) and 6.216 gal (m/s²), respectively; the relative error was −8.679 %. The near-ground blast pressure shock waves 200 and 300 cm from the explosion center were used to analyze the shock wave energy. Figs. 4 and 5 respectively show the blast pressure curves for the points 200 and 300 cm from the explosion center. At 200 cm from the explosion center, the experimental and numerical analysis values for the near-ground blast pressure were 50.297 and 45.826 kPa, respectively; the relative error was −8.890 %. At 300 cm from the explosion center, the experimental and numerical analysis values for the near-ground blast pressure were 22.534 and 20.718 kPa, respectively; the relative error was −8.059 %. These results agree with those of another relevant study [11].

Fig. 3Acceleration curves of explosion experiment and numerical analysis over time, measured perpendicular from the ground 200 cm above the explosive source

Figs. 3, 4, and 5 depict the duration curve and relative error for numerical analysis and experimental results. The reception of surface acceleration and explosion pressure signals in the onsite experiment was triggered by the instrument, and the numerical analysis results presented the duration curve of the particle. Due to the complexity of the explosive environment, material parameters and equations of state of the numerical analysis cannot fully reflect the experimental conditions and complex environmental factors. Although the numerical analysis and experimental results were not consistently accurate, the analyzed duration curve displayed a consistent trend, and the overall analysis results were in agreement with the facts.

The research results demonstrated that a fluid’s dynamic behavior and its interaction with solids could be analyzed using 3D SOLID164 elements together with numerical fluid-solid interaction models and the MMALE algorithm, thereby solving the dynamic analysis problems of material, geometric, and contact nonlinearity. The results indicated that the 3D spatial solid structure model and analysis parameters were reasonable. The results were then used to analyze the dynamic effects of an explosion on soil.

Fig. 4Blast pressure curves of explosion experiment and numerical analysis over time, measured 200 cm from the explosive source

Fig. 5Blast pressure curves of explosion experiment and numerical analysis over time, measured 300 cm from the explosive source

3.2. Shock and seismic waves

A shock wave transmits energy, in the form of pressure, to the ground directly or through the air. When the ground explodes, the blast impacts the air and causes a shock wave, creating air-induced ground motion and triggering directly transmitted ground shocks. Fig. 6 illustrates the transmission of ground explosion shock waves over time. Fig. 6(a-h) depict the free-field shock wave rapidly moving spherically outwards under pressure. In a few milliseconds, the energy attenuated. During the process of pressure transmission, if no reflection wave is encountered, a dynamic pressure similar to wind pressure is formed in the area behind the blast. Following an explosion, the region of influence from the center outward can be divided into the compression zone, fracture zone, and vibration zone. Influenced by the shock wave, high-intensity stress impact is generated. From the explosion’s center, shock, compression, and seismic waves emanate. A seismic wave triggered by an explosion is triggered mainly by the air impact. When the surface of the ground explodes, a compression wave and elastic seismic wave in the soil can simultaneously be measured. The compression wave has abundant energy, which triggers vibration in particles on the ground surface. Soil is a material with low tensile strength. Under the influence of shock waves, particle vibrations accelerate, thereby increasing the shear force in the soil, resulting in particle strain and displacement.

Fig. 6Free-zone shock wave transmission over time in the ground surface explosion

a)$t=$ 249.95 μs

b)$t=$ 499.44 μs

c)$t=$999.39 μs

d)$t=$1,508.10 μs

e)$t=$ 2,009.60 μs

f)$t=$2,499.80 μs

g)$t=$2,999.80 μs

h)$t=$5,000.00 μs

Fig. 7Plastic strain on soil over time

a)$t=$ 249.95 μs

b)$t=$ 499.44 μs

c)$t=$ 999.39 μs

d)$t=$ 1,508.10 μs

e)$t=$ 2,009.60 μs

f)$t=$ 2,499.80 μs

g)$t=$ 2,999.80 μs

i) Experimental plan and numerical model

h)$t=$ 5,000.00 μs

3.3. Material failure and deformations in surrounding soil

Fig. 7 illustrates the soil plastic strain over time. The maximum strain values in the $X$, $Y$, and $Z$ directions were 0.013, 0.0068, and 0.012, respectively. In dynamic reaction analysis, the pressure, stress, and strain on a material are evaluated. The large strain generated by an explosion shock wave influences the internal stress of the soil, leading to soil strain. The blast influences the stress condition of the soil, resulting in particle displacement and deformation. When the stress on the soil exceeds its yield strength, its internal structural stability is affected. Through analysis of material deformation, plastic strain over time, displacement, and increase in shear strain, the stability state of soil was determined. The results revealed that after the yield strength of the soil was exceeded due to the blast, plastic deformation occurred in the material. The plastic strain analysis process revealed that the failure in the soil caused by the blast loading effect began in the area of plastic deformation. Uneven stress distribution reduced the shear strength of the material. Moreover, blast loading triggered a compression wave in the soil, leading to increased strain in the material. When this strain exceeded the yield strength of the soil, displacement occurred, and when this displacement reached a critical level, the soil deformed and broke. Therefore, the blast affected the stress state of the soil, influencing its strength and stability.

Fig. 8 indicates the crater effects caused by ground contact explosions of 113.389 g (0.25 lb) of TNT. Soil can withstand pressure more easily than it can tension and shear forces. When researchers analyze the dynamic reaction of soil against ground explosions, they must consider that a ground explosion forms a hypocenter that transmits energy outward through the soil. The shock wave generated by the explosion is transmitted through soil, and most of the energy is consumed in the deformation and failure of the soil as well as through pressure on the surrounding soil. The blast increases the shear force in the soil and reduces its resistance to shear force. The soil structure is influenced by an uneven distribution of stress, and stress becomes concentrated. When the material can no longer withstand the stress, failure occurs. The analysis results revealed that the explosion caused a sudden change in pressure, resulting in severe deformation and failure in the soil. The shock wave caused failure, starting with local plastic deformation. When the second invariant of stress, $J_2$, is larger than the yield strength of the material, permanent deformation and failure of the material occurs.

Fig. 8Crater caused by ground contact explosion

The stress and strain of soil exhibit nonlinear and plastic characteristics. This study triggered an earthquake with an explosion and used the processes of soil deformation and failure to analyze the stress and strain behavior of the soil. The research results showed that the shock wave triggered vibration in the particles on the surface, which altered the stress state of the soil. When the stress exceeded the yield strength of the soil, the stability of the internal structure was disturbed. Furthermore, uneven stress on the soil caused stress concentration, reducing the soil’s shear strength. Failure originated from local plastic deformation. The subsequent plastic deformation and increased displacement were the crucial conditions prompting soil failure. In addition, before soil failure, a shear failure band was formed. When displacement began to increase substantially, the sudden displacement indicated that the soil had failed. The shock wave energy level near the explosion center was high, and the soil bore a high-energy compression wave. The material failed primarily because of tension and compression. Therefore, shear failure is an essential condition of soil failure. The explosion shock wave at the ground surface affected the material’s stability and caused destruction of the medium. Subject to the direct action of the stress wave, plastic deformation occurred outside of the crater, creating a shear failure zone that inflicted irreversible deformation.