**1. Introduction**

With the aim of damage reduction in person subjected to ballistic impact, polymeric, carbon, and glass fibers are commonly used to develop protective systems. However, during last decades it has been found that a combination of high strength yarns in different directions generates flexible woven fabrics which are light and highly resistant at the same time [1]. Therefore, aramid fibers are one of the most used protection materials nowadays [2], with a growing trend in the industry.

Experimental studies have proved that the ballistic performance of aramid fabrics depends on many factors, such as the projectile geometry [3,4], the impact velocity [3,5–7], the friction between yarns [8–10], the woven structure of the fabric [6,11–13] and material properties [14]. However, most research has found that each property individually does not control the ballistic performance [15,16], rather the combined effect of all of these properties play a key role.

New studies to predict and to reduce the damage are constantly under study [17,18]. The use of surrogate models is common in optimal design problems to approximate the objective functions, as neuronal networks [19,20] or genetic algorithms [21]. However, a big dataset is needed in most cases to obtain good predicted values what can be expensive.

A useful tool for obtaining improved understanding of the behavior of aramid materials under impact loads is the use of numerical models. The ballistic limit, V50, the deformed shape of the woven fabric and the effect of the internal interlayer friction can be studied using the finite element method. However, the geometry of the fabric structure and the failure mechanisms of the yarns make the modelling of 3D fabrics complex [22–30]. To decrease the complexity of the finite element model, new simple models using shell elements [31–35] and truss elements [36] have been successfully implemented.

Duan et al. [22], Rao et al. [23] and Grujicic et al. [29,30] developed 3D numerical models with solid elements representing yarns as homogenous continua. These studies analyzed mainly the friction coefficient between the yarns and the clamping conditions. Modelling results showed that the fabric boundary condition is a primary factor that influenced the friction effect, when only two edges are clamped, fabric reduces the residual velocity of the projectile and absorbs energy more effectively. Fabrics with high stiffness decelerate the projectile relatively rapidly meanwhile fabrics with high-strength yarns need more time to initiate the failure. The material which combines both properties is the most favorable of all the examined ones under the imposed boundary conditions considered in the study.

Despite the accuracy of 3D models, some researchers observed that they have two limitations. The first is that these models do not consider the statistical variability in the fabric geometry and material properties. The second one is the high computational cost due to the realistic representation of the fabric that requires a fine mesh with many elements.

To address the first problem and improve the prediction of the V<sup>50</sup> velocities, Nilakantan et al. [24,26,37] implemented a framework which incorporated the inherent statistical variability in the system. This framework can be used to calculate the probabilistic velocity response (PVR) curve for projectiles with different size and shape and for different clamped conditions. These studies analyzed the projectile geometry [26], the impact location sensitivity [27] and the fabric clamping conditions [28,38]. The findings showed that projectiles with a small impact area resulted in a higher residual velocity when the impact occurs in the gap between the yarns compared to when it occurs at the yarn crossover junction. The large cylindrical projectile, which had the flattest and the largest impact face, showed no sensitivity to impact location. The study also remarked that a circular shaped projectile shows lower dependency on the residual velocity with the projectile impact location than other configurations. This is mainly due to the fact that this configuration is not sensitive to in-plane fabric rotations, giving it a slight advantage over the other clamping configurations.

To address the second problem related to the computational efficiency, researchers presented different solutions. One option was the use of shell elements instead of solid elements [32–35]. Ha–Minh et al. [32,33] studied the effect of the number of elements used in the model (four and eight elements by yarn) and showed, in terms of velocity evolution, similar results in both cases. This is a significant conclusion because the computation time of the model with 8 elements is double than the model with 4 elements.

Chu et al. [35] studied mechanical properties of the yarns. It was concluded that the yarn density does not affect significantly the ballistic performance of the fabric. As opposed to this, it was observed that a high value of longitudinal Young's modulus produces a faster deceleration of the projectile.

The second solution to reduce the computational cost is to modify the original 3D model, generating a hybrid model, which is comprised of zones with different modelling resolutions and different finite element formulation, all coupled together with impedance matching interfaces [31,39–42]. Barauskas et al. [31] developed a model for a plain-woven single-ply. The fabric model presented three different zones: a zone close to the impact area including failure modelling, a second zone where the woven structure did not undergo failure, and a last zone (far from the impact area) implemented as orthotropic membrane. During this research, it was proved that it is not easy to validate this kind of approach due to the many different levels involved in the process.

The original model of Nilakantan et al. [26], was modified in later studies [39–41]. New models which combine solid and shell elements were tested under different distributions of both meso-level and macro-level in the fabric producing a low cost computational model.

Finally, Ha-Minh et al. [42] carried out a study comparing three different hybrid models, that was a modification of their initial macroscopic-model. The numerical analysis showed that the difference among the models was negligible for results such as the evolution of projectile velocity and the global impact fabric behavior.

With the objective of studying the friction coefficient in aramid woven, Das et al. [36] implemented a model with truss elements. The study analyzed the influence of the friction coefficient and the clamping conditions. It was reported that there is a limit for the friction between yarns since the shape of the projectile is negligible and affects in a negative way to the energy dissipation and the failure mechanism.

In the present study, a simplified finite element model is developed using truss elements with the aim of establishing the advantage that 1D element models present over solid element models. A validation analysis proved that the simplified model is efficient in terms of computational cost and accuracy. A parametric study has been carried out using the numerical model, with specific focus on the effect of the impact angle, a parameter that is not usually analyzed in the literature. Results are used to generate a response surface diagram able to predict the ballistic limit for low impact velocities. This methodology will be a useful tool for rapid impact response analysis in the industry.

#### **2. Theoretical Study**

The study of the physical mechanisms by which impact pressure waves reach the chest and cause injury are continuously under study in the research community. It is a goal to determine how the impact pressure waves are transmitted to the thorax or the brain in order to implement effective preventive measures and reduce the exposure risks [43,44].

The energy absorption rate of aramid fabric during impact depends on many variables. Among all of them, it should be mentioned the fiber modulus and the material modulus, which is related to the ability to brake the projectile [22,45].

The breaking energy of the yarn is determined by its characteristics, such as tensile strength, elongation or modulus. These factors affect the transmitting velocity of the stress wave that is generated by the impact of the projectile. The wave must be dispersed rapidly, and the breaking energy must be as high as possible to increase the impact resistance and the ballistic limit (velocity at which the projectile passes through the material) [46].

The characteristics of the material, such as the interlacing geometry of the fibers, the thickness and the number of layers, also affect the distance at which the impact disturbance will have moved in a given time. They also influence the breaking energy of the material (ability to resist breakage due to an external force), which depends on the tensile strength and elongation of the material.

#### **3. Numerical Model**

The numerical models implemented in this study were developed using the finite element explicit code ABAQUS. The following sections describe the methodology used to generate both models.

#### *3.1. Projectile*

All projectiles were implemented as rigid solid objects to reduce the computational time [28]. Four geometric shapes (blunt, hemispherical, conical and truncated conical) was defined to keep the mass constant to 7.5 g for all cases (see Figure 1). The friction coefficient between projectile and material was established to be 0.22, taken from the literature [27]. As a rigid 3D object contained in a 3D space, a projectile has six degrees of freedom for motion. All rotation degrees have been restrained.

The projectile can only traverse in the direction of the impact. The projectile velocity normal to the fabric is specified before each simulation and it ranged between 30–130 m/s based on experimental tests carried out in the same material by Yu et al. [47]. normal to the fabric is specified before each simulation and it ranged between 30–130 m/s based on experimental tests carried out in the same material by Yu et al. [47].

*Materials* **2019**, *12*, x FOR PEER REVIEW 4 of 16

**Figure 1.** Dimensions of the projectiles: (**a**) blunt projectile, (**b**) conical point projectile, (**c**) hemispherical projectile and (**d**) truncated conical point projectile. **Figure 1.** Dimensions of the projectiles: (**a**) blunt projectile, (**b**) conical point projectile, (**c**) hemispherical projectile and (**d**) truncated conical point projectile.

#### *3.2. Specimen 3.2. Specimen*

The specimen was impacted by the projectile in the center. It was fully clamped, resulting an effective area of 101 × 101 mm2 (Figure 2a). Each yarn had density of 1440 kg/m3 and was modelled as a linear-elastic material until failure to capture the behavior of the woven fabric subjected to ballistic impact (Figure 2b) [23,24]. When the ultimate strength was reached, the element was assumed to be damaged and removed from the computational domain. The mechanical properties of the aramid fabric were taken from the literature [36], where no apparent plastic deformation before fracture was reported. Table 1 shows the most important parameters of the yarns in both the directions. The static and dynamic friction coefficient values between yarns are 0.186 and 0.17, respectively. Both coefficients were obtained from semi-analytical model based on yarn pull-out tests The specimen was impacted by the projectile in the center. It was fully clamped, resulting an <sup>e</sup>ffective area of 101 <sup>×</sup> 101 mm<sup>2</sup> (Figure 2a). Each yarn had density of 1440 kg/m<sup>3</sup> and was modelled as a linear-elastic material until failure to capture the behavior of the woven fabric subjected to ballistic impact (Figure 2b) [23,24]. When the ultimate strength was reached, the element was assumed to be damaged and removed from the computational domain. The mechanical properties of the aramid fabric were taken from the literature [36], where no apparent plastic deformation before fracture was reported. Table 1 shows the most important parameters of the yarns in both the directions. The static and dynamic friction coefficient values between yarns are 0.186 and 0.17, respectively. Both coefficients were obtained from semi-analytical model based on yarn pull-out tests [36]. *Materials* **2019**, *12*, x FOR PEER REVIEW 5 of 16 of an approximate cross section of 0.064 mm2 (Figure 3b). Dimensions of the yarns are also given in Figure 3.

On the other hand, the simplified model presents the yarns as 3D one-dimensional truss elements, which use linear interpolation for position and displacement, and have a constant stress. **Figure 2.** Fabric impact test setup with angular impact of the projectiles: (**a**) front view and (**b**) side view. **Figure 2.** Fabric impact test setup with angular impact of the projectiles: (**a**) front view and (**b**) side view.

In this study, the element T3D2 (two-node linear displacement) was chosen. It had three degrees of freedom by node. Truss elements are long, slender structural members that can transmit only axial force and do not transmit moments [48]. The cross section was considered to be made up of trusses

**Yarn Direction 2a (mm) 2b (mm) Warp** 0.75 0.11 **Fill** 0.63 0.12

**4. Validation and Comparison between Models** 

was dissipated.

(**a**) (**b**) **Figure 3.** Yarn geometrical parameters. (**a**) Solid elements model and (**b**) Truss elements model.

The study encompasses evaluation of computational time, projectile geometry effects and ballistic curve prediction (initial and residual velocities). For the flat nose projectile case, with and impact angle of 0°, the ballistic curves obtained with the 1D and 3D models proposed, and the

In the Figure 4, three different regions, typical from this kind of curves, are observed:

1. In the first one (*Vi* < 44 m/s), the negative residual velocities indicate that the projectile rebounds due to the impact. This phase was dominated by fiber elongation. The yarns returned to the original configuration and released most of the elastic energy stored, imparting the projectile almost its original velocity towards the opposite direction. Thus, only a small amount of energy

2. In the second phase (44 m/s < *Vi* < 52 m/s) an abrupt increment from a negative to a positive value of the residual velocity was observed. This corresponds to the rupture failure of the fabric. Some of the yarns surpassed their tensile yield stress and collapsed, dissipating energy in the process.

corresponding experimental values published in [36], are presented in Figure 4.


**Table 1.** Material properties of the Kevlar yarns [36].

*Materials* **2019**, *12*, x FOR PEER REVIEW 5 of 16

of an approximate cross section of 0.064 mm2 (Figure 3b). Dimensions of the yarns are also given in

The difference between both numerical models proposed, presented in Figure 3, is the representation of the yarns. In the complete full model case, the yarn dimensions were taken from [36] and were modelled with 3D solid elements (Figure 3a). The 3D elements were the standard volume elements of Abaqus. This element can be composed of a single homogeneous or heterogeneous material; they are more accurate if not distorted, particularly for quadrilaterals and hexahedra. In the present study the element used was C3D8R (eight-node linear brick) with reduced integration and hourglass control. It has three degrees of freedom by node [48]. (**a**) (**b**) **Figure 2.** Fabric impact test setup with angular impact of the projectiles: (**a**) front view and (**b**) side view.

**Figure 3.** Yarn geometrical parameters. (**a**) Solid elements model and (**b**) Truss elements model. **Figure 3.** Yarn geometrical parameters. (**a**) Solid elements model and (**b**) Truss elements model.

**4. Validation and Comparison between Models**  The study encompasses evaluation of computational time, projectile geometry effects and ballistic curve prediction (initial and residual velocities). For the flat nose projectile case, with and impact angle of 0°, the ballistic curves obtained with the 1D and 3D models proposed, and the corresponding experimental values published in [36], are presented in Figure 4. In the Figure 4, three different regions, typical from this kind of curves, are observed: On the other hand, the simplified model presents the yarns as 3D one-dimensional truss elements, which use linear interpolation for position and displacement, and have a constant stress. In this study, the element T3D2 (two-node linear displacement) was chosen. It had three degrees of freedom by node. Truss elements are long, slender structural members that can transmit only axial force and do not transmit moments [48]. The cross section was considered to be made up of trusses of an approximate cross section of 0.064 mm<sup>2</sup> (Figure 3b). Dimensions of the yarns are also given in Figure 3.

#### 1. In the first one (*Vi* < 44 m/s), the negative residual velocities indicate that the projectile rebounds due to the impact. This phase was dominated by fiber elongation. The yarns returned to the **4. Validation and Comparison between Models**

original configuration and released most of the elastic energy stored, imparting the projectile almost its original velocity towards the opposite direction. Thus, only a small amount of energy was dissipated. 2. In the second phase (44 m/s < *Vi* < 52 m/s) an abrupt increment from a negative to a positive value of the residual velocity was observed. This corresponds to the rupture failure of the fabric. Some The study encompasses evaluation of computational time, projectile geometry effects and ballistic curve prediction (initial and residual velocities). For the flat nose projectile case, with and impact angle of 0◦ , the ballistic curves obtained with the 1D and 3D models proposed, and the corresponding experimental values published in [36], are presented in Figure 4.

of the yarns surpassed their tensile yield stress and collapsed, dissipating energy in the process. In the Figure 4, three different regions, typical from this kind of curves, are observed:


6b.

Some of the yarns surpassed their tensile yield stress and collapsed, dissipating energy in the process. The number of breaking (failed) yarns increased rapidly in a short range of velocities, allowing the projectile to pierce or penetrate the fabric. The number of breaking (failed) yarns increased rapidly in a short range of velocities, allowing the projectile to pierce or penetrate the fabric. *Materials* **2019**, *12*, x FOR PEER REVIEW 6 of 16

*Materials* **2019**, *12*, x FOR PEER REVIEW 6 of 16

3. If the initial velocity is further increased (*V<sup>i</sup>* > 52 m/s), the projectile always pierced the fabric, but the initial velocity was reduced by the dissipated energy. That is because the yarns first suffered an elongation, and then they failed by rupture because of the impact, absorbing some portion of the impact energy from the projectile. 3. If the initial velocity is further increased (*Vi* > 52 m/s), the projectile always pierced the fabric, but the initial velocity was reduced by the dissipated energy. That is because the yarns first suffered an elongation, and then they failed by rupture because of the impact, absorbing some portion of the impact energy from the projectile. The number of breaking (failed) yarns increased rapidly in a short range of velocities, allowing the projectile to pierce or penetrate the fabric. 3. If the initial velocity is further increased (*Vi* > 52 m/s), the projectile always pierced the fabric, but the initial velocity was reduced by the dissipated energy. That is because the yarns first

The results obtained accurately reproduced the experimental trends above and below the ballistic limit. Here, negative residual velocities mean projectile rebound (impact velocity below the ballistic limit). In particular, the ballistic limit error was close to 4%, as can be seen in Table 2. The results obtained accurately reproduced the experimental trends above and below the ballistic limit. Here, negative residual velocities mean projectile rebound (impact velocity below the ballistic limit). In particular, the ballistic limit error was close to 4%, as can be seen in Table 2. suffered an elongation, and then they failed by rupture because of the impact, absorbing some portion of the impact energy from the projectile. The results obtained accurately reproduced the experimental trends above and below the

ballistic limit. Here, negative residual velocities mean projectile rebound (impact velocity below the

**Figure 4.** Calibrated ballistic curves for flat projectile. Experimental data taken from [36]. **Figure 4.** Calibrated ballistic curves for flat projectile. Experimental data taken from [36].

**Table 2.** Ballistic limit results for the truss elements model. **Table 2.** Ballistic limit results for the truss elements model. **Figure 4.** Calibrated ballistic curves for flat projectile. Experimental data taken from [36]. **Table 2.** Ballistic limit results for the truss elements model.


Once the model was calibrated for the flat nose projectile, the hemispherical point projectile case was analyzed. Is this case, the general trends were also quite well predicted, see Figure 5. The three different regions described for the Figure 4 were also observed here. Although the error for the ballistic limit derived from the change of projectile geometry was close to −7% (Table 2), it can be considered acceptable according to the simplifications of the model and the level of engineering accuracy desired. Figure 6 shows the stressed fibers (typical cross section) at different instants of the simulation (Figure 6a). The fiber separation and fiber breakage during perforation is shown in Figure Once the model was calibrated for the flat nose projectile, the hemispherical point projectile case was analyzed. Is this case, the general trends were also quite well predicted, see Figure 5. The three different regions described for the Figure 4 were also observed here. Although the error for the ballistic limit derived from the change of projectile geometry was close to −7% (Table 2), it can be considered acceptable according to the simplifications of the model and the level of engineering accuracy desired. Figure 6 shows the stressed fibers (typical cross section) at different instants of the simulation (Figure 6a). The fiber separation and fiber breakage during perforation is shown in Figure 6b. Once the model was calibrated for the flat nose projectile, the hemispherical point projectile case was analyzed. Is this case, the general trends were also quite well predicted, see Figure 5. The three different regions described for the Figure 4 were also observed here. Although the error for the ballistic limit derived from the change of projectile geometry was close to −7% (Table 2), it can be considered acceptable according to the simplifications of the model and the level of engineering accuracy desired. Figure 6 shows the stressed fibers (typical cross section) at different instants of the simulation (Figure 6a). The fiber separation and fiber breakage during perforation is shown in Figure 6b.

**Figure 5.** Validated ballistic curves for hemispherical point projectile. Experimental data taken from [36]. **Figure 5.** Validated ballistic curves for hemispherical point projectile. Experimental data taken from [36].

**Figure 6.** (**a**) Penetration of the hemispherical point projectile into a fully clamped layer at different time instants and (**b**) detail of the yarns failure. **Figure 6.** (**a**) Penetration of the hemispherical point projectile into a fully clamped layer at different time instants and (**b**) detail of the yarns failure.

In general, it can be observed that both 1D and 3D element models present reasonable accuracy when predicting the residual velocity of the projectile in terms of magnitude and trending. It is possible to observe the overestimation of residual velocity predicted with the solid model; in particular, it is slightly larger than that obtained with the simplified model, which involves that the solid model is conservative. In general, it can be observed that both 1D and 3D element models present reasonable accuracy when predicting the residual velocity of the projectile in terms of magnitude and trending. It is possible to observe the overestimation of residual velocity predicted with the solid model; in particular, it is slightly larger than that obtained with the simplified model, which involves that the solid model is conservative.

Despite an improved velocity estimation by the 3D model, the results obtained by the simplified model (which is around 99% faster, Table 3) lead to consideration it as an efficient tool in terms of computational requirement and accuracy level. Despite an improved velocity estimation by the 3D model, the results obtained by the simplified model (which is around 99% faster, Table 3) lead to consideration it as an efficient tool in terms of computational requirement and accuracy level.


**Hemispherical projectile** 270 s 29,058 s 99.07%

low computational time. Results are discussed in the following sections.

**Table 3.** Comparison of computational time for the case of 60 m/s between 1D and 3D models. **Table 3.** Comparison of computational time for the case of 60 m/s between 1D and 3D models.
