Numerical Material Model for Composite Laminates in High-Velocity Impact Simulation

Abstract

A numerical material model for composite laminate, was developed and integrated into the nonlinear dynamic explicit finite element programs as a material user subroutine. This model coupling nonlinear state of equation (EOS), was a macro-mechanics model, which was used to simulate the major mechanical behaviors of composite laminate under high-velocity impact conditions. The basic theoretical framework of the developed material model was introduced. An inverse flyer plate simulation was conducted, which demonstrated the advantage of the developed model in characterizing the nonlinear shock response. The developed model and its implementation were validated through a classic ballistic impact issue, i.e. projectile impacting on Kevlar29/Phenolic laminate. The failure modes and ballistic limit velocity were analyzed, and a good agreement was achieved when comparing with the analytical and experimental results. The computational capacity of this model, for Kevlar/Epoxy laminates with different architectures, i.e. plain-woven and cross-plied laminates, was further evaluated and the residual velocity curves and damage cone were accurately predicted.

Keywords:
Material model; composite laminate; high-velocity impact; numerical simulation

1 INTRODUCTION

Composite materials, due to their inherently superior mechanical properties such as high-specific strength and stiffness, have been widely used in aerospace, civil, and armor protection applications (Zhang et al., 2008; Tabiei and Ivanov, 2002). To improve the ballistic performance of composite laminate, it is essential to further study its ballistic impact response and better understand its impact/penetration behaviors. Over the past decade, considerable research work has been done by many scholars from different angles or aspects. The efforts to study the ballistic behaviors of composites mainly focus on the ballistic limit, residual velocity, energy absorption and failure modes. Three methods, i.e. experiment, theoretical analysis and numerical simulation are selected to conduct the studies.

Theoretical model can explicitly express the relationship between the ballistic performance of composites and influential parameters such as impact velocity, nose shape of projectiles, material characteristics and boundary conditions. Analytical study is a powerful tool to completely understand the damage and energy dissipation mechanisms of composites (Naik and Doshi, 2005). Analytical investigations for composite materials under low to high impact velocities have been carried out and several theoretical models have been also proposed by many researchers (Wen, 2000; Billon and Robinson, 2001; Naik and Doshi, 2005). Experimental study is another significant method because it is not only closer to the real physical process but also is the base of the other methods. Many experimental studies are reported in the literatures (Zhu et al., 1992; Wang et al., 2005; Babu et al., 2007; Zu et al. 2015), which give insight of the anti-ballistic mechanism of composites, and provide a large amount of data for the future studies. Numerous intervening parameters as mentioned above, make the impact issue of composite materials are rather complicated. Theoretical models usually depend on many assumptions and simplified conditions, which seriously hinder their utility for general engineering problems. Moreover, experiment methods are usually subjected to the technical conditions, and the internal impact process is difficult to observe. Therefore, designing composite armors or shields solely based on the theoretical models and experimental data is uneconomical and impracticable (Silva et al., 2005; Xin and Wen, 2012).

Fortunately, recent advances of computer science and numerical algorithms, offer a possibility of using numerical simulation combined with slight experiment tests for evaluating composite materials. The material model and modeling method for composite play a key role in numerical simulation. Hoof (1999) developed a numerical model for composite laminate. In this model, composite laminates were constructed ply by ply and a tie-break contact algorithm was introduced to simulate the interlaminar failure. Tabiei and Ivanov (2002) presented a micromechanical model for flexible woven fabric. The yarns’ reorientation and fabric architecture were considered in this model. It could simulate the dual behaviors of the flexible fabric, i.e. the trellis mechanism behavior (before yarns locking) and generally anisotropic elastic properties (after yarns packing). Grujicic (2012) developed a fiber-level model for Kevlar^® KM2 ballistic fabric, in which yarn was represented as a bundle of collinear discrete fibers, and each fiber was discretized into a number of three-dimensional beam elements. Using this fiber-level model, the influence, of fiber transverse properties and the friction among fibers, on the ballistic behavior of fabric have been investigated. Anderson (1994) presented a constitutive relationship for anisotropic materials, in which the coupling problem of the volumetric and deviatoric had been overcome. Hayhurst (1999), Clegg (2006) and Riedel (2006) et al. implemented Anderson’ approach and developed an advanced material model for Nextel and Kevlar/epoxy cloth. This model was based on smoothed-particle hydrodynamics (SPH) algorithm and used to predict the hypervelocity impact (HVI) behaviors of Nextel and Kevlar/epoxy cloth.

Composite laminate is fabricated by several single fabrics, and fabric is assembled by a number of yarns in specific patterns such as weaving or stacking. While the yarns are usually made by hundreds of the high-strength fibers (Grujicic et al., 2012). The more meticulous of the numerical model will the more reflect the microstructural of composites, and the more accuracy/fidelity of the numerical results will be obtained. However, it will make the computational scale increase rapidly, and thus result in a dramatic decline in computational efficiency. On the other hand, for the engineering issues, e.g. designing of composite armors or shields, the macroscopic parameters, such as ballistic limit, residual velocity, energy absorption and penetration/perforation behavior, are the most important design references.

Aiming at such problem, a computational macro-mechanical model was developed in this paper, which was based on Anderson’ methodology and had been implemented in the nonlinear dynamic explicit finite element code LS-DYNA as a user-defined material model. The organization of the paper is as follows. The theoretical framework of this model is discussed and an inverse flyer plate simulation is conducted to demonstrate the nonlinear shock response of this developed model in Section 2. To validate this model, a high-velocity penetration problem for Kevlar29/Phenolic laminate has been investigated by numerical simulation in Section 3. The computational capacity of this model for Kevlar/Epoxy laminates with different architectures is further evaluated in Section 4. The main conclusions are summarized in Section 5.

2 MATERIAL MODEL

As mentioned above, composite materials are generally synthesized by two or more materials through chemical method and possess themselves hierarchical structure. To enhance the engineering practicability and improve the computational efficiency, the composite material model developed in this paper is a computational macro-mechanical model, which equates the whole composite laminate into an orthotropic homogeneous material as shown in Figure 1. This model accounting for the nonlinear shock effects, aims to capture the main mechanical behaviors and evaluate the macroscopic parameters of composite laminates under high-velocity impact loading.

Figure 1
Schematic representation of the equivalent conversion for composite laminate.

Comparing with isotropic materials (e.g. metallic materials), the behaviors of anisotropic materials subject to impact seem to be more complex. The phenomena observed in the impact tests for composite laminates are as follows (Hayhurst et al., 1999; Silva et al., 2005):

1. Material anisotropy

2. Shock response

3. Anisotropic strength degradation

4. Coupling of volumetric and deviatoric response

5. Delamination

2.1 Constitutive Model and EOS

Different to isotropic materials, composite laminate is usually orthotropic and the strength effects (constitutive model) and pressure effects (EOS) are strongly coupled together. To deal with this problem, Anderson (1994) have proposed a method, in which the total strain, ( _ij has been split into deviatoric, ${\epsilon }_{ij}^d$ and volumetric, ( _vol components. Therefore, the stress-strain relationship of composite laminate can be expressed as:

$\left[\begin{array}{l}{\sigma }_{11}\\ {\sigma }_{22}\\ {\sigma }_{33}\\ {\sigma }_{23}\\ {\sigma }_{31}\\ {\sigma }_{12}\end{array}\right]=\left[\begin{array}{cccccc}{C}_{11}& C_{12}& C_{13}& 0& 0& 0\\ C_{21}& C_{22}& C_{23}& 0& 0& 0\\ C_{31}& C_{32}& C_{33}& 0& 0& 0\\ 0& 0& 0& C_{44}& 0& 0\\ 0& 0& 0& 0& C_{55}& 0\\ 0& 0& 0& 0& 0& C_{66}\end{array}\right]\left[\begin{array}{c}{\epsilon }_{11}^d+\frac{1}{3}{\epsilon }_{vol}\\ {\epsilon }_{22}^d+\frac{1}{3}{\epsilon }_{vol}\\ {\epsilon }_{33}^d+\frac{1}{3}{\epsilon }_{vol}\\ {\epsilon }_{{}_{23}}^d\\ {\epsilon }_{{}_{31}}^d\\ {\epsilon }_{{}_{12}}^d\end{array}\right]$ (1)

where σ_ij is the stress tensor; C _ij are 9 independent elastic constants.

Pressure is defined as a third of the trace of the stresses, i.e. P = -( σ₁₁ = σ₂₂ + σ₃₃)/3. Therefore, expanding Eq. (1), the expression of pressure can be obtained:

$$\begin{array}{l}P=-\frac{1}{9}\left(C_{11}+C_{22}+C_{33}+2C_{12}+2C_{23}+2C_{31}\right){\epsilon }_{\nu ol}\\ -\frac{1}{3}\left(C_{11}+C_{12}+C_{13}\right){\epsilon }_{11}^d-\frac{1}{3}\left(C_{21}+C_{22}+C_{23}\right){\epsilon }_{22}^d-\frac{1}{3}\left(C_{31}+C_{32}+C_{33}\right){\epsilon }_{33}^d\end{array}$$ (2)

Moreover, the deviatoric stress, $${\sigma }_{ij}^d$$ can be expressed as:

$$\{\begin{array}{c}\begin{array}{l}{\sigma }_{11}^d=\frac{1}{9}\left(2C_{11}-C_{22}-C_{33}+C_{12}-2C_{23}+C_{31}\right){\epsilon }_{vol}\\ \begin{array}{cc}& \end{array}+\frac{1}{3}\left(2C_{11}-C_{12}-C_{13}\right){\epsilon }_{11}^d+\frac{1}{3}\left(2C_{12}-C_{22}-C_{23}\right){\epsilon }_{22}^d+\frac{1}{3}\left(2C_{31}-C_{32}-C_{33}\right){\epsilon }_{33}^d\end{array}\\ \\ \begin{array}{l}{\sigma }_{22}^d=\frac{1}{9}\left(-C_{11}+2C_{22}-C_{33}+C_{12}+C_{23}-2C_{31}\right){\epsilon }_{vol}\\ \begin{array}{cc}& +\frac{1}{3}\left(-C_{11}+2C_{12}-C_{13}\right){\epsilon }_{11}^d+\frac{1}{3}\left(-C_{12}+2C_{22}-C_{23}\right){\epsilon }_{22}^d+\frac{1}{3}\left(-C_{31}+2C_{32}-C_{33}\right){\epsilon }_{33}^d\end{array}\end{array}\\ \\ \begin{array}{l}{\sigma }_{33}^d=\frac{1}{9}\left(-C_{11}-C_{22}+2C_{33}-2C_{12}+C_{23}+C_{31}\right){\epsilon }_{vol}\\ \begin{array}{cc}& +\frac{1}{3}\left(-C_{11}-C_{12}+2C_{13}\right){\epsilon }_{11}^d+\frac{1}{3}\left(-C_{12}-C_{22}+2C_{23}\right){\epsilon }_{22}^d+\frac{1}{3}\left(-C_{31}-C_{32}+2C_{33}\right){\epsilon }_{33}^d\end{array}\end{array}\end{array}$$ .(3)

From the Eq. (2) and Eq. (3), it is can be seen clearly that both the volumetric and deviatoric components of strain have the contributions to the pressure and deviatoric stress. For the composite materials, a significant conclusion can be drawn that the volumetric and deviatoric response are strongly coupled, as that deviatoric strain can result in spherical stress, while volumetric strain can lead to deviatoric stress.

The first term, on the right side of Eq. (2), is a liner EOS for isotropic Hookean materials. However, for orthotropic materials, it represents the volumetric (thermodynamic) response. To involve the nonlinear shock effects, the first term of Eq. (2), i.e. the contribution to pressure from the volumetric strain, is modified by Mie-Grüneisen EOS. Meanwhile, the contribution to pressure from the deviatoric strain is remained as a correction. Thus the expression of pressure can be rewritten as:

$P=P_{EOS}\left(\mu ,e\right)-\frac{1}{3}\left(C_{11}+C_{21}+C_{31}\right){\epsilon }_{11}^d-\frac{1}{3}\left(C_{12}+C_{22}+C_{32}\right){\epsilon }_{22}^d-\frac{1}{3}\left(C_{13}+C_{23}+C_{33}\right){\epsilon }_{33}^d$ (4)

In this study, the polynomial EOS was selected and developed using Fortran language.

$P_{EOS}\left(\mu ,e\right)=\{\begin{array}{c}{A}_1\mu +A_2{\mu }^2+A_3{\mu }^3\left(B_0+B_1\mu \right){\rho }_{0}e\begin{array}{cc}& \mu >0\end{array}\\ T_1\mu +T_2{\mu }^2+B_0{\rho }_{0}e\begin{array}{ccc}\begin{array}{cc}\begin{array}{cc}& \end{array}& \end{array}& & \begin{array}{cc}& \mu <0\end{array}\end{array}\end{array}$ (5)

where, $\mu =\frac{\rho }{{\rho }_0}-1$, A ₁ = T ₁, and the other parameters will be determined by experiments or simulations (Anderson et al., 1994).

2.2 Failure Model

The failure model developed in this paper is an orthotropic failure model, which includes two stages, i.e. initial failure and post-failure response. The failure initiation can be based on any combination of material stress and/or strain, and the material strain criterion is adapted in this paper. Subsequent to failure initiation, stiffness and strength properties of the failed material will be updated based on the direction or modes of the failure (Silva et al., 2005).

$\left[\begin{array}{l}\Delta {\sigma }_{11}\\ \Delta {\sigma }_{22}\\ \Delta {\sigma }_{33}\\ \Delta {\sigma }_{23}\\ \Delta {\sigma }_{31}\\ \Delta {\sigma }_{12}\end{array}\right]=\left[\begin{array}{cccccc}{C}_{11}& C_{12}& C_{13}& 0& 0& 0\\ C_{21}& C_{22}& C_{23}& 0& 0& 0\\ C_{31}& C_{32}& C_{33}& 0& 0& 0\\ 0& 0& 0& \alpha C_{44}& 0& 0\\ 0& 0& 0& 0& \alpha C_{55}& 0\\ 0& 0& 0& 0& 0& \alpha C_{66}\end{array}\right]\left[\begin{array}{c}\Delta {\epsilon }_{11}\\ \Delta {\epsilon }_{22}\\ \Delta {\epsilon }_{33}\\ \Delta {\epsilon }_{23}\\ \Delta {\epsilon }_{31}\\ \Delta {\epsilon }_{12}\end{array}\right]$ (6)

(1) Delamination failure (11-direction)

In this study, the through-thickness direction of composite laminate is defined as the 11-direction, while the in-plane principal directions are 22- and 33-directions. Delamination will be caused by the excessive stresses (or strains) in the 11- or 12-plane direction. If failure is initiated in either of these two modes, the stress in the 11-direction and the corresponding orthotropic stiffness coefficients C _ij is instantaneously set to zero (Tham et al., 2008; Bandaru et al., 2015, 2016):

$$\Delta {\sigma }_{11}=\mathrm{0,}\begin{array}{cc}& C_{ij}\end{array}=C_{ji}=0\begin{array}{cc}& \begin{array}{cc}i=1& 11-direction\end{array}\end{array}$$ (7)

Eq. (7) shows that this failure model is independent of the directions, i.e. failure in 11-direction does not affect 22- and 33-directions. However, the delamination will in practice be associated with a reduction in shear stiffness, and the fractional residual shear stiffness is maintained through the parameter, α (0.0 to 1.0). Usually, the value of the parameter, α is obtained by experimental tests but in the absence of appropriate material data, a nominal value of 20% is typically used (Silva et al., 2005; Tham et al., 2008).

(2) In-plane failure

The 22- and 33-directions are assumed to be in the plane of the composite. If failure is initiated in 22- and 33-directions, the post-failure response is similar to 11-direction,

$$\{\begin{array}{c}\Delta {\sigma }_{22}=\mathrm{0,}\begin{array}{cc}& C_{ij}\end{array}=C_{ji}=0\begin{array}{cc}& \begin{array}{cc}i=2& 22-direction\end{array}\end{array}\\ \Delta {\sigma }_{33}=\mathrm{0,}\begin{array}{cc}& C_{ij}\end{array}=C_{ji}=0\begin{array}{cc}& \begin{array}{cc}i=3& 33-direction\end{array}\end{array}\end{array}$$ (8)

(3) Combined failure

If all the three material directions fail simultaneously, the material stiffness and strength will become isotropic with no stress deviators and tensile material stresses, indicating that the material can withstand only hydrostatic pressure.

2.3 Erosion Model

Due to the orthogonality of composite laminates, some elements in the impact region will be more easily distorted during the simulations, which may lead to numerical instabilities or even calculation failure (Deka et al., 2008). In this study, element effective strain is used as the erosion criterion, similar approach has been adopted in the references (Grujicic et al., 2009; Bandaru et al., 2015, 2016).

$${\epsilon }_{eff}=\frac{\sqrt{2}}{3}\sqrt{{\left({\epsilon }_1-{\epsilon }_2\right)}^2+{\left({\epsilon }_2-{\epsilon }_3\right)}^2+{\left({\epsilon }_3-{\epsilon }_1\right)}^2+6\left({\epsilon }_{12}^2+{\epsilon }_{23}^2+{\epsilon }_{31}^2\right)}$$ (9)

where, ( _i , ( _ij (i=1 to 3, j=1 to 3) is the strain components.

During the impact process, the composite material will undergo non-failure and failure states and the corresponding cell strain will also undergo this two stages. For instance, ( ₁ is the strain in the thickness direction, while ${\epsilon }_1^{failure}$ and ${\epsilon }_1^{erosion}$ are the corresponding failure and erosion strain components, respectively. In the non-failure stage, the ( ₁ will slowly grow with other strain components before it arrives to ${\epsilon }_1^{failure}$. However, after arriving to ${\epsilon }_1^{failure}$, failure occurs, which causes that composite laminates cannot withstand any loading in the thickness direction. Therefore, the ( ₁ will increase dramatically until it arrives to ${\epsilon }_1^{erosion}$ (Figure 2). At this point, the meshes will be extremely distorted and the deletion mechanism will be triggered.

Figure 2
Growth of the strain in the 11-direction of composite laminate.

2.4 Inverse Flyer Plate Simulations

An inverse flyer plate test (IFPT, Figure 3), creating uniaxial compression state under high-velocity impact condition, was conducted by EMI (Ernst Mach Institute) (Hayhurst et al., 1999). The projectile plate consists of a 6 mm thick Kevlar/Epoxy plate and a 4 mm thick C45-steel backing plate whilst a 2 mm thick C45-steel plate is stationary as the witness plate (target plate). The layered projectile, carried by a polymer sabot, is accelerated by means of a single stage gas gun. The projectile’s velocity is measured via trigger pins mounted at the muzzle of the gun barrel. Once the impact event of projectile/target occurs, shock waves will be generated and propagate in the witness plate, and this signal, i.e. free surface velocity history, will be recorded by a laser velocity interferometer VISAR (velocity interferometer system for any reflector) (Lässig T et al., 2015). The corresponding simulation has been implemented in this study, which can be used to verify the shock response of the developed material model (Figure 4).

Figure 3
Experimental setup of IFPT.

Figure 4
Simulation model of IFPT.

For comparison, the standard material model, MAT_2 (MAT_ORTHOTROPIC_ELASTIC) in LS-DYAN has been introduce. The numerical and experimental results are compared in Figure 5.

Figure 5
Comparison of numerical and experimental flyer velocity histories for two material models.

Figure 5 shows that the two simulation results present a dramatic difference. Good match between the developed model case and the experiment is achieved, while the curve obtained from MAT_2 case shows a great deviation. It is encouraging and means that the developed material model is able to simulate the nonlinear shock response, i.e. phenomena (2).

3 VALIDATION OF THE DEVELOPED MATERIAL MODEL

3.1 Numerical Simulation

The theoretical framework of the developed model was discussed and its nonlinear shock response has been also demonstrated in Section 2. In this section, the performance of this model will be further verified in high-velocity impact condition. The high-velocity impact/penetration experiment is usually conducted by one-stage powder/gas gun system, the main components of which are breech, barrel, safety vessel and target chamber (Figure 6). The projectile and its sabot are launched together by the gun and the impact velocity is controlled by the amount of gun powder. In the safety vessel, the projectile and sabot will be separated before arriving the target; the projectile’s velocity is measured by photoelectric diodes and the impact process is recorded by high speed camera.

Figure 6
Schematic illustration of one stage powder gun.

The classic ballistic impact case is selected form the literature (Hoof et al., 1999), in which the composite target consists of 19-plies Kevlar29/Phenolic fabric and its total thickness is 9.5 mm. The projectile is fragment simulating projectile (FSP) based on STANAG 2920 and US MIL-P-46593 (Figure 7). This is a 1.1g FSP that is fabricated using 4340 steel and heated treated to 29 HRC (Hoof et al., 1999; Silva et al., 2005).

Figure 7
Geometrical properties and mesh model of FSP.

Johnson-Cook material model, including strain/strain rate hardening and thermal softening effects, is selected to model the FSP; while the target plate is modeled by the developed material model. The parameters of this two material models are summarized in Table 1 (Ramadhan et al., 2013; Ansari and Chakrabarti, 2016) and Table 2 (Tham et al., 2008; Bandaru et al., 2015), respectively.

Table 1
Johnson-Cook material model parameters for the FSP.

Parameters	4340 steel
Density (g/cm3)	7.83
Bulk modulus (kPa)	1.59E+8
Yield stress (kPa)	7.92E+5
Harding constant (kPa)	5.10E+5
Harding exponent	0.26
Strain rate constant	0.014
Thermal softening exponent	1.03
Melting temperature (K)	1793

Table 2
Composite material model for Kevlar29/Phenolic panels.

Parameters	Kevlar29/Phenolic
Young modulus 11 (kPa)	1.9480E+06
Young modulus 22 (kPa)	1.7989E+07
Young modulus 33 (kPa)	1.7989E+07
Poisson ratio 12	0.0756
Poisson ratio 23	0.08
Poisson ratio 31	0.6981
Shear modulus 12 (kPa)	2.23500E+05
Shear modulus 23 (kPa)	1.85701E+06
Shear modulus 31 (kPa)	2.23500E+05
Tensile failure strain 11	0.02
Tensile failure strain 22	0.06
Tensile failure strain 33	0.06

Brick elements with single integration point, which are more robust than the fully integrated elements, are used to discretize the target and projectile. Each ply of composite laminate is modeled using a single solid element, and thus 19 elements are assigned along the thickness direction. Considering computational accuracy or robust, the cell size (about 0.5 mm) of the projectile remains consistent with the target plate. An eroding surface-to-surface contact algorithm (CONTACT_ERODING_SURFACE_TO_SURFACE) is used to simulate the contact between the projectile and composite laminate. Moreover, stiffness-based hourglass control algorithm has been also introduced to hinder the hourglass mode. The finite element impact model for projectile/target is shown in Figure 8.

Figure 8
Ballistic impact model for projectile/composite laminate.

3.2 Discussion of the Simulation Results

(1) Failure modes

A series of simulations with different impact velocities were conducted. The case with the initial impact velocity of 590 m/s is selected as an example to discuss the simulation results. Figure 9 shows the plot of strain-11 and the red color denotes the failure strain-11 (seen Table 1). The failure mode, in the form of a crucifix, means that cracks occur and evolve along the warp and weft direction on the top surface of the target plate. This cross-shaped damage pattern was also reported in many references (Silva et al., 2005; Tham et al., 2008; Bandaru et al. 2016), which presents the characteristic of anisotropic strength degradation, i.e. phenomena (3) mentioned in Section 2. The damage pattern on the back surface is shown in Figure 10. Due to the effects of shock-wave reflection/composite-laminate bending or both, the stress state on the back surface is tensile. The phenomena of back-face bulging and fiber fracture were observed in the simulation.

Figure 9
Cross-shaped damage on the top surface of target plate.

Figure 10
Back-face bulging and fiber fracture.

During the penetration process, three distinct armor-penetration stages are simulated which correspond to three major failure modes (Cheng, 2003). As shown in Figure 11, in the initial penetration stage, the projectile highly compress the composite plate which lead to shear/cut loading conditions formed at the edge of FSP. Punching shear failure is the dominate failure mode, and the damage of the fibers and matrix is mostly around the rim and underneath the projectile in this stage. With the decrease of the projectile velocity, fibers absorb the projectile’s kinetic energy by stretching, and meanwhile the damage is fully extended to the outer zones. Fibers breakage is the major failure mode in the second stage. Due to the reflection of the compression wave and extensive bulging of the back-face, the final stage is dominated by delamination and extensive stretching of the filaments.

Figure 11
Failure modes in three distinct armor-penetration stages.

These three penetration stages are simulated and shown in Figure 11. The tensile and crimped meshes in the impact zone present the fiber stretching, shearing and crushing. The red color demonstrates the distribution of the failure strains (${\epsilon }_1^{failure}$) indicating the delamination damage of the target plate. In the damage region near the back surface (about 1/3), meshes have been extremely stretched and some even have reached the erosion criteria (${\epsilon }_1^{erosion}$). It just goes to show that more serious delamination failure (or a completely damaged state) has occurred in the final stage. However, the developed model is a macro-mechanical model, and thus the plot of failure strain-11, representing the distribution of the damage, only reflects the macro-mechanical properties or has the average mechanical meaning.

(2) Ballistic limit velocity

Above-mentioned, the failure modes (qualitative analysis) of composite laminates were discussed. To further verify the developed model, the ballistic limit velocity (BLV, quantitative analysis) is also calculated. There are two standard methods to measure the BLV. The first one is the velocity history method and the second one is based on US MIL-STD-662E standard (Bandaru et al., 2016). In the first method, the maximum impact velocity where the projectile can be stopped, is utilized to be defined as the BLV (V _bl ). In the second method, BLV (V ₅₀) is the average of impact velocities with equal number of partial and full penetrations within a small velocity range of 38 m/s. To ensure the numerical reliability and accuracy, the two methods are available and compared in the present study.

The velocity time histories of projectiles with different initial velocities are shown in Figure 12a. The positive residual velocity means that the projectile has fully penetrated the target plate, while the negative represents partial penetration and rebounding. From Figure 12a, it can be observed that for the case of 590 m/s, its residual velocity is -8.03 m/s. It means that the most kinetic energy of the projectile has been absorbed by the target; negative residual velocity means that rebounding is emerged but not very large. Therefore, the impact velocity of 590 m/s can be considered as the BLV (V _bl ). For the case of 610 m/s, its residual velocity is +50.5 m/s which means the target plate has been fully penetrated. Comparing the above two cases, a phenomenon can be observed that slight increase of the initial velocity (that beyond the BLV) leads to a rapid increase of the residual velocity, similar to the report in the reference (Zhu et al., 1992). Figure 12b shows the computational result base on US standard. Six impact simulations have been conducted and among which three cases are full perforation and the other three are partial perforation. The average impact velocity, i.e. V ₅₀ is 586.5 m/s. The BVLs measured by two methods are very close indicating that the numerical result is valid or reliable in turn.

Figure 12
BLVs measured by two methods.

Considerable efforts of building analytic model to predict the BLV for composite laminates have been performed. Based on the assumption that the mean pressure provided by the target consists of two parts, namely quasi-static stress and dynamic stress; Wen (2000) proposed a theoretical model to calculate the BLV for composite laminates:

$$V_{bl}=\frac{\pi \beta \sqrt{{\rho }_t{\sigma }_Y}{D}^{2}T}{4M}\left[1+\sqrt{1+\frac{8M}{\pi {\beta }^2{\rho }_t{D}^{2}T}}\right]$$ (10)

where, M and D are the projectile’s mass and diameter; T, ρ _t and σ_Y are the thickness, density and elastic limit of the composite target, respectively, and β is the geometry parameter of the projectile. In this case, the parameters used to calculate the BLV are listed in Table 3 (Hoof et al, 1999; Wen, 2000).

Table 3
Material and physical properties of the FSP and composite laminate.

σt (g/cm3)	D (mm)	T (mm)	M (g)	σY (Mpa)	β
1.65	5.461	9.5	1.1	194	1.9

The BLVs obtained from simulation (using two methods), analytical equation (Wen, 2000) and experiment (Hoof et al., 1999) are summarized and compared in Table 4. It can be observed that the numerical results either method 1 or method 2 are greatly close to the experiment, but the analytical result is a little conservative. In generally, the errors of the BLVs gotten from the three ways are acceptable, which shows that the developed model is correct or reasonable in terms of predicting the BLV.

Table 4
Comparison of ballistic limit velocity for different methods.

Results	Ballistic limit (m/s)
Experiment	600
Theoretical analysis	566
Simulation method 1	590
Simulation method 2	586.5

4 PREDICTING THE RESIDUAL VELOCITY CURVES OF TWO DIFFERENT KEVLAR/EPOXY LAMINATES

In this section, the computational capability of the developed model, for composite laminates with different architectures, has been evaluated. The residual velocity curves of two different Kevlar/Epoxy laminates, i.e. plain-woven and cross-plied laminate, have been predicted and compared with the experimental results (Wang et al., 2005).

4.1 Numerical Simulation

In the impact tests, composite laminates consist of 75% Kevlar and 25% epoxy (vol. %), which are bonded by chemical method. One kind of laminates is woven from yarns while the other one is stacked in a sequence of o° / 90° (seen Figure 13a and Figure 13b). The thickness of those two kinds of composite laminates is 5 mm with 10 plies, each 0.5-mm thick. The diameter and length of the cylindrical shape projectile are 7.62 mm and 15.24 mm, respectively; the projectile is prepared from 4340 steel and heat-treated to a hardness of HRC 32 ± 1 (Wang et al., 2005).

Figure 13
Composite laminates and finite element impact model.

Not considering the specific structures, both the plain-woven and cross-plied laminates are equated into an orthotropic homogeneous plate, which is modeled by the developed model (seen Figure 13c). Same to the last section, uniform 3D solid elements have been adapted and 20 elements are assigned along the thickness direction. Because a negligible deformation of the projectiles is observed in the test (Wang et al., 2005); the projectile in this simulation is modeled using rigid material model (MAT_20, ρ = 7830kg/m ³, E = 200Gpa, v = 0.3). This is a cost efficient model because that the rigid elements avoid allocating storage to record the history variables during the simulation (Zeng et al., 2005). Table 5 (Xin and Wen, 2012) gives the parameters of composite laminate, and the finite element model is shown in Figure 13.

Table 5
Composite material model for Kevlar/Epoxy panels.

Parameters	Kevlar/Epoxy
Young modulus 11 (kPa)	4.60E+06
Young modulus 22 (kPa)	2.10E+07
Young modulus 33 (kPa)	2.10E+07
Poisson ratio 12	0.14
Poisson ratio 23	0.34
Poisson ratio 31	0.64
Shear modulus 12 (kPa)	1.30E+06
Shear modulus 23 (kPa)	1.30E+06
Shear modulus 31 (kPa)	1.30E+06
Tensile failure strain 11	0.02
Tensile failure strain 22	0.06
Tensile failure strain 33	0.06

4.2 Comparison of the Numerical and Experimental Results

To predict the residual velocities, 12 simulation cases with different initial velocities are conducted; the experimental and numerical results are summarized in Table 6. Figure 14 shows a damage cone in the impact zone at 10 μs for the 553-m/s case that is consistent with Gellert’s model (Gellert et al., 2000). Similar to Section 3, the red color represents the distribution of the failure strains (${\epsilon }_1^{failure}$), which demonstrates the average damage evolution in the transverse direction.

Table 6
The results of experiment and simulation.

Target plate type	Initial Velocity	Residual Velocity		Error(%)
		Experiment	Simulation
Plain-woven laminates	299.4	193.8	188	-2.99
	381	270.3	271	0.26
	422.8	351.8	327	-7.05
	453	372.2	354	-4.89
	553	505	460	-8.91
	Ballistic limit	243	240	-1.23
Cross-plied laminates	246.6	62.6	68.1	8.79
	300.2	190.7	187.4	-1.73
	388.8	302.4	282.6	-6.55
	471.2	416.9	384.4	-7.80
	551.5	498.1	471.5	-5.34
	Ballistic limit	234.2	235	0.34

Figure 14
Damage cone in the impact zone.

The experimental and numerical BLVs for the plain-woven laminate are 243 m/s and 240 m/s, respectively. While for the cross-plied laminate, are 234.2m/s and 235m/s, respectively. Therefore, both of the BLVs, calculated by the developed model, agree well with the experiment. The maximum error is -8.91%, occurring in the 553-m/s case.

The relationship between the residual velocity (V _r ) and initial velocity (V _i ) of the projectile can be expressed as the known Lambert-Jonas penetration model (Lambert and Jonas, 1976)

$$V_r=\{\begin{array}{c}\begin{array}{cccc}\begin{array}{cc}\begin{array}{cc}\begin{array}{cc}\begin{array}{cc}0& \end{array}& \end{array}& \end{array}& \end{array}& V_r\le V_{bl}& & \end{array}\\ a{\left(V_r^p-V_{bl}^p\right)}^{\frac{1}{p}}\begin{array}{c}\begin{array}{cc}& V_r>V_{bl}\begin{array}{cc}\begin{array}{cc}& \end{array}& \end{array}\end{array}\end{array}\end{array}$$ (11)

where, the parameters, i.e. a, p are empirical constants.

Based on the Eq. (11), the residual velocity curves are drawn through fitting the experiment or simulation data (Figure 15). Those curves show some nonlinearity near the BLV; while a linear relationship is exhibited far from it. This shows that composite laminate is sensitive to the variations of initial velocities near the BLV. Moreover, the numerical and experimental curves, both of the different laminates, show a good correlation in the general trend; it means that the developed material model can be suitable for the plain-woven and cross-plied laminates. Furthermore, it is reasonable and feasible to capture the macroscopic mechanical properties of composite laminates through an orthotropic homogeneous model without considering their detailed microscopic structure.

Figure 15a
Comparison of numerical and experimental residual velocity curves: (a) plain woven laminates.

Figure 15b
Comparison of numerical and experimental residual velocity curves: (b) cross-plied laminates.

It is worth noting that the curves of the two cases are in good agreement at low velocity stage. However, deviation is emerged at high velocity stage and it will further enlarge with the increase of initial velocity. The above phenomenon may result from many factors such as model’s approximate errors, strain rate and plastic deformation, etc. To thoroughly understand this deviation, further research is needed and the material model itself has to be improved.

5 CONCLUSIONS

In present study, a computational macro-mechanical model coupling nonlinear EOS is developed and implemented in the nonlinear dynamic explicit finite element code LS-DYNA. The model’s theoretical framework including constitutive model, EOS, failure and erosion model has been introduced. The performance and computational capability of this developed model under high-velocity impact condition are verified and evaluated, and the following conclusions are drawn:

Inverse flyer plate simulations are conducted using the developed model and standard material model (MAT_2). The back-surface velocity history curves of composite laminate is calculated by this two different models and the results are compared with the experiment. Because of the introduction of a polynomial EOS, the developed model can better simulate the shock response of composite laminates than MAT_2 at high-velocity impact condition.

The ballistic performance of the developed model is assessed from qualitative and quantitative analysis. The classic ballistic impact problem, FPS impacting on Kevlar29/Phenolic laminate, is investigated. Using this model, the top-face cross-shaped failure, back-face bulging failure and three typical penetration processes have been simulated. It should be notice that the failure/damage predicted by this model only possesses the macro-mechanical meaning. BVL is also calculated by two numerical methods which agrees well with results obtained from experiment and theoretical analysis.

The computational capacity of this model, for Kevlar/Epoxy laminates with different architectures, i.e. plain-woven and cross-plied laminates, is further evaluated. The residual velocity curves and the damage cone of this two different laminates have been predicted and the results are compared with the experiment and theoretical model. The general trend of the numerical and experimental curves shows a good correlation and the damage cone in the impact zone is consistent with the theoretical model. The study shows that the developed model is a macro-mechanical model, which is suitable for different laminate structures, and the major ballistic impact behaviors of composite laminate with different structures can be captured by this model.

In general, the macroscopic parameters, such as ballistic limit, residual velocity, energy absorption and penetration/perforation behavior, can be accurately calculated by this model. This indicates that equating the whole composite laminates into an orthotropic homogeneous plate, to capture their macroscopic mechanical properties, is reasonable and feasible. In the future work, further research will be conducted to improve this model and node separation technique will be also introduced.

Acknowledgements

This research is supported by National Natural Science Foundation of China (11502010) and the Fundamental Research Funds for the Central Universities of China (YWF-14-YHXY-022).

References

• Anderson, CE, Cox, PA, Johnson, GR, et al., (1994). A constitutive formulation for anisotropic materials suitable for wave propagation computer programs-II. Computational Mechanics 15: 201-223.
• Ansari, MM, Chakrabarti, A., (2016). Impact behavior of FRP composite plate under low to hyper velocity impact. Composites Part B: Engineering 95: 462-474.
• Babu, MG, Velmurugan, R, Gupta, NK., (2007). Heavy mass projectile impact on thin and moderately thick unidirectional fiber/epoxy laminates. Latin American Journal of Solids & Structures 4: 247-265.
• Bandaru, AK, Chavan, VV, Ahmad, S, et al., (2016). Ballistic impact response of Kevlar® reinforced thermoplastic composite armors. International Journal of Impact Engineering 89: 1-13.
• Bandaru, AK, Vetiyatil, L, Ahmad, S., (2015). The effect of hybridization on the ballistic impact behavior of hybrid composite armors. Composites Part B: Engineering 76: 300-319.
• Billon, HH, Robinson, DJ., (2001). Models for the ballistic impact of fabric armour. International Journal of Impact Engineering 25: 411-422.
• Cheng, WL., (2003). High Velocity Impact of Thick Composites. International Journal Impact Engineering 2003: 167-187.
• Clegg, RA, White, DM, Riedel, W, et al., (2006). Hypervelocity impact damage prediction in composites: Part I-material model and characterisation. International Journal of Impact Engineering 33: 190-200.
• Deka, LJ, Bartus, SD, Vaidya, UK., (2008). Damage evolution and energy absorption of E-glass/polypropylene laminates subjected to ballistic impact. Journal of Materials Science 43: 4399-4410.
• Gellert, EP, Cimpoeru, SJ, Woodward, RL., (2000). A study of the effect of target thickness on the ballistic perforation of glass-fibre-reinforced plastic composites. International Journal of Impact Engineering 24: 445-456.
• Grujicic, M, Hariharan, A, Pandurangan, B, et al., (2012). Fiber-Level Modeling of Dynamic Strength of Kevlar® KM2 Ballistic Fabric. Journal of Materials Engineering and Performance 21: 1107-1119.
• Hayhurst, CJ, Hiermaier, SJ, Clegg, RA, et al., (1999). Development of material models for Nextel and Kevlar-epoxy for high pressures and strain rates. International Journal of Impact Engineering 23: 365-376.
• Hoof, JV. (1999). Modelling of impact induced delamination in composite materials, Ph.D. Thesis, Carleton University, Canada.
• Naik, NK, Doshi, AV., (2005). Ballistic Impact Behavior of Thick Composites: Analytical Formulation. AIAA Journal 43: 1525-1536.
• Ramadhan, AA, Abu Talib, AR, Mohd Rafie, AS, et al., (2013). High velocity impact response of Kevlar-29/epoxy and 6061-T6 aluminum laminated panels. Materials & Design 43: 307-321.
• Riedel, W, Nahme, H, White, DM, et al., (2006). Hypervelocity impact damage prediction in composites: Part II-experimental investigations and simulations. International Journal of Impact Engineering 33: 670-680.
• Silva, MAG, Cismaşiu, C, Chiorean, CG., (2005). Numerical simulation of ballistic impact on composite laminates. International Journal of Impact Engineering 31: 289-306.
• Tabiei, A, Ivanov, I., (2002). Computational micro-mechanical model of flexible woven fabric for finite element impact simulation. International Journal for Numerical Methods in Engineering 53: 1259-1276.
• Tham, CY, Tan, VBC, Lee, HP., (2008). Ballistic impact of a KEVLAR® helmet: Experiment and simulations. International Journal of Impact Engineering 35: 304-318.
• Wang, YB, Wang, XJ, Xiu-Zhang, HU, et al., (2005). Experiment study of ballistic resistance of Kevlar laminates. Engineering Mechanics 22: 76-81.
• Wen, HM., (2000). Predicting the penetration and perforation of FRP laminates struck normally by projectiles with different nose shapes. Composite Structures 49: 321-329.
• Xin, SH, Wen, HM., (2012). Numerical study on the perforation of fiber reinforced plastic laminates struck by high velocity projectiles. The Journal of Strain Analysis for Engineering Design 47: 513-523.
• Zeng, T, Fang, D, Lu, T., (2005). Dynamic crashing and impact energy absorption of 3D braided composite tubes. Materials Letters 59: 1491-1496.
• Zhang, GM, Batra, RC, Zheng, J., (2008). Effect of frame size, frame type, and clamping pressure on the ballistic performance of soft body armor. Composites Part B: Engineering 39: 476-489.
• Zhu, G, Goldsmith, W, Dharan, CKH., (1992). Penetration of laminated Kevlar by projectiles-I. Experimental investigation. International Journal of Solids & Structures 29: 399-420.
• Zu, XD, Huang, ZX, Zhai, W, et al., (2015). Performance of kevlar fibre-reinforced rubber composite armour against shaped-charge jet penetration. Latin American Journal of Solids & Structures 12: 507-519.

Publication Dates

• Publication in this collection
  Dec 2017

History

• Received
  08 Feb 2017
• Reviewed
  17 June 2017
• Accepted
  13 July 2017