Parametric Study on Contact Explosion Resistance of Steel Wire Mesh Reinforced Geopolymer Based Ultra-High Performance Concrete Slabs Using Calibrated Continuous Surface Cap Model

: This paper conducts a parametric analysis on the response of geopolymer-based ultra-high-performance concrete (G-UHPC) slabs reinforced with steel wire mesh (SWM) subjected to contact explosions using the validated Continuous Surface Cap (CSC) model. Firstly, based on the available experimental data, the CSC model parameters, which account for the yield surface, damage formulation, kinematic hardening, and strain rate effect, were comprehensively developed for G-UHPC. The modiﬁed CSC model was initially assessed by comparing the quasi-static test results of G-UHPC. Then, the numerical modeling was performed on 200 mm thick SWM-reinforced G-UHPC slabs against 0.4 kg and 1.0 kg TNT contact explosions. The fair agreement between the numerical and experimental data concerning the local damage of the slabs was reported to demonstrate the applicability of the material and structural models. With the validated numerical models, a parametric study was further acted upon to explore the contribution of the variables of SWM, slab thickness, and TNT equivalence on the local damage and energy evolution of G-UHPC slabs subjected to contact blasts. Moreover, based on simulation results from the parametric study, an updated empirical model was derived to evaluate the local damage pattern and internal energy absorption rate of SWM-reinforced G-UHPC slabs.


Introduction
In recent years, hazardous loads induced by various blast accidents have been a serious concern in civil and military constructions [1][2][3][4][5]. Under close-in or contact explosions, normal concrete typically exhibits brittle failure with highly localized damage due to the low material strength and ductility [6][7][8]. To strengthen the resistance of concrete structures against close-in and contact explosions, advanced concrete materials, e.g., engineered cementitious composites (ECC), ultra-high-performance concrete (UHPC), etc., have undergone extensive development and investigation over the last few decades [9][10][11]. Compared to normal concrete, the advanced concrete materials exhibited superior dynamic material strength, toughness, and energy absorption capacity due to the existence of fiber reinforcement, which could effectively restrain the propagation of blast stress waves within the concrete structures [12].
With the global focus on sustainability in civil and military constructions, a geopolymer binder system is deemed to be an alternative to the conventional Portland cement in the concrete mixture [13][14][15] since the manufacture of geopolymer can consume less energy and produce lower levels of carbon dioxide. After continuous development in these decades, the raw materials of geopolymer include natural minerals, e.g., kaolin, metakaolin,

A Brief Overview of CSC Model in LS-DYNA
In LS-DYNA, the CSC model mainly includes the failure surface, damage functions, and strain rate formulations. The failure surface is combined with the shear failure surface and cap hardening surface, in which these two surfaces are continuously and smoothly connected, as presented in Figure 1. The failure surface is characterized by three stress tensor invariants (J 1 , J 2 and J 3 ) and the cap hardening parameter κ, which is given by: where J 1 = 3P denotes the first stress tensor invariant; J 2 = 1 2 S ij S ij denotes the second invariant of the deviatoric stress tensor; J 3 = 1 2 S ij S jk S ki denotes the third stress tensor invariant; P is the pressure; S ij , S jk and S ki are the deviatoric stress tensors, the indexes i, j, k refer to elements of this stress tensor; κ is the value of J 1 at the intersection of the shear damage surface and the cap surface; is the Rubin three-invariant decay coefficient; F c is the abbreviation of hardened cap function; the full formula definition is listed in Equation (5). F f is the abbreviation of the shear surface function; the full formula definition is listed in Equation (2).
Buildings 2022, 12, x FOR PEER REVIEW 3 of 24 subjected to 0.4 kg and 1.0 kg TNT contact explosions. Through comparing the numerical and test data concerning the local damage, the applicability of the calibrated CSC model and structural model was validated. Using the verified numerical models, a parametric study on SWM-reinforced G-UHPC slabs was further carried out concerning the TNT equivalent, slab thickness, and SWM variables (i.e., number of layers, space between adjacent wires, and steel wire diameter). Furthermore, an updated empirical model was proposed to evaluate the local damage pattern and internal energy absorption rate of SWMreinforced G-UHPC slabs against contact explosions. As mentioned above, although several experiments on G-UHPC members with SWM reinforcement subjected to dynamic loads have been reported, the corresponding numerical simulations were urgently required to further interpret the efficiency and effectiveness of SWM parameters. Moreover, rare studies were conducted to adopt the CSC model to simulate the G-UHPC in contact explosions. This paper detailed performed the CSC model calibration, the structural model validation, and the parametric study of SWM variables.

A Brief Overview of CSC Model in LS-DYNA
In LS-DYNA, the CSC model mainly includes the failure surface, damage functions, and strain rate formulations. The failure surface is combined with the shear failure surface and cap hardening surface, in which these two surfaces are continuously and smoothly connected, as presented in Figure 1. The failure surface is characterized by three stress tensor invariants (J 1 , J 2 and J 3 ) and the cap hardening parameter κ, which is given by: (1) where J 1 = 3P denotes the first stress tensor invariant; J 2 = 1 2 S ij S ij denotes the second invariant of the deviatoric stress tensor; J 3 = 1 2 S ij S jk S ki denotes the third stress tensor invariant; P is the pressure; S ij , S jk and S ki are the deviatoric stress tensors, the indexes i, j, k refer to elements of this stress tensor; κ is the value of J 1 at the intersection of the shear damage surface and the cap surface; ℜ is the Rubin three-invariant decay coefficient; Fc is the abbreviation of hardened cap function; the full formula definition is listed in Equation (5). Ff is the abbreviation of the shear surface function; the full formula definition is listed in Equation (2). The shear surface consists of a triaxial compression (TXC) surface, a triaxial extension (TXE) surface, and a triaxial torsion (TOR) surface. Figure 2 presents the triaxial compression and tension meridians, respectively. The controlling equation for the TXC is given by: F f (J 1 ) = αλ exp (-βJ 1 ) + θJ 1 (2) The shear surface consists of a triaxial compression (TXC) surface, a triaxial extension (TXE) surface, and a triaxial torsion (TOR) surface. Figure 2 presents the triaxial compression and tension meridians, respectively. The controlling equation for the TXC is given by: where α, β, λ and θ are the material constants as derived through fitting the data from the triaxial compressive tests. where α, β, λ and θ are the material constants as derived through fitting the data from the triaxial compressive tests.
(a) (b) The TXC intensity ratio of the Rubin scale function at any stress level defines the triaxial torsional and tensile radial functions [40]. The controlling functions for TOR and TXE are given by: where Q 1 is the TOR/TXC intensity ratio and Q 2 is the TXE/TXC intensity ratio; α 1 , β 1 , λ 1 and θ 1 are the TOR input parameters; α 2 , β 2 , λ 2 and θ 2 are the TXE input parameters. The cap hardening surface is a two-part function that is either uniform or elliptical. When the stress state is in the region of the tension or very low confining pressure, the cap hardening surface function is uniform. When the stress state is in the region of the low to high confining pressure, the cap hardening function is elliptical. The cap hardening surface is defined as: where κ 0 denotes the value of J 1 corresponding to the initial intersection between the shear failure and cap hardening surface; X (κ) denotes the intersection between the cap hardening surface and the axis of J 1 , which is given by: where R is the ellipticity ratio of the cap hardening surface. The cap moves to reflect the plastic volume change, in which the cap expansion and contraction are in the light of an isotropic hardening rule [41], which is given by: where ε ν p denotes the plastic volumetric strain; W denotes the maximum plastic volumetric strain; D1 and D2 are the material constants; X 0 is the initial position of the cap if κ = κ 0 . The damage formulations define the strain softening and modulus reduction behavior. The damage criterion is on the basis of the damage energy release rate method as established by Simo et al. [42], which is defined as: The TXC intensity ratio of the Rubin scale function at any stress level defines the triaxial torsional and tensile radial functions [40]. The controlling functions for TOR and TXE are given by: where Q 1 is the TOR/TXC intensity ratio and Q 2 is the TXE/TXC intensity ratio; α 1 , β 1 , λ 1 and θ 1 are the TOR input parameters; α 2 , β 2 , λ 2 and θ 2 are the TXE input parameters. The cap hardening surface is a two-part function that is either uniform or elliptical. When the stress state is in the region of the tension or very low confining pressure, the cap hardening surface function is uniform. When the stress state is in the region of the low to high confining pressure, the cap hardening function is elliptical. The cap hardening surface is defined as: where κ 0 denotes the value of J 1 corresponding to the initial intersection between the shear failure and cap hardening surface; X (κ) denotes the intersection between the cap hardening surface and the axis of J 1 , which is given by: where R is the ellipticity ratio of the cap hardening surface. The cap moves to reflect the plastic volume change, in which the cap expansion and contraction are in the light of an isotropic hardening rule [41], which is given by: where ε p ν denotes the plastic volumetric strain; W denotes the maximum plastic volumetric strain; D 1 and D 2 are the material constants; X 0 is the initial position of the cap if κ = κ 0 .
The damage formulations define the strain softening and modulus reduction behavior. The damage criterion is on the basis of the damage energy release rate method as established by Simo et al. [42], which is defined as: where d denotes a scalar damage parameter that can be converted to a stress tensor σ ij in the absence of damage or a stress tensor σ ij in the presence of damage. The strain rate is accommodated by applying the viscoplastic algorithm, which was extended by Simo et al. [43] on the basis of the Duvaut-Lions formulation. The strength calculation for the strain rate effect is based on the dynamic increase factor (DIF), in which the DIF value equals the ratio of the dynamic strength of the concrete material to quasi-static strength. The relevant strain rate effect equation is defined as: where σ vp ij is viscoplastic stress; σ T ij is elastic stress; σ P ij is plastic stress; γ = ∆t/η 1+∆t/η , η is fluidity coefficient parameter; ∆t is the time step. The DIF for the compressive and tensile strengths is expressed as: where f c is uniaxial compressive strength; f t is direct tensile strength; f c,d is dynamic compressive strength; f t,d is dynamic tensile strength; . ε is strain rate; η c and η t respectively denote the fluidity coefficient for compressive and tensile strengths, and are defined as follows: where η oc and N c are strain rate parameters for compression; η ot and N t are strain rate parameters for tension.

Bulk and Shear Moduli
The bulk modulus K represents the response of material under hydrostatic pressure. The shear modulus G represents the response of material under shear stress. In this study, these two moduli are not independent, and for isotropic materials, they related to Young's modulus E c through the equations: where E c denotes the elastic modulus of G-UHPC, which was determined to be 33 GPa as obtained from the elastic stage in the uniaxial compressive stress-strain curve [31]; ν is the Poisson's ratio of G-UHPC, which was determined to be 0.2 according to the previous review study by Ranjbar and Zhang [44].

Triaxial Compression Surface
The TXC parameters were determined by fitting four types of stress states obtained from the triaxial compression tests. As shown in Figure 2a, the four stress states included A1: triaxial tensile state, B1: biaxial tensile state, C1: uniaxial compression state and D1: triaxial compressive state. The triaxial and biaxial tensile strengths of G-UHPC were set approximately equal to the uniaxial tensile strength [45]. To gain the triaxial compressive strength of G-UHPC at different confining pressure, Khan et al. [46] and Haider et al. [47] carried out the triaxial compression tests on geopolymer concrete with the uniaxial compressive strength up to 90 MPa. Through selecting the triaxial compression test data, as presented in Figure 3, the fitting formula is expressed as: where ; σ 0 = (σ 1 +σ 2 +σ 3 )

3
; σ 1 , σ 2 and σ 3 are the principal stress. is the Poisson's ratio of G-UHPC, which was determined to be 0.2 according to the previous review study by Ranjbar and Zhang [44].

Triaxial Compression Surface
The TXC parameters were determined by fitting four types of stress states obtained from the triaxial compression tests. As shown in Figure 2a, the four stress states included A1: triaxial tensile state, B1: biaxial tensile state, C1: uniaxial compression state and D1: triaxial compressive state. The triaxial and biaxial tensile strengths of G-UHPC were set approximately equal to the uniaxial tensile strength [45]. To gain the triaxial compressive strength of G-UHPC at different confining pressure, Khan et al. [46] and Haider et al. [47] carried out the triaxial compression tests on geopolymer concrete with the uniaxial compressive strength up to 90 MPa. Through selecting the triaxial compression test data, as presented in Figure 3, the fitting formula is expressed as: where ; σ 1 , σ 2 and σ 3 are the principal stress. Triaxial compressive test data [46,47] and fitting curve.
In order to obtain the parameters in the compression meridian, a minimum of five stress states was suggested [48]. Apart from points C1 and D1 in Figure 2a, the other three stress states at various pressure levels (i.e., J 1 = 1.5f c , 3f c and 5f c ) were selected for the triaxial compression states, which could be determined by Equation (17). Based on the aforementioned stress states, the parameters were determined via fitting to the compression failure surface at various strength levels using the least-square method [39]. Figure 4 shows the fitting results, wherein the corresponding parameters of α, β, λ and θ with respect to f c were determined by the following equations: θ = -9 × 10 6 f c 2 + 0.0018f c + 0.0062 (21) Figure 3. Triaxial compressive test data [46,47] and fitting curve.
In order to obtain the parameters in the compression meridian, a minimum of five stress states was suggested [48]. Apart from points C1 and D1 in Figure 2a, the other three stress states at various pressure levels (i.e., J 1 = 1.5 f c , 3 f c and 5 f c ) were selected for the triaxial compression states, which could be determined by Equation (17). Based on the aforementioned stress states, the parameters were determined via fitting to the compression failure surface at various strength levels using the least-square method [39]. Figure 4 shows the fitting results, wherein the corresponding parameters of α, β, λ and θ with respect to f c were determined by the following equations:

Cap Hardening Surface
The cap hardening surface parameters consist of initial cap position (X0), cap aspect ratio (R), maximum effective plastic volume strain (W), linear shape parameter (D1), and quadratic shape parameter (D2). The cap parameters could be determined based on the

Cap Hardening Surface
The cap hardening surface parameters consist of initial cap position (X 0 ), cap aspect ratio (R), maximum effective plastic volume strain (W), linear shape parameter (D 1 ), and quadratic shape parameter (D 2 ). The cap parameters could be determined based on the hydrostatic pressure-volumetric strain curve obtained from the hydrostatic compression tests. Since limited hydrostatic tests were conducted on G-UHPC, this study adopted the available test data on PC-UHPC to assume the cap parameters of the CSC model. This assumption was in view of the similarity between these two types of UHPC materials in terms of density, Poisson's ratio, and porosity [49,50]. Based on hydrostatic compression data by Williams et al. [51] and Xu et al. [52], the relationship between X 0 and f c was fitted, which is shown in Figure 5a. The cap ellipticity R can be derived based on the previous study by Jiang et al. [41], the expressions of X 0 and R are as follows: hydrostatic pressure-volumetric strain curve obtained from the hydrostatic compression tests. Since limited hydrostatic tests were conducted on G-UHPC, this study adopted the available test data on PC-UHPC to assume the cap parameters of the CSC model. This assumption was in view of the similarity between these two types of UHPC materials in terms of density, Poisson's ratio, and porosity [49,50]. Based on hydrostatic compression data by Williams et al. [51] and Xu et al. [52], the relationship between X 0 and f c was fitted, which is shown in Figure 5a. The cap ellipticity R can be derived based on the previous study by Jiang et al. [41], the expressions of X 0 and R are as follows: The maximum effective plastic volume strain W was assumed to be the porosity of G-UHPC [38]. In the current study, the value of W was determined to be 0.0445 for G-UHPC in accordance with the proposed value by Williams et al. [51]. After fitting with the test data using Equation (8), as presented in Figure 5b, the values of D1 and D2 were determined to be 3.368 × 10 − ³ MPa −1 and 1.134 × 10 −6 MPa −2 , respectively. The default values of 1.0 and 0 were respectively taken for N H and C H as suggested in the original CSC model [35].

Damage
The model parameters governing the damage evolution mainly include the compressive fracture energy (GFC), tension fracture energy (GFT), shear fracture energy (GFS), ductile shape softening (B), brittle shape softening (D), shear-to-compression transition parameter (pwrc), shear-to-tension transition parameter (pwrt) and pressure softening parameter (pmod). As recommended in the original CSC model [35], the default values of 0, 1.0 and 5.0 were respectively set for pmod, pwrt and pwrc, and GFS was assumed to be equal to GFT. The damage parameters GFT and D govern the strain-softening behavior under tension. The damage parameters GFC and B govern the strain-softening behavior under compression. The fracture energy (GFT and GFC) in the CSC model is sensitive to the element size [38]. To determine GFC, B, GFS, and D, the single-element model under uniaxial compression and the full-scale model under four-point bending were built, as presented in Figures 6a and 7a, respectively. The details about the single-element simulation setup could be found in the study by Guo et al. [37], and the four-point bending simulation setup, as well as the loading scheme, could be found in the study by Wei et al. [39]. The element size adopted in the models was 10 mm, which was consistent with that in the contact explosion simulations for the slabs. Figure 6b presents the influence of GFC and B on the uniaxial compressive stress-strain curve. The compressive strain hardening The maximum effective plastic volume strain W was assumed to be the porosity of G-UHPC [38]. In the current study, the value of W was determined to be 0.0445 for G-UHPC in accordance with the proposed value by Williams et al. [51]. After fitting with the test data using Equation (8), as presented in Figure 5b, the values of D 1 and D 2 were determined to be 3.368 × 10 −3 MPa −1 and 1.134 × 10 −6 MPa −2 , respectively. The default values of 1.0 and 0 were respectively taken for N H and C H as suggested in the original CSC model [35].

Damage
The model parameters governing the damage evolution mainly include the compressive fracture energy (GFC), tension fracture energy (GFT), shear fracture energy (GFS), ductile shape softening (B), brittle shape softening (D), shear-to-compression transition parameter (pwrc), shear-to-tension transition parameter (pwrt) and pressure softening parameter (pmod). As recommended in the original CSC model [35], the default values of 0, 1.0 and 5.0 were respectively set for pmod, pwrt and pwrc, and GFS was assumed to be equal to GFT. The damage parameters GFT and D govern the strain-softening behavior under tension. The damage parameters GFC and B govern the strain-softening behavior under compression. The fracture energy (GFT and GFC) in the CSC model is sensitive to the element size [38]. To determine GFC, B, GFS, and D, the single-element model under uniaxial compression and the full-scale model under four-point bending were built, as presented in Figures 6a and 7a, respectively. The details about the single-element simulation setup could be found in the study by Guo et al. [37], and the four-point bending simulation setup, as well as the loading scheme, could be found in the study by Wei et al. [39]. The element size adopted in the models was 10 mm, which was consistent with that in the contact explosion simulations for the slabs. Figure 6b presents the influence of GFC and B on the uniaxial compressive stress-strain curve. The compressive strain hardening behavior of the CSC model increases with the increase in parameters GFC and B. After fitting with the test results as reported in Ref. [31], GFC and B were determined to be 700 Pa·m and 10, respectively. Figure 7b presents the effect of GFT and D on the force-displacement curve, and the tensile strain hardening behavior of the CSC model also increases with the increase in parameters GFT and D. After fitting with the test result as reported in Ref. [24], GFT and D were taken as 40 Pa·m and 5000, respectively. behavior of the CSC model increases with the increase in parameters GFC and B. After fitting with the test results as reported in Ref. [31], GFC and B were determined to be 700 Pa·m and 10, respectively. Figure 7b presents the effect of GFT and D on the force-displacement curve, and the tensile strain hardening behavior of the CSC model also increases with the increase in parameters GFT and D. After fitting with the test result as reported in Ref. [24], GFT and D were taken as 40 Pa·m and 5000, respectively.  DIF c = 1.0 for ε < ε sc 1.0 + 0.01ε 0.88 for ε ≥ ε sc (24) where ε sc = 1.2 × 10 5 s 1 . For the tensile strength:  behavior of the CSC model increases with the increase in parameters GFC and B. After fitting with the test results as reported in Ref. [31], GFC and B were determined to be 700 Pa·m and 10, respectively. Figure 7b presents the effect of GFT and D on the force-displacement curve, and the tensile strain hardening behavior of the CSC model also increases with the increase in parameters GFT and D. After fitting with the test result as reported in Ref. [24], GFT and D were taken as 40 Pa·m and 5000, respectively.  DIF c = 1.0 for ε < ε sc 1.0 + 0.01ε 0.88 for ε ≥ ε sc (24) where ε sc = 1.2 × 10 5 s 1 . For the tensile strength:  In compression:

Strain Rate Effect
In tension: overt = 0.64f t ε 0.55 (31) srate and repow were taken as the default value of 1.0, as recommended in the original CSC (a) (b)

Validation of Numerical Models
To validate the numerical models, contact explosion tests were numerically simulated on the multi-layered SWM-reinforced G-UHPC slabs with the utilization of the calibrated CSC model. The following contents briefly introduce the contact explosion tests, finite element modeling, material models excluding concrete, and the comparison between the numerical and experimental results. For the compressive strength:

Contact Explosion Tests
where . ε sc = 1.2 × 10 −5 s −1 . For the tensile strength: where . ε st = 1.0 × 10 −6 s −1 . To prevent overprediction of the DIF values at high strain rates, the LS-DYNA user manual [32] suggests interposing a cut-off value for the DIF-strain rate curves. In the current study, referring to the recommendation from the LS-DYNA user manual [32], the cut-off values were set at a strain rate of 300 s −1 in both compression and tension. Then, the active input key parameters, including η oc , η ot , N c , N t , overc and overt (Equations (11)-(14)), were derived by fitting data on the proposed curve made based on Equations (24) and (25), which are given by: In compression: In tension: srate and repow were taken as the default value of 1.0, as recommended in the original CSC model [35].

Validation of Numerical Models
To validate the numerical models, contact explosion tests were numerically simulated on the multi-layered SWM-reinforced G-UHPC slabs with the utilization of the calibrated CSC model. The following contents briefly introduce the contact explosion tests, finite element modeling, material models excluding concrete, and the comparison between the numerical and experimental results.

Contact Explosion Tests
In the contact explosion tests, as reported in Refs [31,64], two 200 mm × 1500 mm × 1500 mm G-UHPC slabs reinforced with multi-layered SWM were tested. The first slab was reinforced by 10-layer SWM reinforcement (labeled as G-UHPC-10SWM) and subjected to 0.4 kg of TNT explosive. The second slab was reinforced by 20-layer SWM reinforcement (labeled as G-UHPC-20SWM) and subjected to 1.0 kg of TNT explosive. Further, 1 mm diameter 304 stainless steel wires with a 10 mm interval for SWM. The SWM layers were designed to be evenly distributed, and the SWM layers next to each other were spaced the same. The concrete cover depth was 10 mm. The yield strength and elastic modulus of the reinforcement were 500 MPa and 200 GPa, respectively. The 28-day uniaxial compressive strength of G-UHPC without steel fibers was 90 MPa. Figure 9 presents the schematic diagram of the slab, TNT explosive, and SWM reinforcement. In the tests, the SWM-reinforced G-UHPC slabs were simply supported by the square steel frame with a side length of 1500 mm located upon four trapezoidal concrete piers [64]. The setup of the contact explosion tests can be seen in Figure 10. In the contact explosion tests, as reported in Refs [31,64], two 200 mm × 1500 mm × 1500 mm G-UHPC slabs reinforced with multi-layered SWM were tested. The first slab was reinforced by 10-layer SWM reinforcement (labeled as G-UHPC-10SWM) and subjected to 0.4 kg of TNT explosive. The second slab was reinforced by 20-layer SWM reinforcement (labeled as G-UHPC-20SWM) and subjected to 1.0 kg of TNT explosive. Further, 1 mm diameter 304 stainless steel wires with a 10 mm interval for SWM. The SWM layers were designed to be evenly distributed, and the SWM layers next to each other were spaced the same. The concrete cover depth was 10 mm. The yield strength and elastic modulus of the reinforcement were 500 MPa and 200 GPa, respectively. The 28-day uniaxial compressive strength of G-UHPC without steel fibers was 90 MPa. Figure 9 presents the schematic diagram of the slab, TNT explosive, and SWM reinforcement. In the tests, the SWM-reinforced G-UHPC slabs were simply supported by the square steel frame with a side length of 1500 mm located upon four trapezoidal concrete piers [64]. The setup of the contact explosion tests can be seen in Figure 10.   Figure 11 shows the 3D finite element models of G-UHPC-10SWM and G-UHPC-20SWM subjected to TNT contact explosions through the multi-material ALE algorithm. The G-UHPC slab and steel frame were modeled utilizing a single-point integration algorithm and the eight-node hexahedron Lagrangian elements, wherein a 10 mm mesh size was chosen for the slab after a convergence test. The SWM was modeled by the Hughes-Liu beam elements with a mesh size of 10 mm and the cross-sectional integration algorithm. The TNT explosive and air were modeled by the Euler elements with a mesh size   Figure 11 shows the 3D finite element models of G-UHPC-10SWM and G-UHPC-20SWM subjected to TNT contact explosions through the multi-material ALE algorithm. The G-UHPC slab and steel frame were modeled utilizing a single-point integration algorithm and the eight-node hexahedron Lagrangian elements, wherein a 10 mm mesh size was chosen for the slab after a convergence test. The SWM was modeled by the Hughes-Liu beam elements with a mesh size of 10 mm and the cross-sectional integration algorithm. The TNT explosive and air were modeled by the Euler elements with a mesh size  Figure 11 shows the 3D finite element models of G-UHPC-10SWM and G-UHPC-20SWM subjected to TNT contact explosions through the multi-material ALE algorithm. The G-UHPC slab and steel frame were modeled utilizing a single-point integration algorithm and the eight-node hexahedron Lagrangian elements, wherein a 10 mm mesh size was chosen for the slab after a convergence test. The SWM was modeled by the Hughes-Liu beam elements with a mesh size of 10 mm and the cross-sectional integration algorithm. The TNT explosive and air were modeled by the Euler elements with a mesh size of 10 mm after a convergence test. CONTACT_AUTOMATIC_SURFACE_TO_SURFACE was adopted to define the contact between the concrete slab and steel frame, in which the static and dynamic frictional coefficients were both set to 0.5. The blast wave propagation, the interaction between the blast wave and the G-UHPC slab, as well as the full constraint between the slab and SWM reinforcement, were defined through CONSTRAINED_LAGRANGE_IN_SOLID. Non-reflecting boundary conditions were set on the surrounding sides of the air domain to avoid reflected stress waves.

Finite Element Modelling
Buildings 2022, 12, x FOR PEER REVIEW 13 of 24 of 10 mm after a convergence test. CONTACT_AUTOMATIC_SURFACE_TO_SURFACE was adopted to define the contact between the concrete slab and steel frame, in which the static and dynamic frictional coefficients were both set to 0.5. The blast wave propagation, the interaction between the blast wave and the G-UHPC slab, as well as the full constraint between the slab and SWM reinforcement, were defined through CONSTRAINED_LA-GRANGE_IN_SOLID. Non-reflecting boundary conditions were set on the surrounding sides of the air domain to avoid reflected stress waves.

G-UHPC
In accordance with the calibration procedures as introduced in Section 3, the CSC model parameters were determined for G-UHPC, which are listed in Table 1. It was generally accepted that a key problem with Lagrangian elements under high-rate loading is the large distortion that causes computational overflow [65]. To solve this problem, an erosion algorithm activated by MAT_ADD_EROSION was adopted to remove large distortion elements once the erosion criterion was reached. In the current study, the

G-UHPC
In accordance with the calibration procedures as introduced in Section 3, the CSC model parameters were determined for G-UHPC, which are listed in Table 1. It was generally accepted that a key problem with Lagrangian elements under high-rate loading is the large distortion that causes computational overflow [65]. To solve this problem, an erosion algorithm activated by MAT_ADD_EROSION was adopted to remove large distortion elements once the erosion criterion was reached. In the current study, the maximum principal strain of 0.1 was ultimately used as the erosion criterion for G-UHPC via several trial simulations.

TNT Explosive and Air
The TNT explosive was modeled using MAT_HIGH_EXPLOSIVE_BURN. EOS_JWL was employed for the TNT explosive to simulate the detonation procedure. The pressure in EOS is defined as: where V denotes the relative volume; E denotes the specific internal energy; ω, A, B 0 , R 1 and R 2 are the active input parameters. Table 2 lists the model and EOS parameters for the TNT explosive. The air was modeled by MAT_NULL, along with EOS_LINEAR_POLYNOMIAL. The pressure P in EOS was defined by: where C 0 , C 1 , C 2 , C 3 , C 4 , C 5 and C 6 are user-defined constants; µ = 1/V − 1, in which V is the relative volume; E 0 is the internal energy per initial volume. Table 2 lists the model and EOS parameters for the air.

SWM and Square Steel Frame
MAT_PIECEWISE_LINEAR_PLASTICITY was employed for SWM. The strain rate effect on SWM was considered via the incorporation of the curve of DIF and strain rate. The DIF values could be calculated via the equation as proposed by Malvar et al. [66], which is given by: where . ε is the strain rate; α = 0.074 − 0.04 f y /60, in which f y is the yield strength of SWM. Due to the negligible deformation after the explosion, the square steel frame was treated as a rigid body. Thus MAT_RIGID was adopted to model the square steel frame. Table 3 lists the model parameters for SWM and square steel frame.

Figures 12 and 13
show the numerical results of the local damage for G-UHPC-10SWM subjected to 0.4 kg TNT contact explosion and G-UHPC-20SWM subjected to 1.0 kg TNT contact explosion, respectively, wherein the effective plastic strain ranging from 0.3 (no damage) to 1 (complete damage) represents the local damage to the slab. As observed from the results of physical tests, the front, rear, and side faces of the slab suffered evident damage along with cracks induced by the stress waves. The steel frame to support the slab also exacerbated the damage close to the edges, especially on the rear face of the slab, as shown in Figure 13. The numerical results were then compared with the test results concerning the crack distribution, crater, and scabbing damage, from which fair agreement between the numerical and test results was exhibited. Therefore, the material and numerical model used in this study are suitable for reproducing the local damage of G-UHPC slabs with SWM reinforcement under contact explosions.
served from the results of physical tests, the front, rear, and side faces of the slab suffered evident damage along with cracks induced by the stress waves. The steel frame to support the slab also exacerbated the damage close to the edges, especially on the rear face of the slab, as shown in Figure 13. The numerical results were then compared with the test results concerning the crack distribution, crater, and scabbing damage, from which fair agreement between the numerical and test results was exhibited. Therefore, the material and numerical model used in this study are suitable for reproducing the local damage of G-UHPC slabs with SWM reinforcement under contact explosions.  Additionally, the energy evolutions of SWM under two blast scenarios, i.e., 0.4 kg and 1.0 kg, were obtained from the numerical simulations, which are shown in Figure 14. Achieved from the numerical results, the total energies induced by 0.4 kg and 1.0 kg TNT explosives were 1750 kJ and 4380 kJ, respectively. To quantitatively evaluate the effect of Additionally, the energy evolutions of SWM under two blast scenarios, i.e., 0.4 kg and 1.0 kg, were obtained from the numerical simulations, which are shown in Figure 14. Achieved from the numerical results, the total energies induced by 0.4 kg and 1.0 kg TNT explosives were 1750 kJ and 4380 kJ, respectively. To quantitatively evaluate the effect of the SWM volumetric content on the energy absorption capacity under contact explosions, the internal energy absorption rate of SWM (E * r ) was defined, which is given by: where E s and E m denote the internal energies absorbed by the SWM reinforcement and G-UHPC slab, respectively. As shown in Figure 15,

Parametric Study
To comprehensively understand the contribution of the SWM reinforcement on the contact explosion resistance of G-UHPC slabs, a parametric study was further performed by using the validated finite element models. Table 4 lists the parametric study matrix for SWM-reinforced G-UHPC slabs subjected to contact explosions, in which the SWM layer number (N), space between adjacent wires per layer (L), steel wire diameter (d), TNT equivalent (W) and slab thickness (T) were included for investigation. The yield strength of SWM and the uniaxial compressive strength of G-UHPC were kept as 500 MPa, and 90

Parametric Study
To comprehensively understand the contribution of the SWM reinforcement on the contact explosion resistance of G-UHPC slabs, a parametric study was further performed by using the validated finite element models. Table 4 lists the parametric study matrix for SWM-reinforced G-UHPC slabs subjected to contact explosions, in which the SWM layer number (N), space between adjacent wires per layer (L), steel wire diameter (d), TNT equivalent (W) and slab thickness (T) were included for investigation. The yield strength

Parametric Study
To comprehensively understand the contribution of the SWM reinforcement on the contact explosion resistance of G-UHPC slabs, a parametric study was further performed by using the validated finite element models. Table 4 lists the parametric study matrix for SWM-reinforced G-UHPC slabs subjected to contact explosions, in which the SWM layer number (N), space between adjacent wires per layer (L), steel wire diameter (d), TNT equivalent (W) and slab thickness (T) were included for investigation. The yield strength of SWM and the uniaxial compressive strength of G-UHPC were kept as 500 MPa, and 90 MPa, respectively, and 125 test cases were performed for each case (W-T-N 1 -10-1.0, W-T-N 2 -L-1.0, and W-T-N 2 -10-d) listed in Table 4.   Table 5 lists the classifications of the local damage for the G-UHPC slabs after contact explosions [67], wherein three damage levels were included concerning the spall depth. Based on the damage classifications, the result on the local damage of SWM reinforced G-UHPC slabs obtained from the parametric study was normalized, which is defined as: where H 1 is the crater depth on the front face; H 2 is the scabbing depth on the rear face; H 3 is the tunnel depth. Table 5. Local damage classifications of concrete slabs after contact explosions [67].

Damage Level Damage Description Damage Scheme
Mild A very shallow spalling to a third of the slab thickness and no change in the slab thickness to a few noticeable fissures  Table 5 lists the classifications of the local damage for the G-UHPC slabs after contact explosions [67], wherein three damage levels were included concerning the spall depth. Table 5. Local damage classifications of concrete slabs after contact explosions [67].

Damage Level Damage Description Damage Scheme
Mild A very shallow spalling to a third of the slab thickness and no change in the slab thickness to a few noticeable fissures Moderate From more than a third to more than two thirds of the slab thickness spalling Severe From just over two third of the slab thickness spall to breach Figure 16 shows the influence of N, L, and d on the relationship between the normalized damage d* and T/W 1/3 . It can be observed that the increase in N and d can decrease the normalized damage of G-UHPC slabs subjected to contact explosions under the identical T/W 1/3 , but the increase in L will exacerbate the normalized damage. Based on the normalized results from the parametric study (as shown in Figure 17), a series of empirical equations concerning T, W, N, L, and d were proposed to evaluate the local damage (D*) of SWM-reinforced G-UHPC slabs subjected to contact explosions, which are given by: (37) For the mild damage level:  Table 5 lists the classifications of the local damage for the G-UHPC slabs after contact explosions [67], wherein three damage levels were included concerning the spall depth. Table 5. Local damage classifications of concrete slabs after contact explosions [67].

Damage Level Damage Description Damage Scheme
Mild A very shallow spalling to a third of the slab thickness and no change in the slab thickness to a few noticeable fissures Moderate From more than a third to more than two thirds of the slab thickness spalling Severe From just over two third of the slab thickness spall to breach Figure 16 shows the influence of N, L, and d on the relationship between the normalized damage d* and T/W 1/3 . It can be observed that the increase in N and d can decrease the normalized damage of G-UHPC slabs subjected to contact explosions under the identical T/W 1/3 , but the increase in L will exacerbate the normalized damage. Based on the normalized results from the parametric study (as shown in Figure 17), a series of empirical equations concerning T, W, N, L, and d were proposed to evaluate the local damage (D*) of SWM-reinforced G-UHPC slabs subjected to contact explosions, which are given by: (37) For the mild damage level:  Table 5 lists the classifications of the local damage for the G-UHPC slabs after contact explosions [67], wherein three damage levels were included concerning the spall depth. Table 5. Local damage classifications of concrete slabs after contact explosions [67].

Damage Level Damage Description Damage Scheme
Mild A very shallow spalling to a third of the slab thickness and no change in the slab thickness to a few noticeable fissures Moderate From more than a third to more than two thirds of the slab thickness spalling Severe From just over two third of the slab thickness spall to breach Figure 16 shows the influence of N, L, and d on the relationship between the normalized damage d* and T/W 1/3 . It can be observed that the increase in N and d can decrease the normalized damage of G-UHPC slabs subjected to contact explosions under the identical T/W 1/3 , but the increase in L will exacerbate the normalized damage. Based on the normalized results from the parametric study (as shown in Figure 17), a series of empirical equations concerning T, W, N, L, and d were proposed to evaluate the local damage (D*) of SWM-reinforced G-UHPC slabs subjected to contact explosions, which are given by: (37) For the mild damage level:  Table 5 lists the classifications of the local damage for the G-UHPC slabs after contact explosions [67], wherein three damage levels were included concerning the spall depth. Table 5. Local damage classifications of concrete slabs after contact explosions [67].

Damage Level Damage Description Damage Scheme
Mild A very shallow spalling to a third of the slab thickness and no change in the slab thickness to a few noticeable fissures Moderate From more than a third to more than two thirds of the slab thickness spalling Severe From just over two third of the slab thickness spall to breach Figure 16 shows the influence of N, L, and d on the relationship between the normalized damage d* and T/W 1/3 . It can be observed that the increase in N and d can decrease the normalized damage of G-UHPC slabs subjected to contact explosions under the identical T/W 1/3 , but the increase in L will exacerbate the normalized damage. Based on the normalized results from the parametric study (as shown in Figure 17), a series of empirical equations concerning T, W, N, L, and d were proposed to evaluate the local damage (D*) of SWM-reinforced G-UHPC slabs subjected to contact explosions, which are given by: (37) For the mild damage level:  Table 5 lists the classifications of the local damage for the G-UHPC slabs after contact explosions [67], wherein three damage levels were included concerning the spall depth. Table 5. Local damage classifications of concrete slabs after contact explosions [67].

Damage Level Damage Description Damage Scheme
Mild A very shallow spalling to a third of the slab thickness and no change in the slab thickness to a few noticeable fissures Moderate From more than a third to more than two thirds of the slab thickness spalling Severe From just over two third of the slab thickness spall to breach Figure 16 shows the influence of N, L, and d on the relationship between the normalized damage d* and T/W 1/3 . It can be observed that the increase in N and d can decrease the normalized damage of G-UHPC slabs subjected to contact explosions under the identical T/W 1/3 , but the increase in L will exacerbate the normalized damage. Based on the normalized results from the parametric study (as shown in Figure 17), a series of empirical equations concerning T, W, N, L, and d were proposed to evaluate the local damage (D*) of SWM-reinforced G-UHPC slabs subjected to contact explosions, which are given by: For the mild damage level:  Table 5 lists the classifications of the local damage for the G-UHPC slabs after contact explosions [67], wherein three damage levels were included concerning the spall depth. Table 5. Local damage classifications of concrete slabs after contact explosions [67].

Damage Level Damage Description Damage Scheme
Mild A very shallow spalling to a third of the slab thickness and no change in the slab thickness to a few noticeable fissures Moderate From more than a third to more than two thirds of the slab thickness spalling Severe From just over two third of the slab thickness spall to breach Figure 16 shows the influence of N, L, and d on the relationship between the normalized damage d* and T/W 1/3 . It can be observed that the increase in N and d can decrease the normalized damage of G-UHPC slabs subjected to contact explosions under the identical T/W 1/3 , but the increase in L will exacerbate the normalized damage. Based on the normalized results from the parametric study (as shown in Figure 17), a series of empirical equations concerning T, W, N, L, and d were proposed to evaluate the local damage (D*) of SWM-reinforced G-UHPC slabs subjected to contact explosions, which are given by: For the mild damage level:  Table 5 lists the classifications of the local damage for the G-UHPC slabs after contact explosions [67], wherein three damage levels were included concerning the spall depth. Table 5. Local damage classifications of concrete slabs after contact explosions [67].

Damage Level Damage Description Damage Scheme
Mild A very shallow spalling to a third of the slab thickness and no change in the slab thickness to a few noticeable fissures Moderate From more than a third to more than two thirds of the slab thickness spalling Severe From just over two third of the slab thickness spall to breach Figure 16 shows the influence of N, L, and d on the relationship between the normalized damage d* and T/W 1/3 . It can be observed that the increase in N and d can decrease the normalized damage of G-UHPC slabs subjected to contact explosions under the identical T/W 1/3 , but the increase in L will exacerbate the normalized damage. Based on the normalized results from the parametric study (as shown in Figure 17), a series of empirical equations concerning T, W, N, L, and d were proposed to evaluate the local damage (D*) of SWM-reinforced G-UHPC slabs subjected to contact explosions, which are given by: L 0.18 (37) For the mild damage level:  Table 5 lists the classifications of the local damage for the G-UHPC slabs after contact explosions [67], wherein three damage levels were included concerning the spall depth. Table 5. Local damage classifications of concrete slabs after contact explosions [67].

Damage Level Damage Description Damage Scheme
Mild A very shallow spalling to a third of the slab thickness and no change in the slab thickness to a few noticeable fissures Moderate From more than a third to more than two thirds of the slab thickness spalling Severe From just over two third of the slab thickness spall to breach Figure 16 shows the influence of N, L, and d on the relationship between the normalized damage d* and T/W 1/3 . It can be observed that the increase in N and d can decrease the normalized damage of G-UHPC slabs subjected to contact explosions under the identical T/W 1/3 , but the increase in L will exacerbate the normalized damage. Based on the normalized results from the parametric study (as shown in Figure 17), a series of empirical equations concerning T, W, N, L, and d were proposed to evaluate the local damage (D*) of SWM-reinforced G-UHPC slabs subjected to contact explosions, which are given by: L 0.18 (37) For the mild damage level: For the severe damage level: D* < 1.08 (38c) Figure 16 shows the influence of N, L, and d on the relationship between the normalized damage d* and T/W 1/3 . It can be observed that the increase in N and d can decrease the normalized damage of G-UHPC slabs subjected to contact explosions under the identical T/W 1/3 , but the increase in L will exacerbate the normalized damage. Based on the normalized results from the parametric study (as shown in Figure 17), a series of empirical equations concerning T, W, N, L, and d were proposed to evaluate the local damage (D * ) of SWM-reinforced G-UHPC slabs subjected to contact explosions, which are given by: For the mild damage level : D * ≥ 1.70 (38a) For the moderate damage level : 1.08 ≤ D * < 1.70 (38b) For the severe damage level : D * < 1.08 (38c) where the unit system is 'm-kg'. What is also noted is that the empirical equations were valid for 0.10 m ≤ T ≤ 0.30 m and 0.2 kg ≤ W ≤ 2.4 kg. The increase in N, L, and d can improve the volumetric content of SWM (V s ) in the G-UHPC slab. Based on the numerical results from the parametric study, the relationship among E * r , V s and T/W 1/3 was obtained, which is shown in Figure 18. An empirical equation concerning V s and T/W 1/3 was then proposed to predict E * r under various contact explosion scenarios, which is given by: Moderate From more than a third to more than two thirds of the slab thickness spalling Severe From just over two third of the slab thickness spall to breach Figure 16 shows the influence of N, L, and d on the relationship between the normalized damage d* and T/W 1/3 . It can be observed that the increase in N and d can decrease the normalized damage of G-UHPC slabs subjected to contact explosions under the identical T/W 1/3 , but the increase in L will exacerbate the normalized damage. Based on the normalized results from the parametric study (as shown in Figure 17), a series of empirical equations concerning T, W, N, L, and d were proposed to evaluate the local damage (D*) of SWM-reinforced G-UHPC slabs subjected to contact explosions, which are given by:  among E r * , V s and T/W 1/3 was obtained, which is shown in Figure 18. An empirical equation concerning V s and T/W 1/3 was then proposed to predict E r * under various contact explosion scenarios, which is given by:

Conclusions
In this study, the dynamic response of G-UHPC slabs with SWM reinforcement against contact explosions was numerically performed. To reproduce the existing tests, the Continuous Surface Cap (CSC) model for G-UHPC and the structural model for contact explosion tests in nonlinear finite element software LS-DYNA have been validated. The results noted that structural models based on the developed CSC model reasonably capture the local damage of experimental specimens, which highlighted that the CSC Figure 18. Nonlinear fitting plot of internal energy absorption rate for SWM.

Conclusions
In this study, the dynamic response of G-UHPC slabs with SWM reinforcement against contact explosions was numerically performed. To reproduce the existing tests, the Continuous Surface Cap (CSC) model for G-UHPC and the structural model for contact explosion tests in nonlinear finite element software LS-DYNA have been validated. The results noted that structural models based on the developed CSC model reasonably capture the local damage of experimental specimens, which highlighted that the CSC model could be adopted to simulate G-UHPC with SWM reinforcement in contact explosion tests. With the validated finite element models, a parametric study was further carried out to explore the influence of SWM layer number, space between adjacent wires, wire diameter, TNT equivalent, and slab thickness on the local damage levels of G-UHPC slabs with SWM reinforcement under contact explosions. Parametric analysis results show the influence of SWM variables on local damage; based on the results from the parametric study, a series of empirical equations concerning the aforementioned variables were established. These empirical equations can be used to identify the damage mode of G-UHPC slabs with SWM reinforcement under contact explosion and predict the internal energy absorption rate of SWM. The proposed empirical equations could provide a general guideline for designing G-UHPC with SWM reinforcement against contact explosions.