Methodology for Predicting the Structural Response of RPC-Filled Steel Tubular Columns under Blast Loading

Abstract

Reactive powder concrete-filled steel tube (RPC-FST) is a critical supporting component of large-span, high-rise, and heavy-load structures. The collapse of RPC-FST may occur under explosive load. Therefore, investigation of the dynamic response is essential for understanding the mechanisms of collapse. In this article, the numerical simulation of reactive powder concrete (RPC) adopted the modified Holmquist–Johnson–Cook (HJC) model and the modified Karagozian and Case (K&C) model. The dynamic response of RPC-FST columns under explosive load is analyzed based on arbitrary Lagrange–Euler (ALE) method. The proposed model is verified by experimental results. Results show that the modified HJC model and modified K&C model can be applied to simulate the dynamic response of RPC-FST columns under explosive load. As compared with the modified HJC model, the modified K&C model has more accurate results. This phenomenon mainly accounts for the lack of ultimate strain of RPC (EFMIN). To analyze the reliability sensitivity of RPC-FST, an efficient probability analysis method is proposed based on the Kriging model and Monte Carlo simulation (MCS). The proposed method considers five nonlinear factors, including weight and distance of TNT, height and section diameter of RPC-FST, and steel tube thickness. Finally, the sensitivity of each factor is evaluated. Results show that TNT weight greatly influences the reliability of the RPC-FST, followed by TNT distance, RPC-FST height, RPC-FST section diameter, and steel tube thickness. In addition, the RPC-FST dynamic response analysis method based on the Kriging model and MCS can improve the calculation efficiency by more than 200 times compared with the ALE method.

1. Introduction

Terrorist activities increased considerably over the past two decades. According to a survey, explosions account for 50% of terrorist attacks, which has raised deep concerns regarding the safety of civilian and military infrastructure [1,2,3]. Therefore, the evaluation of reactive powder concrete-filled steel tube (RPC-FST) performance and the improvement of RPC-FST reliability under explosive loads are crucial in construction.

The explosion process is a typical nonlinear transient dynamic issue. The methods of applying explosion load are directly related to the accuracy of the numerical analysis. In general, the methods of applying explosion load can be roughly divided into three types: arbitrary Lagrange–Euler (ALE), load blast enhanced (LBE), and pressure-time history method. Among these, the LBE-based method defines explosion load based on the empirical blast loading function, which is derived from experimental explosion data [4,5]. Compared with the LBE-based method, the pressure-time history method defines the explosion load according to an exponentially decaying function (such as Friedlander’s equation), which can be replaced with an equivalent triangular function when the negative phase effect can be ignored [6,7]. The LBE-based method and pressure-time history method have high computational efficiency. Nevertheless, these simplified methods lack a profound understanding of the functional mechanism between air and RPC-FST. Contrary to the aforementioned methods, the ALE-based method can accurately restore local details of blast conditions considering fluid–solid interaction (FSI) theory. Based on this method, Thai [8] assessed the residual strength of fiber-reinforced concrete columns under explosion load. Kostopoulos [9] investigated the blast resistance of a composite foam-core sacrificial cladding for steel-reinforced concrete structures. Wu [10] investigated the dynamic response of ultra-high-performance cement-based composite-filled steel tube (UHPCC-FST) under close-range explosion.

Besides the methods of applying explosion load, the material model is also a significant factor in precisely simulating the dynamic response of a structure under explosion load. The LS-DYNA program offers many dynamic constitutive models for simulating concrete materials, such as the Holmquist–Johnson–Cook (HJC) model and the Karagozian and Case (K&C) model. The HJC model is mainly suitable for material under large strain, high strain rate, and high pressure [11]. Based on the HJC material model, Zhu [12] established a theoretical model for shaped charge jets’ penetration into concrete targets. Kristoffersen [13] investigated the ballistic perforation resistance of concrete slabs. Wan [14] calibrated parameters to investigate the blast resistance of ultra-high-performance concrete (UHPCs) slabs. It should be noted that the modified HJC model parameters are suitable for UHPCs slabs with 10.5% porosity of coarse aggregate. However, the applicability of HJC material model parameters modified for RPC materials with low porosity without coarse aggregate needs to be deeply investigated.

Different from the HJC material model, the K&C material model has great advantages in describing the material under damage evolution, restraint effect, and shear expansion [15]. Based on the K&C material model, Wang [16] estimated the residual axial load-bearing capacity of UHPCC-FST specimens subjected to contact explosion. Liu [17] investigated the damage evolution of reinforced concrete piers with carbon-fiber-reinforced polymer under contact explosion. Zhang [18] calibrated the model parameters of UHPCs with high-velocity projectile impact experimental data. Meanwhile, reactive powder concrete (RPC) is a low-porosity material due to its dense microstructure [19]. However, the calculation accuracies of the HJC material model and K&C material model have not been deeply investigated for RPC with low porosity. Therefore, the application scopes of the K&C material model and HJC material model require more discussion.

There are still many challenges in calculating the dynamic response of RPC-FST columns under explosion loads. The correct description of the RPC material model under steel tube constraint has research value. In addition, the influence of different types of fiber and different fiber contents in RPC on the dynamic response of an RPC-FST column under explosion load is still worth discussing. Although most scholars have investigated the mechanical properties of RPC-FST with experiments [20,21,22], the changes in the mechanical properties of RPC-FST columns under explosion loads are not clear. Therefore, dynamic response analysis of RPC-FST columns under explosion load is still in the exploratory stage.

In practical scenarios, it is impossible to accurately describe and deterministically control the nonlinear factors of the explosion process. Therefore, the characteristic parameters for estimating the stability of the RPC-FST are random. For example, manufacturing error causes randomness in the dimensions of the RPC-FST column. Moreover, the randomness of weight and distance of TNT cause uncertainty of the explosion load. Therefore, it is critical to establish the probability analysis for the explosion process. Hussein [23] investigated the reliability of composite wood-sand-wood blast walls based on Monte Carlo simulation (MCS). Song [24] presented the reliability analysis method for a steel frame structure under explosion load based on Bayesian theory. Ding [25] established an effective reliability evaluation framework to predict the failure risk of steel frame structures under explosion load. Shi [26] predicted the damage of reinforced concrete wall panels under various threats based on MCS. Momeni [27] presented an improved calculation method, based on the MCS method and finite element approach, to evaluate the minimum safe scaled distance for steel columns under dynamic blast loads. Beyond any dispute, the application of these above methods can provide guidance for the probability analysis of the explosion process. Nevertheless, the reliability sensitivity analysis of the RPC-FST column based on MCS still requires further investigation.

In this article, the critical axial deformation is employed as the stability threshold of RPC-FST. Meanwhile, the application scopes of the K&C material model and HJC material model are theoretically illustrated with numerical calculation. To improve the computational efficiency of sensitivity assessment, the surrogate Kriging model is adopted to replace the modified K&C model. The rest of this paper is organized as follows: In Section 2, the explosion load is applied based on the ALE method. In addition, the modified K&C material model for low-porosity RPC materials is established. In Section 3, it is verified that the modified K&C model has advantages over the modified HJC model in material properties of low-porosity RPC. In addition, a numerical example establishes the basis for the Kriging model. In Section 4, the RPC-FST sensitivity analysis of the random explosion process is proposed. The sensitivity of the output with respect to each factor is evaluated.

2. RPC-FST Model Based on ALE Method

As the core component of blast-resistant structures, RPC-FST plays a significant role in protecting buildings in explosively hazardous areas from fatal damage. Therefore, the dynamic response of RPC-FST under explosive load has attracted wide attention. In the process of dynamic response analysis, the robustness of the model directly affects the accuracy of numerical calculation. In this section, the K&C model is modified based on the low porosity of RPC material. In addition, the finite element model (FEM) of RPC-FST under explosion load is established based on the ALE method, as shown in Figure 1.

2.1. Fluid–Solid Interaction of RPC-FST

The approach for solving FSI problems is the ALE formulation for the fluid domain and the Lagrangian formulation for the structure domain. This approach is called the ALE method, and it has been intensively used for problems involving small and large structure displacements with no topological changes in the structure. Therefore, the ALE method from LS-DYNA provides a possibility to model multi-phase highly dynamic problems. The nonlinear dynamic analysis program LS-DYNA is used to analyze the dynamic response of RPC-FST columns under explosion load. A numerical model composed of four parts is established. The TNT and air domain are treated as ALE parts, while the RPC-FST columns are treated as Lagrangian parts.

2.2. Establishing the Model

The Jones–Wilkins–Lee equation of state (JWL EOS) defines the pressure as a function of the relative volume V and internal energy E₀. TNT has been established as material by using the JWL EOS method. It can be expressed as:

$$P=A(1-\frac{\omega }{R_{1}V})e^{-R_{1}V}+B(1-\frac{\omega }{R_{2}V})e^{-R_{2}V}+\frac{\omega E_0}{V}$$

(1)

where A, B, R₁, R₂, and ω are the parameters related to the TNT type.

The air domain is established with the equation of state, which represents the relationship between the element pressure and the internal energy given in Equation (2).

$$P=C_0+C_1\mu +C_2{\mu }^2+C_3{\mu }^3+(C_4+C_5\mu +C_6{\mu }^2)E$$

(2)

where P is the pressure, C_i (i = 1, 2, …,6) values are constant, and E is the internal energy per unit volume.

The plastic kinematic model is an elastic-plastic model with kinematic and isotropic hardening. Therefore, the plastic kinematic model is selected for steel tube materials under impact load. The strain rate effect on steel material can be incorporated by the Cowper and Symonds (CS) model (Equation (3)), which scales the yield stress with the factor.

$$DIF=1+{(\frac{\dot{\epsilon }}{C})}^{1/p}$$

(3)

where $\dot{\epsilon }$ is the strain rate. C and p denote strain rate parameters for the CS strain rate model. In this study, C and p are set as 6844 s⁻¹ and 3.91, respectively. This combination accurately predicts the dynamic response of concrete-filled steel tubular columns, as verified by tensile experiments [28]. The constants mentioned above will be explained in detail later in Section 3.

RPC has extremely high compressive strength and fracture energy [29,30,31]. There is a lack of accurate model parameters to describe the characteristics of RPC materials in RPC-FST columns. The K&C model is considered a promising dynamic constitutive model [15], and thus can be used to reveal the material properties under impact and explosion loads. The K&C model is modified according to the characteristics of RPC without coarse aggregate.

The K&C model is specifically used to calculate the concrete structural response under blast and impact loadings. It is comprised of two parts—volume response and deviatoric response.

In volume response, the volume change of materials under different pressures is observed by applying stress. The equation of state *EOS_TABULATED_COMPRESSION correlates the pressure p and the volumetric strain ɛ_V. In deviatoric response, three failure surface strength models—yield failure surface Δσ_y, maximum failure surface Δσ_m, and residual failure surface Δσ_r—have to be characterized via a_0y, a_1y, a_2y, a₀, a₁, a₂, a_1f, and a_2f. The deviatoric stresses remain elastic during loading or reloading until the stress achieves the yield failure surface. Then, the deviatoric stress further increases until the maximum failure surface is reached. Beyond this stage, the response can be perfectly plastic or softened to the residual failure surface. The damage function captures the hardening and softening behavior of three failure surfaces. According to the characteristics of 5% porosity of RPC material in RPC-FST, the modified K&C model process is shown in Figure 2.

The dynamic increase factor (DIF) relationship proposed by Hou [32] is adopted, which considers the strain rate effect of RPC. The DIF of the RPC compressive strength can be expressed as in Equation (4).

$$DIF=a{\left(\frac{\dot{\epsilon }}{{\dot{\epsilon }}_s}\right)}^b\left[\dot{\epsilon }\right]\le \dot{\epsilon }\le 317s^{-1}$$

(4)

where ($\dot{\epsilon }$) is the critical strain rate of DIF. If the strain rate is lower than ($\dot{\epsilon }$), the strain rate effect on the compressive strength can be neglected. Moreover, (${\dot{\epsilon }}_s$) = 1 s⁻¹ is used to make indexes dimensionless, where ($\dot{\epsilon }$), a, and b are the constant coefficients for RPC. a and b are set as 0.291 and 0.310, respectively.

It is worth noting that detonation products are transmitted in the air with high-frequency shock waves, and their physical process and dynamic response are exceedingly complex. To reveal this process accurately, the determination of air domain size is essential. According to TNT distance, TNT size, RPC-FST, and other model parameters in the practical project, the high-frequency shock wave shows a large amplitude, and the behavior of RPC-FST displacement is geometrically nonlinear. Therefore, the size of the air domain was chosen as 3000 mm × 1000 mm × 8000 mm. The non-reflection boundary keyword is used for air boundary conditions to prevent wave reflection from the boundary. The mesh size denotes 20 mm, and the termination time denotes 50 ms.

3. Model Validation and Numerical Examples

In this section, it will be revealed that the modified K&C model has advantages over the modified HJC model for the description of low-porosity RPC material properties in RPC-FST under explosion load. Modified parameters of the K&C model, damage function η(λ), and equation of state are shown in

Table 1. Modified K&C model parameters.

Parameter	Value	Parameter	Value
Mass density R0 (g/mm3)	0.002	Poisson’s ratio PR	0.2
Uniaxial tensile strength T (MPa)	9	Failure surface parameter a0	3.3 × 107
Failure surface parameter a1	0.45	Failure surface parameter a2	7.4 × 10−10
Compressive damage parameter B1	1.6	Tensile damage scaling exponent B2	1.4
Damage scaling coefficient B3	1.15	Initial yield surface coefficient a0y	2.5 × 107
Initial yield surface coefficient a1y	0.63	Initial yield surface coefficient a2y	2.3 × 10−9
Residual failure surface coefficient a1f	0.44	Residual failure surface coefficient a2f	1.1 × 10−9
Aggregate diameter LOCWID (mm)	0.03	Post peak dilatancy decay EDROP	1
Stretch factor Sλ	100	Fractional dilatancy OMEGA	0.5
Unit conversion for length RSIZE	0.04	Unit conversion for stress UCF	145

,

Table 2. Modified damage function η(λ).

Scale Factor λ				Damage Function η
Parameter	Value	Parameter	Value	Parameter	Value	Parameter	Value
λ1	0	λ8	3.2 × 10−4	η1	0	η8	0.85
λ2	8.0 × 10−6	λ9	5.2 × 10−4	η2	0.85	η9	0.67
λ3	2.4 × 10−5	λ10	5.7 × 10−4	η3	0.97	η10	0.57
λ4	4.0 × 10−5	λ11	1	η4	0.99	η11	0.37
λ5	5.6 × 10−5	λ12	1.0 × 101	η5	1	η12	0.37
λ6	7.2 × 10−5	λ13	1.0 × 1010	η6	0.99	η13	0.37
λ7	8.8 × 10−5			η7	0.97

and

Table 3. Modified equation of state.

Volumetric Strain		Pressure		Bulk Unloading Modulus
Parameter	Value	Parameter	Value (MPa)	Parameter	Value (MPa)
εv1	0	C1	0	K1	20,540
εv2	0.0015	C2	30.81	K2	20,541
εv3	0.002	C3	39.43	K3	20,830
εv4	0.0025	C4	49.03	K4	21,870
εv5	0.0029	C5	73.05	K5	26,020
εv6	0.0043	C6	94.94	K6	30,190
εv7	0.0101	C7	103.32	K7	34,340
εv8	0.0305	C8	196.3	K8	37,480
εv9	0.0513	C9	296.08	K9	84,330
εv10	0.0726	C10	420.06	K10	102,700

[18]. The modified HJC model is shown in Appendix A. Air equation of state parameters are shown in

Table 4. Air equation of state parameters.

Parameter	C0	C1	C2	C3	C4	C5	C6	E (J/mm3)	V
Value	0	0	0	0	0.4	0.4	0	0.3	1

. The parameters of TNT are shown in

Table 5. TNT parameters.

Parameter	Value	Parameter	Value
TNT density R0 (g/mm3)	1.63 × 10−3	Detonation velocity D (m/s)	6.93 × 103
Parameter A	3.71 × 105	Parameter B	3.23 × 103
Parameter R1	4.15	Parameter R2	0.95
Parameter w	0.3	Detonation energy E0 (J/mm3)	7000
Initial relative volume V	1	Chapman-Jouget pressure PCJ	2.70 × 104

. The material properties of the steel tube are shown in

Table 6. Material properties of steel tube.

Parameter	Value	Parameter	Value
Mass density R (g/mm3)	7.85 × 10−3	Young’s modulus E (MPa)	2.05 × 105
Poisson’s ratio PR	0.3	Yield stress SIGY (MPa)	318.5
Hardening parameter BETA	0	Strain rate parameter SRC	6844
Effective plastic strain FS	0.25	Strain rate parameter SRP	3.91

.

3.1. Model Competition

In order to profoundly investigate the accuracy of the modified K&C model and the modified HJC model, we decided to verify the applicability of the two modified material models to RPC materials according to the explosion experiment conducted by Wang [33] in the recent report. The compressive strength of RPC is 135.1 MPa, and steel tube yield strength is 318.5 MPa. The top of the column limits the displacement in the x, y, and z directions, and the bottom end is a fixed-end constraint. Model parameters are shown in

Table 7. Model parameters.

Specimen	TNT Weight (kg)	TNT Height (mm)	TNT Distance (mm)	Steel Tube Thickness (mm)
C-15	4	250	221	6
C-16	4	950	221	6
C-17	7	955	230	6
C-18	10	945	275	6

.

The modified K&C model, HJC model, modified HJC model, and experiment comparison are shown in Figure 3. Among them, the material parameters of the HJC material model are determined according to Wan [14]. It should be pointed out that the residual displacement results of the experiments are used for the verification. The results show that the modified K&C model parameters have better accuracy in calculating dynamic response than the modified HJC model parameters and parameters of the modified HJC model proposed by Wan [14]. The displacement time history curve of the modified K&C model converges to the experimental value faster through vibration. The increase of RPC damage degree in RPC-FST under near-range explosion load is the main reason for the inaccurate numerical calculation.

The structural response under explosion loading is influenced by the constitutive models of the material, which will bring uncertainty to structural protection. The result of residual displacement with the HJC model is obviously lower than the experimental value. This phenomenon is mainly because the parameters are modified for RPC containing coarse aggregate [14]. The porosity of RPC containing coarse aggregate is generally 10.5%, while the porosity of RPC without coarse aggregate is generally 5%, which leads to a certain deviation in the expression of the equation of state under this working condition. In addition, determining the parameters of the damage evolution equation based only on steel fiber content needs to be further discussed. Therefore, the dynamic response analysis of the modified HJC model is carried out. The results show that the modified HJC model has obvious advantages over the HJC model. This phenomenon is mainly because the modified HJC model considers the influence of porosity on the equation of state. Furthermore, the effects of the water/cement ratio, steel fiber content, and steel fiber length on the damage evolution equation are considered. However, the advantages of the modified K&C model are still obvious. Through analysis, it should be pointed out that the slope of the assumed failure surface is characterized by uncertainty, the values of ultimate strains ɛ_x obtained by the uniaxial compression test are too large, and the damage parameter D is small (Figure A2). Therefore, the deviation of the damage evolution equation leads to the calculation results of the HJC model being lower than the experiment results.

From Figure 4, it is found that the dynamic responses of RPC-FST calculated based on the modified K&C model are consistent with the experimental results. From the experiment, the steel tube surface of the crater exhibited many tiny scars. The results show that some elements fail on the blast surface, which is consistent with the experimental phenomenon. Due to some solid elements of steel tube reaching ultimate strength after bearing the explosion load and being unable to continue bearing the load, the result is solid element failure.

3.2. Numerical Example

According to the model comparison, it is found that the modified K&C model is more suitable for the RPC material model in RPC-FST. Therefore, the RPC material model is established by using the modified K&C model. In addition, a group of samples after the Latin hypercube sampling plan is selected as the structural parameters of RPC-FST. This process establishes the basis for the reliability sensitivity analysis of RPC-FST in Section 4. The dynamic response analysis of the RPC-FST column under explosive load is carried out based on FSI theory. The bottom and top of the RPC-FST columns are restrained in the x, y, and z directions, whereas all the other degrees of freedom can move freely.

The ALE method can evaluate a blast scenario with more detail, including explosive detonation. Therefore, the velocity of the shock wave and the displacement of RPC-FST are extracted. Figure 5 illustrates the functional mechanism between air and RPC-FST.

The evolutions of RPC-FST dynamic response under explosion load can be roughly divided into four phases: process A, process B, process C, and process D, as shown in Figure 5a. Meanwhile, the details of each process are shown in Figure 5b. The results show that the shock wave did not reach the blast surface of RPC-FST at 0.5 ms. The maximum displacement of the RPC-FST occurred at 2.5 ms. In addition, the maximum displacement (process C) in the RPC-FST displacement time history curve is selected to establish the Kriging surrogate model.

4. Sensitivity Analysis of RPC-FST

Irreversible RPC-FST damage will emerge during the explosion process. Meanwhile, the accumulation of the damage process will decrease the durability of the RPC-FST [22]. Furthermore, the characteristic parameters for estimating the damage of the RPC-FST are random due to manufacturing and measurement errors [34,35]. Therefore, the quantification and propagation of randomness of the explosion process are critical in RPC-FST design. In addition, the low calculation efficiency of the modified K&C model under explosion load is considered. The surrogate Kriging model is adopted to replace the modified K&C model for sensitivity analysis of RPC-FST under random explosion load.

For probabilistic analysis, the reliability of RPC-FST under explosion load can be defined as a multi-dimensional integral:

$$P_s={\displaystyle \int \cdots {\displaystyle {\int }_{G(x)>0}{f}_X(x)}}dx$$

(5)

where

$$G\left(x\right)=\left[\left(\frac{H}{2}\right)\times \tan \left(\theta \right)\right]-\epsilon \left(x\right)$$

(6)

where G(x) is the limit state function, which can be established based on Section 2; x = (H, W, t, L, D)^T is a vector of independent random variables used to calculate the dynamic response of the RPC-FST system; X is a random parameter vector corresponding to x; H is the height of RPC-FST; W is the weight of TNT; t describes the thickness of steel tube; L represents the distance between TNT and RPC-FST; D is the section diameter of RPC-FST; and θ denotes the deflection angle of RPC-FST at the upper damage limit [36]. According to TM5-1300 [3], the damage levels of concrete-filled steel tube columns can be divided into low damage, moderate damage, and high damage (

Table 8. Damage level of concrete-filled steel tube column.

Damage Level	Low Damage	Moderate Damage	High Damage
θ	0–2°	2–5°	5–12°

). Therefore, the failure probability is calculated for the low damage of a 0–2° column, as shown in Figure 6. In order to calculate the sensitivity of the allowable deflection angle between 0 and 2°, we decided to choose 1.4° as the allowable deflection angle. ε(x) is the maximum displacement based on FSI theory. However, due to the high dimensionality and complicated integrand of Equation (5), both analytical and direct numerical methods cannot be applied directly.

4.1. Monte Carlo Simulation (MCS) Method

MCS is an accurate technique to solve complex multi-dimensional integrals [37]. The problem solution is transformed into the expectation of the probability model. Meanwhile, statistical analysis is adopted to evaluate the probability approximation of the multi-dimensional integral [38]. The conversion of Equation (5) based on MCS can be expressed as:

$$P_f={\displaystyle \underset{F}{\int }{I}_F\left(x\right)f_X\left(x\right)}dx=E\left(I_F\left(x\right)\right)$$

(7)

where F is the failure domain of RPC-FST, and I_F(x) is the indication function of the failure domain, which can be obtained by:

$$I_F\left(X\right)=\left\{\begin{array}{l}0\hspace{1em}G\left(X\right)>0\\ 1\hspace{1em}G\left(X\right)\le 0\end{array}\right.$$

(8)

Therefore, the failure probability of RPC-FST can be calculated by:

$$P_f=\frac{1}{N}{\displaystyle \sum _{j=1}^{N}E\left(I_F\left(x_j\right)\right)}\approx \frac{1}{N}{\displaystyle \sum _{j=1}^N{I}_F\left(x_j\right)}$$

(9)

4.2. Sensitivity Analysis Based on Kriging

To ensure the accuracy of the reliability sensitivity analysis for RPC-FST based on MCS, numerous repeatability calculations should be completed. However, the limit state function Equation (6) is highly complex and implicit. Therefore, the reliability sensitivity analysis for RPC-FST is limited by its high computational cost. To improve the calculation efficiency of the limit state function, the surrogate model Kriging is adopted to establish an approximation form for G(x).

The surrogate Kriging model was created from an initial number of computational RPC-FST dynamic response calculations in a design of experiments (DoE) sampling plan. There exist several DoE techniques, and in this work, the Latin hypercube sampling (LHS) plan was used. The LHS plan divides each design parameter into N equally sized intervals, where the same value of a parameter can only occur once [39]. TNT weight, TNT distance, steel tube thickness, RPC-FST diameter, and RPC-FST height are five design parameters. An example of a randomly generated LHS plan can be seen in Figure 7.

The basic idea of the Kriging method includes two parts: a regression model used to represent the global trend and a stochastic process used to represent the local behavior [40]. Therefore, the Kriging model can be expressed as Equation (10).

$$g_k(X)={\displaystyle \sum _{i=1}^p{f}_i(X){\beta }_i+z(X)}$$

(10)

where f(X) = {f₁(X), f₂(X), …f_p(X)}^T is the vector of regression basis functions; β = {β₁, β₂, …β_p}^T is the vector of the regression coefficients; p is the number of the basic functions in the regression model; and z(X) is assumed to be a Gaussian stationary process with a mean of zero and a standard deviation of σ. The covariance matrix can be expressed as:

$$\mathrm{c}\mathrm{o}\mathrm{v}[z(X^{(i)}),z(X^{(j)})]={\sigma }^2[R(X^{(i)},X^{(j)},\theta )]$$

(11)

where R(X⁽ⁱ⁾, X^(j)) is the spatial auto-correlation function between input samples X⁽ⁱ⁾ and X^(j), which coordinates the smoothness of the Gaussian model. It can be expressed as:

$$R[X^{(i)},X^{(j)},\theta ]=\exp (-{{\displaystyle \sum _{k=1}^m{\theta }_k\left|X_k^{(i)}-X_k^{(j)}\right|}}^2)$$

(12)

where θ = {θ₁, θ₂ … θ_k}^T is a hyperparameters vector defining the auto-correlation function. The hyperparameters vector θ can be estimated as follows using the maximum likelihood method:

$$\max F(\theta )=-\frac{m\ln ({\widehat{\sigma }}^2)+\ln (R)}{2}\hspace{1em}{\theta }_k\ge 0$$

(13)

where

$$\left\{\begin{array}{l}\begin{array}{l}{\widehat{\sigma }}^2={(g-F_k\widehat{\beta })}^T{R}^{-1}(g-F_k\widehat{\beta })/m\hfill \\ \widehat{\beta }={(F_k{}^T{R}^{-1}{F}_k)}^{-1}{F}_k{}^T{R}^{-1}g\hfill \end{array}\hfill \end{array}\right.$$

(14)

where ${\widehat{\sigma }}^2$ is an estimate of σ; $\widehat{\beta }$ is an estimate of β; F_k is a matrix, which gathers the regression functions based on training points; and g is the response vector corresponding to training points. The predicted response vector of the Kriging model at a given unknown point X can be expressed as:

$$g_k(X)=f^T(X)\widehat{\beta }+r^T(X)R^{-1}(g-F\widehat{\beta })$$

(15)

where r(X) is the correlation vector between the training and predicting points. It can be expressed as:

$$r^T(X)=\left\{R(X,X^{(1)}),R(X,X^{(2)})\dots R(X,X^{(m)})\right\}$$

(16)

The predictive mean and variance for the Kriging model can be performed as follows:

$$\begin{array}{l}{\mu }_{gk}(X)=f^T(X)\widehat{\beta }+r^T(X)R^{-1}(g-F\widehat{\beta })\hfill \\ {\sigma }^2{}_{gk}(X)={\widehat{\sigma }}^2\left\{1-r^T(X)R^{-1}r(X)+{[F^T{R}^{-1}r(X)-f(X)]}^T{(F^T{R}^{-1}F)}^{-1}[F^T{R}^{-1}r(X)-f(X)]\right\}\hfill \end{array}$$

(17)

To make the output Kriging surrogate model more accurately predict the displacement of the RPC-FST column, the relative error-index β_R (Equation (18)) is used to evaluate the model [41].

$${\beta }_R=\frac{g_k(X_L)-g_L(X_L)}{g_L(X_L)}$$

(18)

where X_L is a variable matrix composed of N_L samples generated by the LHS plan, g_k(X_L) is the output by the updated Kriging surrogate model, and g_L(X_L) is the output by the theory of FSI analysis.

Based on the MCS method, the reliability can be expressed as:

$$\frac{\partial P_f}{\partial {\theta }_{X_i}^{(k)}}={\displaystyle \int \cdots {\displaystyle {\int }_{R^n}{I}_F(x)\frac{\partial f_X(x)}{\partial {\theta }_{X_i}^{(k)}}}}\frac{1}{f_X(x)}{f}_X(x)dx=E[\frac{I_F(x)}{f_x(x)}\frac{\partial f_X(x)}{\partial {\theta }_{X_i}^{(k)}}]$$

(19)

The reliability sensitivity is obtained as Equation (20).

$${\beta }_S=\frac{\partial {\widehat{P}}_f}{\partial {\theta }_{X_i}^{(k)}}=\frac{1}{N_{MC}}{\displaystyle \sum _{j=1}^{N_{MC}}\frac{I_F(X)}{f_X(X)}\frac{\partial f_X(X)}{\partial {\theta }_{X_i}^{(k)}}}$$

(20)

The sensitivity coefficient is calculated as follows:

$$\mathrm{c}\mathrm{o}\mathrm{v}({\beta }_S)\approx \sqrt{\frac{\frac{1}{N_{MC}-1}\left(\frac{1}{N_{MC}}{\displaystyle \sum _{j=1}^{N_{MC}}{[\frac{I_{{}_F}(X)}{f_X(X)}\frac{\partial f_X(X)}{\partial {\theta }_{X_i}^{(k)}}]}^2-{({\beta }_S)}^2}\right)}{\left|{\beta }_S\right|}}$$

(21)

The sensitivity of the five input variables of steel tube thickness, RPC-FST column height, TNT weight, TNT distance, RPC-FST section diameter to the displacement error risk of RPC-FST columns, and sensitivity coefficient are evaluated by Equations (20) and (21). Due to the units of five random variables not being unified and the corresponding orders of magnitude of the five input variables being different, the sensitivity standardization calculation is required. After standardized distribution parameters μ_Xi and σ_Xi, the corresponding sensitivity can be expressed as Equation (22) [42,43].

$$\begin{array}{l}{S}_{{\mu }_{X_i}}=\frac{\partial P_f}{\partial {\mu }_{X_i}}\frac{{\sigma }_{X_i}}{P_f}\\ S_{{\sigma }_{X_i}}=\frac{\partial P_f}{\partial {\sigma }_{X_i}}\frac{{\sigma }_{X_i}}{P_f}\end{array}$$

(22)

The reliability sensitivity gradient of the i-th random variable to P_f is:

$$S_i=\sqrt{S_{{\mu }_{X_i}}^2+S_{{\sigma }_{X_i}}^2}$$

(23)

The sensitivity factor of the RPC-FST system can be obtained by standardizing the gradient of five parameters, as given in Equation (24).

$${\lambda }_i=\frac{S_i}{{\displaystyle \sum _{k=1}^{N_{MC}}{S}_k}}\times 100\%$$

(24)

where N_MC is the number of sample pools.

4.3. Sensitivity Analysis of the RPC-FST

Based on a numerical example provided in Section 3, 20 sample groups are selected to establish the Kriging surrogate model. The accuracy of the Kriging surrogate model needs to be verified, so 10 groups of data obtained by numerical calculation are used as verification groups and standard values. Absolute error E_A and relative error-index β_R are critical parameters for evaluating the response surface. The design of the experiment is shown in

Table 9. Design of experiment.

Variable	Lower Bound	Upper Bound	Unit
TNT weight W	150	300	kg
Distance of TNT L	1000	2000	mm
Diameter of RPC-FST D	300	500	mm
Thickness of steel tube t	10	20	mm
Height of RPC-FST H	2000	3000	mm

. The Bayesian method can update the coefficient of variation based on samples. It should be pointed out that coefficients of variation can be obtained in large-scale tests [44]. The distribution of random variables is shown in

Table 10. Distribution of random variables.

Variable	Distribution	Mean	Variable Coefficient	Unit
W	normal	225	0.01	kg
L	normal	1500	0.001	mm
D	normal	400	0.01	mm
t	normal	15	0.1	mm
H	normal	2500	0.001	mm

. These numerical examples are tested on a personal computer with Intel i7-9700 and 16 GB memory. The analysis structure of the relative error-index is shown in

Table 11. Testing results based on testing samples.

N	W (kg)	L (m)	D (mm)	T (mm)	H (mm)	gL(XL) (mm)	gk(XL) (mm)	EA	βR
1	278.98	1.57	318.33	18.90	2516.67	16.00	14.91	−1.09	−0.07
2	149.32	1.82	338.67	12.50	2333.67	18.40	17.99	−0.41	−0.02
3	154.46	1.85	446.75	15.67	2798.33	2.73	3.00	0.26	0.10
4	245.45	1.17	371.67	15.70	2450.00	21.01	19.43	−1.58	−0.08
5	255.15	1.50	335.00	15.30	2216.67	18.72	19.82	1.10	0.06
6	224.02	1.23	461.67	17.90	2950.00	7.16	6.76	−0.40	−0.06
7	175.75	1.37	358.33	14.90	2416.67	18.16	17.74	−0.41	−0.02
8	166.05	1.10	421.67	18.10	2016.67	10.93	9.47	−1.46	−0.13
9	236.95	1.90	368.33	19.10	2650.00	5.68	6.50	0.82	0.15
10	240.18	1.17	495.00	14.30	2183.33	20.19	19.21	−0.98	−0.05

. The data shows that the error fluctuates within an acceptable range. It confirmed that the derived Kriging surrogate model could predict the displacement of the RPC-FST accurately. The Kriging surrogate model constructed by the limit state function is used for calculation. Furthermore, 10⁶ MCS is used to calculate the reliability sensitivity prognosis of the RPC-FST system [45].

From Figure 8a, it is obvious that the overall maximum displacement of RPC-FST has a greater correlation with the five variables of steel tube thickness, RPC-FST column height, TNT weight, TNT distance, and RPC-FST section diameter. However, the positive and negative correlations of the five variables with RPC-FST maximum displacement are different. The maximum displacement risk of RPC-FST increases with the increase of TNT weight and RPC-FST column height (positive correlation) and increases with the decrease of RPC-FST section diameter, steel tube thickness, and TNT distance (negative correlation). This can be concluded from Figure 8b.

To compare the efficiency of the three models (the modified HJC model, the modified K&C model, and the Kriging surrogate model) in calculating the maximum displacement, we selected four sets of data for analysis, as shown in Table 11. The results show that the updated Kriging surrogate model in this paper can significantly improve the sensitivity analysis efficiency of the RPC-FST system under explosion load, as shown in

Table 12. Efficiency competitions of three models.

Testing Sample	Total CPU Time (s)
	Modified HJC Model	Modified K&C Model	Kriging Surrogate Model
1	9562	15,423	67
2	15,447	19,661	70
3	11,212	15,223	62
4	11,329	18,195	69

.

5. Conclusions

The displacement prediction model of RPC-FST under random explosion load is established based on fluid–solid interaction (FSI) theory. The main conclusions are as follows:

(1) The modified K&C model has higher accuracy than the modified HJC model in calculating the dynamic response of RPC-FST columns under explosion load. This is mainly because the previous HJC model is applicable to RPC materials with coarse aggregate and porosity of 10.5%. Although the HJC model is calibrated by modifying the porosity and related parameters, a lack of test for the ultimate strain of RPC (EFMIN) caused the calculated deflection to be lower than the test deflection in the modified HJC model.

(2) It takes a long time to use the ALE method to carry out the numerical calculations under RPC-FST explosion loading. The updated Kriging surrogate model can significantly shorten the required time for the FSI theoretical analysis by more than 200 times. Besides, the RPC-FST dynamic response analysis based on the updated Kriging surrogate model has high accuracy.

(3) Based on the analysis of the five random parameters by the Kriging surrogate model, the TNT weight has the most significant impact on the risk of maximum displacement of RPC-FST compared with other parameters.