An investigation of ballistic response of reinforced and sandwich concrete panels using computational techniques

Abstract

Structural performance of nuclear containment structures and power plant facilities is of critical importance for public safety. The performance of concrete in a high-speed hard projectile impact is a complex problem due to a combination of multiple failure modes including brittle tensile fracture, crushing, and spalling. In this study, reinforced concrete (RC) and steel-concrete-steel sandwich (SCSS) panels are investigated under high-speed hard projectile impact. Two modeling techniques, smoothed particle hydrodynamics (SPH) and conventional finite element (FE) analysis with element erosion are used. Penetration depth and global deformation are compared between doubly RC and SCSS panels in order to identify the advantages of the presence of steel plates over the reinforcement layers. A parametric analysis of the front and rear plate thicknesses of the SCSS configuration showed that the SCSS panel with a thick front plate has the best performance in controlling the hard projectile. While a thick rear plate is effective in the case of a large and soft projectile as the plate reduces the rear deformation. The effects of the impact angle and impact velocity are also considered. It was observed that the impact angle for the flat nose missile is critical and the front steel plate is effective in minimizing penetration depth.

Keywords

concrete panels / projectile impact / finite element modeling / smoothed particle hydrodynamics / strain rate effect

Introduction

As of 2018, there are 99 operating nuclear power plants and spent fuel storage facilities in 30 states of the United States of America. The potential consequences of beyond-design events are of critical importance as was seen in the unfortunate incidents of the Fukushima Power Plant in Japan and the World Trade Center in New York. Steel-concrete-steel sandwich (SCSS) and reinforced concrete (RC) panels have been widely used in nuclear power plant facilities as prefabricated elements for modular construction [¹,²]. Various nuclear-related structures including containment buildings and dry storage systems are built using SCSS panels due to the low cost of fabrication and high levels of safety [³,⁴]. To fabricate SCSS panels, plain concrete is placed between two steel plates. One of the main advantages of the SCSS panels over conventional RC counterparts is the lack of scabbing at the protected faces of the panel in the case of a projectile impact. This phenomenon happens due to wave propagation through the concrete even if the projectile has not penetrated completely [⁵]. The application of steel plates or carbon fiber sheets on the protected face minimizes the risk of fragmentation. Furthermore, the increased ductility of the panels with these materials absorbs the kinetic energy of a projectile. The main failure modes in the case of hard projectile impact, which has limited deformation or considered rigid, are penetration, scabbing, and spalling. A high-speed hard projectile might penetrate through the containment building and reactor walls or the overpack of a dry storage system and cause radioactive leakage. The impact may also cause a tip-over incident [⁶,⁷]. This study focuses on the structural performance of SCSS and RC panels under hard projectile where scabbing and perforation is the main concern. A well-performing panel is expected to have a minimum intrusion of the projectile with a minimum amount of scabbing.

Studies on the protection of nuclear power plant facilities from projectile impact started in the 1950s. Several full-scale experiments were performed with different types of projectiles including fighter jets and aircraft engines [⁸]. Hanifehzadeh et al. [⁹] studied aircraft engine impact on a composite storage container. It was found that SCSS configuration is more effective than RC with a single reforcing steel layer in controlling the penetration of a soft projectile while the rear plate thickness controls the global deformation. Sohel and Liew [¹⁰] studied sandwich panels made with lightweight concrete impacted by a low-velocity projectile through an experimental program. A novel J-hook was proposed for efficient composite behavior. It was observed that one-fourteenth the deflection of the span length causes a 30% reduction in flexural strength. Abdel-Kader and Fouda [¹¹] studied the influence of reinforcement ratio by testing 500 mm × 500 mm × 100 mm panels. A hard projectile with 23 mm diameter and a velocity from 200 to 430 m/s was used. It was found that the presence of a rear steel plate improves the perforation resistance more effectively than a front steel plate.

Bruhl et al. [¹²] proposed a single degree-of-freedom model to calculate the maximum displacement of a sandwich panel subjected to missile impact. The load was modeled as a concentrated load at the center of the panel, which is monotonically increased with a loading diameter-to-span length ratio of 0.25. Conventional finite element (FE) modeling was used to calibrate the model and it was concluded that under the fixed boundary conditions, a plastic collapse causes radial yield lines and concrete cracking, while in the simply supported case, a diagonal cracking pattern is observed. It was also indicated that unlike the yield strength of the steel and the reinforcement ratio, the compressive strength of concrete does not have a major impact on the structural response of the panels. Heckötter and Vepsä [¹³] performed large-scale experimental and numerical analysis on RC slabs impacted by a hard and a soft missile and reported the global and local failure modes. The conventional FE analysis and smoothed particle hydrodynamics (SPH) in predicting the structural response were compared. The results showed that both methods could simulate the impact problem with acceptable accuracy; however, both methods underestimated the residual missile velocity. For impact velocities less than 100 m/s, the predicted results were close; while, for higher velocities, the SPH results were closer to the experimental data.

Different computational methods for modeling of fracture is available in the literature including FEs, extended FEs and meshless methods. There are several issues associated with the mesh-based analysis including loss of accuracy under large deformations, lack of arbitrary crack growth paths, and discontinuity. Mesh generation for a complicated geometry is another limitation of grid-based methods. Under extreme deformation of elements, erosion algorithms could partially solve the problem. However, with element erosion, mass conservation principle is violated. Furthermore, the erosion criterion is solely based on experimental data. Several mesh-free methods have been developed recently to model fracture propagation in two and three dimensions [¹⁴–¹⁷]. The main concept of the mesh-free methods is to provide a numerical solution for the integral or partial differential equations with different boundary conditions for a set of arbitrarily distributed nodes or particles. SPH, which was introduced by Lucy [¹⁸] and extended by Gingold and Monaghan [¹⁹], is currently known as one of the efficient techniques to solve applied mechanics problems. SPH involves a set of particles, which have defined material properties and interact with each other within the range obtained from weight or smoothing functions [²⁰]. An approximation is performed to obtain the value of a parameter in an arbitrary particle, by summing the contributions from a set of neighboring particles. The radius of the neighboring area is referred to as the smoothing length. Since the connectivity between particles is calculated and the connections can change during the analysis, large deformations could easily be handled. It is assumed that the nodes are distributed uniformly in the volume. If the mass and density of a specific part are known, the weight and volume of the particles are obtained by dividing the total volume and weight by the number of particles.

The cracking particles method (CPM) was introduced by Rabczuk and Belytschko [²¹] where a discontinuity along a line in 2D or a plane in 3D is introduced at each crack lacation. The cracking criterion is checked independently at each particle and the collection of the cracked particles is approximated as the crack path. The advantage of the method is its computational efficiency and no existance of spurious instabilities when used with Lagrangian kernels. Several studies have been published recently based on the element-free Galerkin (EFG) method for dynamic fracture modeling including cohesive cracks and fragmentation. Rabczuk and Zi [²²] presented the extended EFG method for dynamic fracture with cohesive cracks and to handle the presence of multiple cracks. In EFG, the jump in the displacement field of the particles defines the cracks. Rabczuk and Eibl [²³] used the Lagrangian mesh-free method to model the dynamic failure and fragmentation of concrete under blast and impact load based on the model by Randles and Libersky [²⁴]. A viscous damage model was proposed that account for the strain rate effects. The method was verified, and the crack pattern was compared with experimental data. Alternative methods such as peridynamics to model dynamic fracture and impact loading have been developed [²⁵]. The fracture behavior is defined based on the breakage of the bonds between particles.

In this study, the structural behavior of two different panels: RC and SCSS against high-speed hard projectile impact is investigated. The penetration depth and rear deformation are considered for comparison of the structural behavior at different impact angles and velocities. The penetration depth is important in the case of a high-speed sharp projectile impact; while, the rear deformation is an indicator of the global energy absorption of the system. The material models for steel and concrete are initially calibrated and verified with the experimental data available in the literature. Four well-known empirical equations are also used for comparison purposes. The strain rate effect for steel and concrete are included in the material constitutive formulations for a more accurate prediction of the experimental results. First, the conventional FE and SPH techniques are compared and the advantages and disadvantages of each model are discussed. The most effective configuration for RC and SCSS is identified. Next, a parametric study was performed to determine the optimal plate thicknesses for SCSS panels. It was reported in previous studies that the reinforcement ratio has no significant effect on the performance of RC panels [²⁶,²⁷]. Therefore, here, only one model is developed for the RC panel and the parametric study was performed only on the SCSS panel. The main difference between this study and previous work is the modeling of highly confined conditions of concrete provided by the thick steel plate and implementation of the strain rate effects for concrete simultaneously. The results presented here advance the understanding of the RC and SCSS panel behavior under ballistic loads and may be used in the design of protective structures in the future.

Modeling

Abaqus/Explicit [²⁸] is used for the dynamic analysis in this paper. The explicit dynamic integration method (also known as the forward Euler or central difference algorithm) has the ability to analyze problems involving high-speed, large nonlinearity, post-buckling, and collapse [²⁹]. Another advantage of Abaqus/Explicit is that it requires much less disk space and memory than Abaqus/Standard for the same simulation. In the explicit analysis at the end of each time step, the system matrices are updated, and the new system of equations is solved without iteration. If the increments are small enough, accurate results are obtained, otherwise, the solution diverges. The explicit algorithm uses the central difference method to integrate the equation of motion explicitly through time. During each increment, the initial kinematic conditions are used to calculate the kinematic conditions for the next step. The analysis is conducted through many small load increments according to

(1){{\displaystyle \mathit{M}\ddot{\mathit{u}}|}}_t={{\displaystyle \left(\mathit{P}-\mathit{I}\right)|}}_t,

where \mathit{u}, \dot{\mathit{u}}, and \ddot{\mathit{u}} are displacement, velocity, and acceleration vectors, respectively, M is the mass matrix, P is the applied load vector, and I is the internal force vector. Using a diagonal mass matrix for efficiency in the explicit procedure, the nodal accelerations can be obtained at any given time, t, using

(2){{\displaystyle \ddot{\mathit{u}}|}}_t={\mathit{M}}^{-\mathrm{1}}·{{\displaystyle \left(\mathit{P}-\mathit{I}\right)|}}_t.

The velocity and displacement are obtained using the central difference integration rule

(3){{\displaystyle \dot{\mathit{u}}|}}_{t+\frac{\mathrm{\Delta }t}{\mathrm{2}}}={{\displaystyle \dot{\mathit{u}}|}}_{t-\frac{\mathrm{\Delta }t}{\mathrm{2}}}+\left({\displaystyle \frac{{{\displaystyle \mathrm{\Delta }t|}}_{t+\mathrm{\Delta }t}+{{\displaystyle \mathrm{\Delta }t|}}_t}{\mathrm{2}}}\right){{\displaystyle \ddot{\mathit{u}}|}}_t,

(4){{\displaystyle \mathit{u}|}}_{t+\mathrm{\Delta }t}={{\displaystyle \mathit{u}|}}_t+{{\displaystyle \mathrm{\Delta }t|}}_{t+\mathrm{\Delta }t}·{{\displaystyle \dot{\mathit{u}}|}}_{t+\frac{\mathrm{\Delta }t}{\mathrm{2}}},

where Δt is the time step. Unlike the implicit solution, no iteration for internal force equilibrium is required for displacement, velocity and acceleration in the explicit solution. The central difference integration rule is conditionally stable, that is, the solutions become unstable and diverge rapidly if the time increment is large. The main disadvantage of explicit analysis is the stability criterion which usually requires very small time-steps. The time step is based on the smallest transit time required to cross the smallest element in the model by a dilatational wave. Therefore, both accuracy and computational efficiency should be considered in defining the element size. The time increment must satisfy

(5)\mathrm{\Delta }t\le {\displaystyle \frac{\mathrm{2}}{{{\displaystyle \omega }}_{\max }}}\left(\sqrt{\mathrm{1}+{{{\displaystyle \xi }}_{\max }}^{\mathrm{2}}}-{{\displaystyle \xi }}_{\max }\right),

where {{\displaystyle \omega }}_{\max } is the element maximum eigenvalue and {{\displaystyle \xi }}_{\max } is the fraction of critical damping in the mode with the highest frequency. By default, Abaqus/Explicit considers bulk viscosity damping of 1.2 to minimize high frequency oscillations. An estimate of stable time increment is given by

(6)\mathrm{\Delta }t=\min \left({\displaystyle \frac{L^{\mathrm{e}}}{{{\displaystyle c}}_{\mathrm{d}}}}\right),

where L^e is the characteristic length and c_d is the dilation wave speed of the material. For further details regarding the analysis method, one is referred to the user manual of the software [²⁹]. The geometry, boundary, and loading conditions, meshing, the nonlinear behavior of materials and the consideration of the strain rate effects are provided below.

Geometry and FE model

The geometry for the projectile impact studies here is selected from a previously published experimental study in Ref. [³⁰] such that the computer models may be validated. Lee et al. [³⁰] used a 50 kg steel cylinder with a yield strength of 950 MPa and ultimate tensile strength of 1.1 GPa as a projectile impacting a 2 m× 2 m SCSS panel as shown in Fig. 1(a). The projectile had 155 mm diameter with an impact velocity of 314 m/s. As shown in Fig. 1(b), in addition to the SCSS configuration in Ref. [³⁰], an RC panel is also considered in this study.

Considering the symmetry in the geometry, only one-quarter of the panels was modeled. The dimensions of the parts and the concrete strength in the experiment are provided in Table 1. The type and the total number of the elements used for each component in modeling are summarized in Table 2. C3D8R, which is an 8-noded brick element, was used in the Lagrangian formulation for all parts except for the reinforcement of the RC panels. The reinforcement was modeled with T3D2, i.e., 2-noded truss elements. The truss elements were coupled with brick elements of the panel using the “EMBEDDED ELEMENTS” function of Abaqus [²⁸] in which the translational degrees of freedom of the truss nodes are constrained to the neighboring nodes of the host brick elements [²⁹]. By partitioning the concrete core, a finer mesh was created for the impact area to achieve more accurate results (see Fig. 2). Since the element conversion in the SPH simulations is computationally costly, element conversion was only assigned to the impact area.

To capture the perforation of the missile, an element deletion function is necessary. A user-subroutine was implemented to delete the failed elements with respect to equivalent plastic strain in compression, which is defined as

(7){{\displaystyle \overline{\mathit{ϵ}}}}^{\mathrm{p}\mathrm{l}}={{\displaystyle \int }}_{\mathrm{0}}^t\sqrt{{\displaystyle \frac{\mathrm{2}}{\mathrm{3}}}{\dot{\mathit{ϵ}}}^{\mathrm{p}\mathrm{l}}:{\dot{\mathit{ϵ}}}^{\mathrm{p}\mathrm{l}}\mathrm{d}t},

where {\dot{\mathit{ϵ}}}^{\mathrm{p}\mathrm{l}} is the rate of plastic flow. The penetration depth is highly sensitive to the element erosion limit in compression. Therefore, multiple analysis with different limits were performed to obtain the proper penetration depth to match the experimental data. A value of 0.225 for the equivalent compressive plastic strain was found to provide acceptable results.

The front and rear plates were tied to the concrete, which means that the displacement in all three directions is identical between the nodes of the two parts on the shared surface. Due to the higher rigidity of the concrete core, the font plate has a higher plastic deformation during the impact. Therefore, a finer mesh was defined for the plate, 5 mm versus 20 mm, as compared to the concrete. The hard contact was used between the projectile and the front plate in the normal direction, and friction with allowed separation was used after contact. Hard contact means that no intrusion to the slave surface is allowed for the master nodes. Penalty method in tangential direction with a kinetic friction coefficient of 0.6, for steel-to-steel contact, was applied.

Constitutive models with strain rate effect

Concrete

The strain rate effect for concrete becomes important in the case of blast and impact loading. It has been frequently reported that a noticeable increase in strength, called dynamic strength, is observed when concrete is subjected to a high strain rate [³¹–³³]. The main reasons for this increase in concrete strength are the viscosity, inertia force, and the lateral confinement [³⁴]. To consider this phenomenon, a dynamic increase factor (DIF), defined as the ratio of the strength under high strain rate load to the quasi-static strength in uniaxial loading, is applied here to compressive and tensile strength. The proposed empirical equations are normally a function of strain rate only. However, the actual material behavior is a result of various interacting factors including concrete quality, aggregate size, curing and moisture conditions, age, and confinement and boundary conditions [³⁵]. Therefore, a validation of the model against experimental data are essential. The equations recommended for DIF are based on experimental data obtained from drop-hammer, servo-hydraulic, and split Hopkinson tests [³⁵]. The well-known CEB-FIP model [³⁶] is used in this study, which defines DIF in compression as

(8)

(9)DIF=\gamma {{\displaystyle \left({\displaystyle \frac{{{\displaystyle \dot{ϵ}}}_{\mathrm{d}}}{{{\displaystyle \dot{ϵ}}}_{\mathrm{s}}}}\right)}}^{\mathrm{1}/\mathrm{3}},\hspace{0.17em}\mathrm{f}\mathrm{o}\mathrm{r}\hspace{0.17em}\dot{ϵ}>{\mathrm{30}\hspace{0.17em}\mathrm{s}}^{-\mathrm{1}},

(10)\log \gamma =\mathrm{6.156}a-\mathrm{2},

(11)\log \gamma =\mathrm{6.156}a-\mathrm{2},

where {{\displaystyle \dot{ϵ}}}_{\mathrm{d}}, {{\displaystyle \dot{ϵ}}}_{\mathrm{s}} are the dynamic and quasi-static strain rate, respectively, the latter of which is taken as 30 × 10⁻⁶ s⁻¹. It is important to note that in Eqs. (1)–(10), the strain components are the equivalent strain measures. The model is a bilinear relation where the slope changes at 30 s⁻¹. Similarly, for tension the DIF is given as

(12)DIF={{\displaystyle \left({\displaystyle \frac{{{\displaystyle \dot{ϵ}}}_{\mathrm{d}}}{{{\displaystyle \dot{ϵ}}}_{\mathrm{s}}}}\right)}}^{\mathrm{1.016}\delta },\hspace{0.17em}\mathrm{f}\mathrm{o}\mathrm{r}\hspace{0.17em}\dot{ϵ}\le {\mathrm{30}\hspace{0.17em}\mathrm{s}}^{-\mathrm{1}},

(13)DIF=\beta {{\displaystyle \left({\displaystyle \frac{{{\displaystyle \dot{ϵ}}}_{\mathrm{d}}}{{{\displaystyle \dot{ϵ}}}_{\mathrm{s}}}}\right)}}^{\mathrm{1}/\mathrm{3}},\hspace{0.17em}\mathrm{f}\mathrm{o}\mathrm{r}\hspace{0.17em}\dot{ϵ}>{\mathrm{30}\hspace{0.17em}\mathrm{s}}^{-\mathrm{1}},

(14)\log \beta =\mathrm{7.11}\delta -\mathrm{2.33},

(15)\delta =\mathrm{1}/\left(\mathrm{10}+{\displaystyle \frac{{{\displaystyle \mathrm{6}f}}_{\mathrm{c}}^{\text{'}}}{{{\displaystyle f}}_{\mathrm{0}}^{\text{'}}}}\right),

where f′₀ is the reference strength (taken as 10 MPa), and f′_c is the static compressive strength of concrete. Equations (8), (9), (12), and (13) are applicable to strain rates from 3×10⁻⁶ to 300 s⁻¹. Unified Facilities Criteria (UFC) 3-340-02 [³⁷] published by the US Department of Defense provids guidelines for estimating the average strain rate, {{\displaystyle \dot{ϵ}}}_{\mathrm{c}}, in the dynamic loading scenarios as

(16){{\displaystyle \dot{ϵ}}}_{\mathrm{c}}={\displaystyle \frac{{{\displaystyle ϵ}}_{\mathrm{c}\mathrm{0}}}{{{\displaystyle t}}_{\mathrm{e}}}},

where {{\displaystyle ϵ}}_{\mathrm{c}\mathrm{0}} is the strain at the maximum stress level for concrete (taken equal to 0.002) and t_e is the time to reach yielding. From the FE model, t_e was estimated as 6.1 ×10⁻⁴ s, which results in a strain rate of 3.3 s⁻¹. Consequently, the DIF for concrete in compression and tension were calculated as 1.48 and 1.69, respectively, for 28 MPa concrete.

The Concrete Damage Plasticity (CDP) model available in Abaqus [²⁸] was used for modeling the concrete. CDP is based on the concept of isotropic damage elasticity in conjunction with isotropic tensile and compressive plasticity with the inclusion of strain hardening in compression. The model considers the material degradation under load reversal, which is required in a wave propagation simulation. However, there are no criteria defined for the removal of elements in tensile cracking or compressive crushing, and consequently scabbing and spalling could not be captured by the conventional FE approach. The damage parameter is assigned separately for tension and compression as a scalar value. The parameters used in the CDP model are given in Table 3. The default values in Abaqus [²⁹], which are obtained from experimental data in Refs. [³⁸–⁴¹] were used in the model. A detailed description of the model and the related parameters are available in Jankowiak and Lodygowsky [⁴²] and Abaqus User Manual [²⁹].

The same material properties as in the experiments performed by Lee et al. [³⁰] were used for concrete as shown in Table 4. The Hognestad parabola [⁴³], given by

(17){{\displaystyle \sigma }}_{\mathrm{c}}={{\displaystyle f}}_{\mathrm{c}}^{\text{'}}\left[\mathrm{2}\left({\displaystyle \frac{{{\displaystyle ϵ}}_{\mathrm{c}}}{{{\displaystyle ϵ}}_{\mathrm{c}\mathrm{0}}}}\right)-{{\displaystyle \left({\displaystyle \frac{{{\displaystyle ϵ}}_{\mathrm{c}}}{{{\displaystyle ϵ}}_{\mathrm{c}\mathrm{0}}}}\right)}}^{\mathrm{2}}\right],

where σ_c and ε_c are the compressive stress and strain, respectively, are used to represent the concrete stress-strain behavior in compression. In tension, up to tensile strength, a linear behavior was assumed, and the softening part of the curve was considered using the model proposed by Wang and Hsu [⁴⁴] as

(18){{\displaystyle \sigma }}_{\mathrm{t}}={{\displaystyle E}}_{\mathrm{c}}{{\displaystyle ϵ}}_{\mathrm{t}},\hspace{0.17em}\mathrm{i}\mathrm{f}\hspace{0.17em}{{\displaystyle ϵ}}_{\mathrm{t}}\le {{\displaystyle ϵ}}_{\mathrm{c}\mathrm{r}},

(19){{\displaystyle \sigma }}_{\mathrm{t}}={{\displaystyle f}}_{\mathrm{c}\mathrm{r}}{{\displaystyle \left({\displaystyle \frac{{{\displaystyle ϵ}}_{\mathrm{c}\mathrm{r}}}{{{\displaystyle ϵ}}_{\mathrm{t}}}}\right)}}^{\mathrm{0.4}},\hspace{0.17em}\mathrm{i}\mathrm{f}\hspace{0.17em}{{\displaystyle ϵ}}_{\mathrm{t}}>{{\displaystyle ϵ}}_{\mathrm{c}\mathrm{r}},

where σ_t and ε_t are the tensile stress and strain at any given point, respectively, E_c is the modulus of elasticity of concrete, f_cr is the cracking stress and ε_cr is the cracking strain of concrete. Equations (17), (18), and (19) are plotted in Fig. 3 for 23 MPa concrete for static and dynamic conditions.

Additional concrete properties used in the FE model are presented in Table 4. The tensile strength f_t, and the modulus of elasticity E_c, were calculated based on the compressive strength according to the equations in Table 4. The methodology for calibration of the model based on standard material tests is provided in Hanifehzadeh et al. [⁴¹].

Steel

The yield and ultimate strength of steel are also sensitive to the strain rate. For mild steel, the yield stress and strain as well as the ultimate stress increase with increasing strain rate; however, the modulus of elasticity remains more or less constant [³¹,⁴⁵]. To consider the effect of strain rate on the mechanical properties of steel, Cowper and Symonds [⁴⁶] equation that defined DIF as

(20)DIF=\mathrm{1.0}+{{\displaystyle \left({\displaystyle \frac{\dot{ϵ}}{C}}\right)}}^{\frac{\mathrm{1}}{q}}

where \dot{ϵ} is the strain rate, C and q are constants assumed as 40.4 and 5 for mild steel. Malvar and Crawford [⁴⁷] performed a comprehensive experimental study on the strain rate effect on steel reinforcement. The following equation was proposed

(21)DIF={{\displaystyle \left({\displaystyle \frac{\dot{ϵ}}{{\mathrm{10}}^{-\mathrm{4}}}}\right)}}^{{{\displaystyle \alpha }}_{\mathrm{s}}},

where for the yield stress, α_s = α_fy is found as

(22){{\displaystyle a}}_{\mathrm{s}}=\mathrm{0.074}-\mathrm{0.040}{\displaystyle \frac{{{\displaystyle f}}_{\mathrm{y}}}{\mathrm{60}}},

where f_y is the yield strength of steel (f_y is in ksi and 1 ksi= 6.895 MPa). Equation (21) is valid for a yield stress range from 290 to 710 MPa and a strain range from 10⁻⁴ to 225 s⁻¹. The UFC 3-340-02 [³⁷] recommendation for calculation of the average strain rate for steel, {{\displaystyle \dot{ϵ}}}_{\mathrm{s}}, is

(23){{\displaystyle \dot{ϵ}}}_{\mathrm{s}}={\displaystyle \frac{{{\displaystyle f}}_{\mathrm{d}\mathrm{y}}}{{{\displaystyle E}}_{\mathrm{s}}{{\displaystyle t}}_{\mathrm{e}}}},

where f_dy is the dynamic strength and E_s is the modulus of elasticity of steel, and f_dy is obtained by static strength times the DIF. In this study, using 6.1×10⁻⁴ for t_e, DIF was calculated as 1.65 for the steel plates and 1.45 for the reinforcement. Note that in the Eq. (23), an initial value for DIF should be assumed to calculate {{\displaystyle \dot{ϵ}}}_{\mathrm{s}}. Then, the DIF value is calculated from Eq. (20) and compared with the initial assumption. If the value is not close, the new value of DIF is used as the next estimate and iterations are performed until the DIF value converges. A bilinear stress strain curve was assumed in this study and adjusted by multiplying with the DIF. Quasi-static material properties of the ASTM A36 [⁴⁸] and ASTM A615 [⁴⁹] steel were used for the plates and reinforcement, respectively, and the dynamic yield strength obtained according to the formulation presented above are summarized in Table 5 and shown in Fig. 4.

Analysis approach

Conventional FE modeling is widely used to study impact problems. Despite sufficient accuracy in most cases, element conversion problems are observed under very large deformations. To address this issue, a user subroutine, VUSDFLD, is a viable solution that can be used to define the criteria to remove elements with excessive deformations. There are other issues associated with element removal including conservation of mass and energy that need to be addressed. Considering a brittle panel, when the impact energy is applied on the free surface, high-velocity impact wave travels through the medium and reflects from the rear face. This traveling wave causes periodic compression and tension stresses in the material. For brittle materials including concrete and ceramics, low tensile strength is the main cause of the failure. Due to the impact wave, cracks parallel to the free surface form within the volume and fragmentation and spalling occurs on the rear face.

While the conventional FE formulation is unable to model the scabbing and fragmentation, the mesh-free SPH method is suitable for brittle materials experiencing high deformation and consequently fragmentation. The SPH model with automatic element conversion feature has recently become available in Abaqus\Explicit [²⁸], which can handle fragmentation problems. In the SPH method, initially, the conversion criteria based on stress or strain should be defined. The conversion of brick elements (C3D8R) to SPH particles occurs during the analysis when the predefined criteria are met. A strain-based criterion, which takes the absolute value of the maximum principal strain, was used here for the SPH model and one particle was generated per parent element. A threshold value of 0.3 was used for conversion based on the Ref. [⁵⁰]. The number of particles per parent element can be defined by the user; however, increasing the number of particles exponentially increases the cost of computation. When the elements are converted to particles, the particles behave like a free flying object similar to fragmentation in concrete. Unlike the element deletion method, the particles still provide resistance to projectile penetration, leading to a more realistic simulation. Further, the mass and the energy of the system are conserved. Both modeling approaches were investigated in this study and the results are compared in terms of penetration depth and global deformation.

The material strength and properties as described above (including the DIF) are exactly the same for both the conventional and the SPH models. To control the excessive element deformation in the FE model, elements were removed from the mesh when their equivalent compressive plastic strain reached 0.225 using VUSDFLD user subroutine. This value is close to the value, 0.2, which was reported by Rodríguez et al. [⁵¹].

Mesh sensitivity analysis and hourglassing control

The sensitivity of the impact simulations to the mesh size was studied for the reference case and the results are shown in Table 6. These results were used to determine the optimum mesh size. As shown in Table 6, three different mesh sizes were generated for the concrete, and the penetration depth and the rear deformation were calculated. A very fine mesh with 5 mm sides and four elements through the depth is considered for the front plate as this part of model is expected to have the highest element deformation due to impact. Since already a fine mesh was used for the front plate, the mesh sensitivity was only performed on the concrete. The minimum mesh size was used in a 350 mm × 350 mm area of the concrete, which is directly impacted by the projectile (see Fig. 2). This corresponds to an element size of 20 mm.

Under certain loading conditions, linear reduced-integration elements can experience a pattern of nonphysical deformation called hourglassing. Hourglass or zero energy mode happens when deformation occurs in the element with zero deformation energy. In other words, a lack of resistance against deformation occurs in the element. This issue is expected in models with coarse meshes. A built-in hourglass control algorithm known as the integral viscoelastic approach exists to prevent this problem without introducing excessive constraints on the element’s physical response. q being an hourglass mode magnitude and Q being the force (or moment) conjugate to q, the integral viscoelastic approach is defined as

(24)Q={{\displaystyle \int }}_{\mathrm{0}}^{t}sK\left(t-t^{\prime }\right){\displaystyle \frac{\mathrm{d}q}{{\mathrm{d}t}^{\prime }}}{\mathrm{d}t}^{\prime },

where K is the hourglass stiffness, and s is a dimensionless scaling factor. If the hourglassing is still observed in the model, the scaling factor can be increased by the user. However; using values larger than the default values can produce an excessively stiff response and sometimes result in instability. In addition to Abaqus [²⁸] built-in hourglass control, it has been recommended to refine the mesh to avoid hourglassing. A more quantitative approach is to look at the history of the artificial strain energy, which is primarily the energy dissipated to control hourglassing deformation. The artificial strain energy may be compared to the internal strain energy. It is necessary to keep the ratio below 5% to ensure that hourglassing does not occur and it does not have a significant impact on the results [²⁹]. Figure 5 shows this ratio for one of the cases analyzed later in this paper. As shown in the figure, the ratio remains below 3.5% during the analysis indicating that hourglassing is not a problem.

Validation of the FE and SPH models

The conventional FE and the SPH methods are validated for the SCSS configuration using the experimental data. The results in terms of penetration depth and rear deformation at the center are compared. The results of the conventional FE simulations (when the velocity of the projectile becomes zero) are presented in Fig. 6. In the case of a soft projectile, a flexural failure is dominant; however, here due to the large depth of the panel compared to the other two dimensions and the high impact velocity, the failure mode is punching shear rather than flexural failure. Scabbing is prevented by the rear plate. Note that the blunt projectile which has the lowest ability to penetrate compared to a spherical or conical nose is used in this study, as had been the case in the experiments reported by Lee et al. [³⁰].

The reaction force time history at the reaction column, see Fig. 2, is shown in Fig. 7. The maximum reaction is developed at t = 1.9 ms with a maximum value of 15.41 MN when the front plate is ruptured. Note that the values shown in Fig. 7 are multiplied by four since a quarter of the problem is modeled due to symmetry. The significant initial resistance is because of the composite reaction of the front plate and concrete core and also the inertia force. Once the front plate is penetrated, the reaction drops to zero and then becomes negative due to rebound and elastic deformation of the panel.

One challenge in using SPH method is the proper application of the boundary conditions. Once the elements are converted to particles, the initial boundary conditions are not applied anymore, and particles fly freely in all directions. The reason is that the symmetry boundary conditions are initially applied to the nodes on the symmetry planes and it is not possible to define the symmetry boundary conditions to the particles in the volume. The simplest solution is to develop a full model instead of a quarter model. However, considering the high computational cost of analysis for the SPH technique, an alternative solution was developed by placing two rigid plates on the symmetry planes (XZ and YZ) as shown in Fig. 8(a). These planes do not have any interaction with the other parts while they act as boundaries, preventing particles moving in those directions. The results from the SPH model are shown in Fig. 8(b). The simulations were performed using a computer with a 3.4 GHz i7 Quad-core processor. The analysis time for the SPH model was 4.7 h, which was more than three times longer than that needed for the conventional FE analysis. In this study, an element deletion limit in tension was not defined in the conventional FE model as it was found to introduce a large error in the predicting the penetration depth. This issue is caused by the impact wave traveling through the concrete block and resulting in compressive and tensile strains. Since the wave travels faster than the projectile, it results in the deletion of elements before the projectile reaches to them. Therefore, no resistance is provided by the elements to the projectile. This issue is one of the disadvantages of the element deletion approach in conventional FE method. On the other hand, the particles still have resistance against penetration after conversion in the SPH model. The rear deformation in the SPH model was 41 mm, which was very close to the experiment. Comparing Figs. 6(b) and 8(b) where the projectile velocity is zero, shows higher penetration in the conventional FE case where penetration could be directly estimated by measuring the location of the projectile. While in the SPH model, the penetration was calculated based on the depth of the last converted element. In other words, the depth of the crushed region is considered as the penetration depth.

As shown in Fig. 8(b), since there is no element removal in the SPH model and due to a lack of compressibility in the particles, the missile is not penetrated through the concrete core. This issue is due to the confinement provided by the steel plates. The simulation results are compared with the experimental data in Table 7. The results showed that the SPH model can predict the structural response more accurately than the conventional FE method. However, the main disadvantage of the SPH method is the high cost of computation.

The iso-surface equivalent plastic strain is shown in Fig. 9 for three different time steps. The figure confirmed that considering a 350 mm€×€350 mm impact region with a fine mesh, see Fig. 2, is a reasonable assumption. No significant plastic deformation is observed outside of this region.

Analytical solution

Multiple empirical models are available in the literature, mostly obtained from curve fitting of projectile shooting experiments. Therefore, the ability of the empirical models to predict the projectile impact response is dependent on the range of parameters considered in the tests. Some of the most popular models are compared here against computer models. A detailed description of these empirical models is provided by Kennedy [⁵²]. The models used in this section estimate the penetration depth with input parameters in US Customary Units. The notation for the input parameters and the conversion of the US Customary units to Standard International (SI) Units is provided in Table 8.

Modified National Defense Research Committee

One of the most well-known empirical models was proposed by Modified National Defense Research Committee (NDRC) [⁵²] for estimating a rigid projectile penetration into a thick concrete target according to

(25)X=\sqrt{{{\displaystyle \mathrm{4}KNWd\left({\displaystyle \frac{v}{\mathrm{1000}d}}\right)}}^{\mathrm{1.8}},}\mathrm{f}\mathrm{o}\mathrm{r}X/d\le \mathrm{2.0},

(26)X={{\displaystyle KNW\left({\displaystyle \frac{v}{\mathrm{1000}d}}\right)}}^{\mathrm{1.8}}+d,\mathrm{f}\mathrm{o}\mathrm{r}X/d>\mathrm{2.0},

where K is the concrete penetrability factor obtained from 180/(f€′_c€)^0.5, and N is the shape factor of the missile, which is 0.72 for a flat-nose projectile. The NDRC equations were obtained based on experimental data with an impact velocity higher than 152 m/s and a projectile diameter ranging from 25 to 406 mm.

Kar formula

Kar [⁵³] proposed a model based on the NDRC model to be applied to nuclear power plants. The modified equation, given as follows, also accounts for the missile material

(27)\begin{array}{c}{\displaystyle \frac{X}{d}}={\mathrm{2}G}^{\mathrm{0.5}},\mathrm{f}\mathrm{o}\mathrm{r}\hspace{0.17em}G\le \mathrm{1},\hfill \\ {\displaystyle \frac{X}{d}}=G+\mathrm{1},\mathrm{f}\mathrm{o}\mathrm{r}\hspace{0.17em}G\le \mathrm{1},\hfill \end{array}

(28)G=\mathrm{3.8}\times {\mathrm{10}}^{-\mathrm{5}}{{\displaystyle \left({\displaystyle \frac{E}{{{\displaystyle E}}_{\mathrm{s}}}}\right)}}^{\mathrm{1.25}}{\displaystyle \frac{N\times M}{d\sqrt{{{\displaystyle f}}_{\mathrm{c}}}}}\times {{\displaystyle \left({\displaystyle \frac{{{\displaystyle V}}_{\mathrm{0}}}{d}}\right)}}^{\mathrm{1.8}},

where E and E_s are Young’s modulus of mild steel and the projectile, respectively.

Army Corps of Engineers

US Army Corps of Engineers (ACE) proposed the following equation based on the NDRC formula

(29)X=({\displaystyle \frac{\mathrm{282.6}}{\sqrt{{{\displaystyle f}}_{\mathrm{c}}^{\text{'}}}}}\hspace{0.17em}{\displaystyle \frac{W}{{{\displaystyle d}}^{\mathrm{3}}}}{d}^{\mathrm{0.215}}{{\displaystyle \left({\displaystyle \frac{V}{\mathrm{1000}}}\right)}}^{\mathrm{1.5}}+\mathrm{0.5}\times d.

This equation applies to a velocity range of 200 to 1000 m/s. The reinforcement percentage was low in their test; however, since the reinforcement did not have a significant impact on the penetration depth, no reinforcement parameter was included in the proposed model [⁵²].

Bechtel Power Corporation

Bechtel Power Corporation proposed the following equation for estimating the penetration of a solid steel projectile in concrete [⁵⁴]

(30)X=\mathrm{15.5}\left({\displaystyle \frac{W^{\mathrm{0.4}}{V}^{\mathrm{0.5}}}{{{\displaystyle f}}_{\mathrm{c}}^{\text{'}\mathrm{0.5}}{d}^{\mathrm{1.2}}}}\right).

Comparison of analytical models with FE simulations

A separate conventional FE model of a plain concrete block was developed with the same dimensions as described in Table 1 to compare with the empirical equations presented above. The results of this comparison are shown in Fig. 10. Note that the empirical equations are a function of limited parameters and valid for the ranges of those parameters they are calibrated for. Therefore, the empirical equations could be used to approximately estimate the penetration depth and a detailed FE analysis is required for more accurate results. From the results shown in Fig. 10, it could be concluded that plain concrete of 672 mm thickness can stop a projectile with an impact velocity of up to 300 m/s.

Parametric analysis

Plate thickness

A parametric analysis was performed to determine the optimal design of SCSS panels in terms of front and rear steel plate thicknesses. The optimal design is assumed to be achieved when the penetration depth and rear deformation are minimized. The same conventional FE model is used throughout this section for faster computation. For this purpose, the front and the rear plate thicknesses are changed, and the analyses are repeated. Three different models were developed based on the plate thicknesses shown in Table 9. Case 1 in Table 9 is the same as in the experiments described in Section 2. In the subsequent cases, the front plate thickness is increased, and the rear plate thickness is reduced to have a constant steel weight in the SCSS panel. The compressive strength of concrete for the analyses was as 23 MPa, identical to that in the experiments.

As shown in Table 9, increasing the front plate thickness significantly reduces the penetration depth of the projectile. However, more deformation is observed in the rear plate. The rear deformation is a function of the composite interaction between the concrete block and the rear plate. Since there is no reinforcement in the SCSS configuration, the tensile force developed by the flexural deformation of the panel is solely carried by the rear plate. Using a thick rear plate is beneficial in the case of the soft projectile such as an aircraft or an aircraft engine when the penetration is not a concern. While in the case of a small diameter hard projectile, specifically with a sharp nose, a thicker front plate is necessary since this is the most effective element in the structure to reduce the velocity of the projectile. The rear plate displacement, projectile velocity, and penetration depth are plotted in Figs. 11, 12, and 13, respectively. A comparison of these results shows that increasing the front plate thickness from 19 to 25 mm reduces the penetration depth by 56%. The rear deformation, however, is increased by 14%. By switching the rear and front plate, cases 1 and 3, the penetration is reduced by 71% and the rear deformation increased by 28%. This indicates that the thickness of the front plate mainly affects the penetration depth.

Impact velocity

In the second part of the parametric analysis, the effect of impact velocity on the penetration depth of the SCSS panel is investigated. The kinetic energy is a function of the mass of the projectile and the impact velocity; therefore, both parameters are indirectly covered in this section. The loading scenarios are summarized in Table 10. Note that the kinetic energy was obtained from, E= 1/2 mv², with the assumed velocities. In this section, the front plate thickness was taken 19 mm, while the same for the rear plate was 32 mm, which corresponds to Case 1 in Fig. 11. The diameter of the projectile was kept constant at 155 mm.

As shown in Table 10, for velocities higher than 400 m/s, the projectile penetrates through the concrete core. For velocities higher than 600 m/s, the rear plate is also perforated. In other words, a kinetic energy of 6250 kJ corresponding to a 600 m/s impact velocity is the failure limit for the studied composite panel. Compared the 300 m/s rupture velocity of the plain concrete core obtained in Section 4, by placing the steel plate, the rupture velocity of the panel is doubled. Based on these results, the composite interaction of concrete core and the rear plate is the main factor in reducing the rear deformation. According to Table 10, in both cases of 3 and 4, the concrete core is perforated. In case 3, the projectile only touches the rear plate, while in case 4, the rear plate is also completely perforated. This indicates that the rear plate can reduce the velocity of the projectile up to 100 m/s with about 50 mm deformation before rupture. The penetration level for the two impact velocities is shown in Fig. 14.

Figure 15 shows the residual velocity of the projectile as a function of the initial velocity. A linear trend is observed for the residual versus initial velocity. This indicates that the protection that the panel provides varies with the initial velocity of the projectile and it is not a constant. In other words, the ability of the panel to reduce the projectile velocity increases as the initial velocity increases.

Impact angle

The effect of the impact angle is studied in this section. It is very unlikely that a projectile will hit a panel in a perpendicular direction; while, most of the empirical equations and experiments are based on perpendicular impact. Therefore, to investigate this issue, additional simulations at different impact angles as shown in Fig. 16 are performed here. These results could also be used for structures with curved surfaces such as containment buildings and dry storage systems. In this section, Case 1 in Table 9 is used for the geometry with an impact velocity of €314 m/s. Different impact angles were considered for the projectile and the penetration depth as well as the rear deformation were measured for all the cases. Considering the nature of the problem, half of the panel was modeled and proper boundary conditions were applied. A summary of the results is presented in Table 11. Note that the penetration depth was measured in the perpendicular direction and the rear deformation was taken as the maximum deformation in the rear plate, which is not necessarily at the center of the panel. Figures 16 and 17 show the condition of the panels for Case 5 and Fig. 18 shows the front plate for Case 6 where the missile only scratches the plate and reflects back.

Result in Table 11 indicate that an impact angle of 75° has a higher penetration depth compared to the perpendicular impact. This is due to the stress concentrations at the corner of the missile, which causes an easier rupture of the front plate. On the other hand, for less than 30 °, the missile tends to reflect from the steel surface. This behavior becomes more dominant at lower angles as shown in Figs. 18 and 19. Only a 74.7 mm perforation is observed for a 30-degree impact angle. The trend of the penetration depth and the maximum rear deformation are presented in Fig. 19. Both the penetration depth and the rear deformation have the same behavior where the maximum value occurs around 75°. NDRC equation considers a 0.72 reduction coefficient on the penetration depth for a flat nose missile. Since the probability of a perpendicular impact in reality is very low, the application of this coefficient for design is an underestimation.

Projectile diameter

In the last part of the parametric study, the diameter of the projectile is considered for damage evaluation of the panel. In this section, the diameter is the only variable of the problem, and the mass and velocity and consequently, kinetic energy remains the same. The purpose of this approach is to focus on the stress concentration in the impact zone, which is the most important parameter for the penetration depth in addition to the kinetic energy. The material stiffness also remains the same so that no major energy dissipation occurs due to the projectile. The diameter of the projectile was increased to 300 from 155 mm in the reference case. As shown in Fig. 20, the missile rebounded from the panel face and despite considerable cracking and crushing in the concrete core, no penetration was observed. The minimum velocity to perforate through the front plate was estimated as 950 m/s for 300 mm diameter. The time history of the rear deformation at the center and a summary of the results are provided in Fig. 21 and Table 12, respectively.

The results in this section show that when the diameter increased by 100%, the nature of the problem completely changes from punching shear to flexure. Although the limited deformation in the projectile is classified as hard impact, the size of the impact area determines the failure mode. The results showed that the penetration depth is highly affected by the impacted area, which is a function of the nose shape and the diameter. The rebound velocity defined as the velocity of the projectile in the reverse direction after the impact increases with increasing impact velocity up to 400 m/s and then remains constant. This is due to the plastic damage in the concrete core at higher impact velocities, which absorbs most of the energy. It has been shown in the previous section that the impact angle is a critical parameter and the reduction coefficient might not yield conservative results for a flat nose projectile. Based on the findings in this and the previous sections, it was found that multiple simulations with different nose shapes and impact angles for safety evaluation are not necessary, only a perpendicular impact with a sharp nose projectile having the minimum expected diameter would provide the worst-case scenario.

Conclusions

In this study, using nonlinear inelastic analysis, the structural response of RC and SCSS panels under high-velocity impact loading of a hard projectile is studied. Two main parameters including penetration depth and global deformation of the panel are investigated. Additionally, two modeling techniques: conventional FE and SPH were used and compared with experimental results to investigate the competency of each technique. Finally, the most effective configuration was identified via a parametric analysis. The following conclusions are made from the findings:

1) Due to a lack of failure criterion, the CDP model is unable to simulate high-speed perforation and it is only suitable for soft projectile impact. To simulate the hard impact where a deep perforation is expected, application of a user subroutine is shown to mitigate this issue.

2) The results of the SPH model were found to be closer to the experimental results compared to the FE method. However, SPH method needs significantly more computational time.

3) The results showed that by changing the steel plate thicknesses in the SCSS configuration, it is possible to reduce the penetration depth up to 71% with the same amount of steel used in the panels, although, this change increases the rear deformation by 28%.

4) The SCSS configuration was found to be more effective in terms of reducing the penetration depth and increasing the cohesiveness of the structure impacted by projectiles. The limit of rupture velocity of the panel increased up by 100% when concrete is sandwiched between steel plates. Conventional RC panel is not suitable to stop a high-velocity hard projectile.

5) To increase safety, it is recommended to consider a thicker front plate compared to the rear plate to minimize the penetration depth. Such an approach shifts the failure of the protective panels from punching shear to flexural yielding, which is a more desirable mode of failure due to higher amount of energy dissipation.

6) The results from the impact angle analyses showed that, for a flat nose missile, with a slight change in the impact angle from the perpendicular direction, the penetration could increase by 34%. This indicates that the reduction coefficient commonly considered for the flat nose missile in the empirical equations could underestimate the penetration depth. The worst-case impact angle for a flat-nose missile is between 75° and 90°.

7) The two most important parameters in the hard projectile impact analysis were found to be the missile diameter and potential stress concentrations due to impact orientation.

References

[1] Nie J G, Hu H S, Fan J S, Tao M X, Li S Y, Liu F J. Experimental study on seismic behavior of high-strength concrete filled double-steel-plate composite walls. Journal of Constructional Steel Research, 2013, 88: 206–219

[2] Takeuchi M, Narikawa M, Matsuo I, Hara K, Usami S. Study on a concrete filled structure for nuclear power plants. Nuclear Engineering and Design, 1998, 179(2): 209–223

[3] Sohel K M A, Richard Liew J Y, Koh C G. Numerical modelling of lightweight steel-concrete-steel sandwich composite beams subjected to impact. Thin-walled Structures, 2015, 94: 135–146

[4] Ebad Sichani M, Padgett J E, Bisadi V. Probabilistic seismic analysis of concrete dry cask structures. Structural Safety, 2018, 73: 87–98

[5] Jiang H, Chorzepa M G. Aircraft impact analysis of nuclear safety-related concrete structures: A review. Engineering Failure Analysis, 2014, 46: 118–133

[6] Huang C C, Wu T Y. A study on dynamic impact of vertical concrete cask tip-over using explicit finite element analysis procedures. Annals of Nuclear Energy, 2009, 36(2): 213–221

[7] Ebad Sichani M, Hanifehzadeh M, Padgett J E, Gencturk B. Probabilistic analysis of vertical concrete dry casks subjected to tip-over and aging effects. Nuclear Engineering and Design, 2019, 343: 232–247

[8] Riera J D. On the stress analysis of structures subjected to aircraft impact forces. Nuclear Engineering and Design, 1968, 8(4): 415–426

[9] Hanifehzadeh M, Gencturk B, Mousavi R. A numerical study of spent nuclear fuel dry storage systems under extreme impact loading. Engineering Structures, 2018, 161(1): 68–81

[10] Sohel K, Liew J R. Behavior of steel-concrete-steel sandwich slabs subject to impact load. Journal of Constructional Steel Research, 2014, 100: 163–175

[11] Abdel-Kader M, Fouda A. Effect of reinforcement on the response of concrete panels to impact of hard projectiles. International Journal of Impact Engineering, 2014, 63: 1–17

[12] Bruhl J C, Varma A H, Kim J M. Static resistance function for steel-plate composite (SC) walls subject to impactive loading. Nuclear Engineering and Design, 2015, 295: 843–859

[13] Heckötter C, Vepsä A. Experimental investigation and numerical analyses of reinforced concrete structures subjected to external missile impact. Progress in Nuclear Energy, 2015, 84: 56–67

[14] Sudhir Sastry Y B, Budarapu P R, Krishna Y, Devaraj S. Studies on ballistic impact of the composite panels. Theoretical and Applied Fracture Mechanics, 2014, 72: 2–12

[15] Budarapu P R, Gracie R, Yang S W, Zhuang X, Rabczuk T. Efficient coarse graining in multiscale modeling of fracture. Theoretical and Applied Fracture Mechanics, 2014, 69: 126–143

[16] Wu Y, Wang D, Wu C T. Three dimensional fragmentation simulation of concrete structures with a nodally regularized meshfree method. Theoretical and Applied Fracture Mechanics, 2014, 72: 89–99

[17] Wu Y, Wang D, Wu C T, Zhang H. A direct displacement smoothing meshfree particle formulation for impact failure modeling. International Journal of Impact Engineering, 2016, 87: 169–185

[18] Lucy L. A numerical approach to the testing of the fission hypothesis. Astronomical Journal, 1977, 82: 1013–1024

[19] Gingold R A, Monaghan J J. Smoothed particle hydrodynamics: Theory and application to non-spherical stars. Monthly Notices of the Royal Astronomical Society, 1977, 181(3): 375–389

[20] Liu G R, Liu M B. Smoothed Particle Hydrodynamics: A Meshfree Particle Method. Singapore: World Scientific, 2003

[21] Rabczuk T, Belytschko T. Cracking particles: A simplified meshfree method for arbitrary evolving cracks. International Journal for Numerical Methods in Engineering, 2004, 61(13): 2316–2343

[22] Rabczuk T, Zi G. A meshfree method based on the local partition of unity for cohesive cracks. Computational Mechanics, 2007, 39(6): 743–760

[23] Rabczuk T, Eibl J. Modelling dynamic failure of concrete with meshfree methods. International Journal of Impact Engineering, 2006, 32(11): 1878–1897

[24] Randles P, Libersky L. Smoothed particle hydrodynamics: Some recent improvements and applications. Computer Methods in Applied Mechanics and Engineering, 1996, 139(1–4): 375–408

[25] Ren H, Zhuang X, Cai Y, Rabczuk T. Dual-horizon peridynamics. International Journal for Numerical Methods in Engineering, 2016, 108(12): 1451–1476

[26] Hughes G. Hard missile impact on reinforced concrete. Nuclear Engineering and Design, 1984, 77(1): 23–35

[27] Dancygier A, Yankelevsky D. High strength concrete response to hard projectile impact. International Journal of Impact Engineering, 1996, 18(6): 583–599

[28] Dassault Systemes. ABAQUS, 6.14, Simulia Corp., Providence, RI, USA, 2016

[29] Hibbitt H D, Karlsson B I, Sorensen E P. ABAQUS User’s & Theory Manuals, 6.14, Dassault Systèmes Simulia Corp., Providence, RI, USA, 2013

[30] Lee S, Cho S S, Jeon J E, Kim K Y, Seo K S. Impact analyses and tests of concrete overpacks of spent nuclear fuel storage casks. Nuclear Engineering and Technology, 2014, 46(1): 73–80

[31] Malvar L J. Review of static and dynamic properties of steel reinforcing bars. Materials Journal, 1998, 95(5): 609–616

[32] Fu H, Erki M, Seckin M. Review of effects of loading rate on concrete in compression. Journal of Structural Engineering, 1991, 117(12): 3645–3659

[33] Ross C A, Tedesco J W, Kuennen S T. Effects of strain rate on concrete strength. Materials Journal, 1995, 92(1): 37–47

[34] Leppänen J. Concrete Structures Subjected to Fragment Impacts. Goteborg: Chalmers University of Technology, 2004

[35] Li Q, Meng H. About the dynamic strength enhancement of concrete-like materials in a split Hopkinson pressure bar test. International Journal of Solids and Structures, 2003, 40(2): 343–360

[36] CEB-FIP. Design of Concrete Structures. Euro-International Committee for Concrete (CEB). Lausanne, Switzerland, 1993

[37] UFC. Unified facilities criteria (UFC), structures to resist the effects of accidental explosions, UFC 3-340-02. U.S. Department of Defence, 2008

[38] Kupfer H, Hilsdorf H K, Rusch H. Behavior of concrete under biaxial stresses. ACI Journal Proceedings, 1969, 66(8): 656–666

[39] Kupfer H B, Gerstle K H. Behavior of concrete under biaxial stresses. Journal of the Engineering Mechanics Division, 1973, 99(4): 853–866

[40] Lee J, Fenves G L. Plastic-damage model for cyclic loading of concrete structures. Journal of Engineering Mechanics, 1998, 124(8): 892–900

[41] Hanifehzadeh M, Gencturk B, Willam K. Dynamic structural response of reinforced concrete dry storage casks subjected to impact considering material degradation. Nuclear Engineering and Design, 2017, 325: 192–204

[42] Jankowiak T, Lodygowski T. Identification of parameters of concrete damage plasticity constitutive model. Foundations of civil and environmental engineering, 2005, 6(1):53–69

[43] Hognestad E. Study of Combined Bending and Axial Load in Reinforced Concrete Members. Urbana: University of Illinois at Urbana Champaign, 1951

[44] Wang T, Hsu T T. Nonlinear finite element analysis of concrete structures using new constitutive models. Computers & Structures, 2001, 79(32): 2781–2791

[45] Harris H G, Sabnis G. Structural modeling and experimental techniques. CRC Press, Taylor & Francis Group, 1999

[46] Cowper G R, Symonds P S. Strain-hardening and Strain-rate Effects in the Impact Loading of Cantilever Beams. DTIC Document, 1957

[47] Malvar L J, Crawford J E. Dynamic increase factors for steel reinforcing bars. In: 28th DDESB Seminar, Orlando, USA, 1998

[48] ASTM. Standard Specification for Carbon Structural Steel, A36/A36M. West Conshohocken: ASTM International, 2014

[49] ASTM. Standard Specification for Deformed and Plain Carbon-Steel Bars for Concrete Reinforcement, ASTM A615/A615M-16. West Conshohocken: ASTM International, 2016

[50] Kawamoto Y, Stepan J. Analytical Study of Reinforced Concrete Slab Subjected to Soft Missile Impact. SMiRT-23, Manchester, UK, 2015

[51] Rodríguez Soler J, Martinez Cutillas F J, Marti Rodriguez J. Concrete constitutive model, calibration and applications. In: SIMULIA Community Conference, Vienna, Austria, 2013

[52] Kennedy R. A review of procedures for the analysis and design of concrete structures to resist missile impact effects. Nuclear Engineering and Design, 1976, 37(2): 183–203

[53] Kar A K. Local effects of tornado-generated missiles. Journal of the Structural Division, 1978, 104(5): 809–816

[54] Linderman R, Rotz J, Yeh G. Design of structures for missile impact. San Francisco: Bechtel Power Corp., 1974

RIGHTS & PERMISSIONS

Higher Education Press and Springer-Verlag GmbH Germany, part of Springer Nature