FEM-SPH Numerical Simulation of Impact Loading on Floating Laminates

Abstract

The study of dynamic events such as impact and hydroelastic slamming is of great importance in determining the structural integrity of naval or maritime structures, particularly those made of composite materials. This topic has been investigated by numerous researchers using analytical, experimental, and numerical approaches. In this study, we propose using a hybrid numerical model combining smoothed-particle hydrodynamics (SPH) and the finite element method (FEM) to investigate the impact of external objects on floating laminates. The results show a good agreement with the available experimental data regarding the impact dynamic and some limitations in the damage determination.

1. Introduction

Marine structures constructed using composite laminates present a unique set of considerations regarding impact loads and slamming. Composite laminates offer several advantages, such as high strength-to-weight/stiffness-to-weight ratios and corrosion resistance, making them an attractive choice for marine applications. Rubino et al. wrote a comprehensive review about using composites in marine applications [1], while recently, Crupi et al. focused on green composites in marine engineering [2].

However, the response of composite laminates to dynamic loading, including impacts and slamming, requires careful analysis and design. When subjected to impact loads, they generally exhibit a complex response influenced by various factors. The energy absorbed during an impact event depends on the laminate thickness, fiber orientation, stacking sequence, and the impactor properties [3].

Delamination, the separation of adjacent layers in a laminate, is a common failure mode under impact loading. Therefore, analyzing the susceptibility of composite laminates to delamination and considering interlaminar fracture toughness properties are crucial aspects of impact load analysis. By incorporating design strategies such as interlaminar reinforcement or toughening layers, the resistance to delamination can be enhanced, ensuring the structural integrity of marine composite laminates [4].

Another essential consideration for marine structures made of composite laminates is their ability to absorb and dissipate impact energy. Composite laminates exhibit good energy absorption characteristics, which can help mitigate the effects of impact loads. Optimizing the stacking sequence and tailoring the laminate layup to distribute and dissipate impact energy throughout the structure are effective design approaches. In relatively recent years, Sutherland made a pretty extensive review, divided into four papers, of the impact testing of marine composite materials [5].

In the case of slamming, composite laminates face the challenge of rapid loading rates and high-pressure impulses. The response of composite laminates to slamming is influenced by factors such as the laminate configuration, thickness, stiffness, and damping properties. Slamming can induce significant structural vibrations, fatigue, and potential damage to the laminate and the underlying support structure. For some structures, the response to this kind of dynamic loading can be the base of the design process, as shown in [6].

A comprehensive analysis of composite laminates’ dynamic behavior is necessary to address these challenges. Advanced computational tools and modeling techniques are employed to simulate and predict the response of composite laminates under slamming conditions, such as the combined smoothed-particle hydrodynamics (SPH) and finite element method (FEM) [7,8,9,10,11].

The scope of this paper is to study a specific case of impact loading that may involve marine structure, i.e., the foreign object impact on a floating laminate.

Korobkin studied the problem [12], presenting asymptotic and numerical analyses of the unsteady hydroelastic behavior of a floating plate due to given external loads, considering some limit cases and specific hypotheses. In particular, concentrated loads are not considered in the paper.

To study the problem of landing airplanes on a floating pontoon, an analytical solution is obtained by [13] and compared with experimental drop test results.

In Ref. [14], the problem of a sinusoidal pressure pulse on a water-backed plate is solved in closed form. In addition, material non-linearity is introduced in [15], where the elastoplastic behavior of a floating plate under dynamic loading is studied by an FEM-boundary element method (BEM) model.

Following up Ref. [14], in Ref. [16], the membrane stiffening effect in fully constrained water-backed 2D beams is added, and a parametric study is conducted by varying the flexural stiffness. An FEM numerical validation of the analytical model is provided.

The same problem is faced experimentally in [17], where a digital image correlation (DIC) procedure is applied to a composite plate and particle image velocimetry (PIV) to the fluid below the impacted plate.

A set of valuable results of impact experiments in different conditions (air-backed and water-backed clamped and floating) is provided in [18], which will be used to validate the numerical model presented henceforth.

The anisotropy of the composite material under dynamic loading is considered by [19] using a moving element method (MEM)-FEM hybrid model, which, on the other side, is not dealing with foreign object impact loading.

Recently, a specific problem involving the impact of a blunt object on ice was solved analytically in [20], while the effect of water-backing as opposed to air-backing was experimentally evaluated in the case of carbon fiber reinforced plastics subjected to explosions [21].

As the literature review shows, when comparing the impact on a laminate to the case when it is floating on water, several essential differences arise due to the presence of the water environment.

One significant distinction lies in the load transmission. When a laminate is impacted in a dry environment, the load generated by the impact is primarily transmitted through the laminate structure itself. However, when the laminate is floating on water, the water acts as a medium for load transmission. The water can absorb and distribute the impact forces, potentially reducing the severity of the impact on the laminate.

The presence of water also introduces a damping effect. Water possesses damping characteristics that can help absorb and dissipate energy during an impact event. As a result, the water dampens the vibrations and reduces the peak stresses experienced by the laminate, thereby mitigating the damage caused by the impact.

Additionally, the interaction between the laminate and water introduces hydrodynamic effects. When an impact occurs on a floating laminate, water displacement, wave generation, and hydrodynamic forces come into play. These effects can influence the overall response of the laminate to the impact, affecting the distribution of stresses and the deformation patterns.

Another important consideration is the buoyancy and restraint provided by the water. The buoyant force exerted by the water on the floating laminate partially counteracts the downward force from the impact. This buoyancy can help reduce the deformation and prevent immediate failure, providing additional support and restraint to the laminate during the impact event.

These distinctions highlight the need to account for the water environment when assessing the impact on a floating laminate. Researchers must consider the unique dynamics and interactions introduced by the presence of water to evaluate the damage accurately, predict the response of the laminate, and design appropriate mitigation strategies.

By understanding these differences, naval engineers can make informed decisions when designing and optimizing floating composite laminates for marine structures, offshore platforms, and watercraft applications. Incorporating the effects of water load transmission, damping, hydrodynamics, and buoyancy into their analyses and design processes ensures the integrity, safety, and performance of the laminate in a floating environment subjected to impacts.

Recently, SPH-FEM hybrid models have been used in other types of applications in the field of impact between solid bodies, such as the simulation of fragmentation of a concrete slab under blast conditions [22,23] or the hypervelocity impact on shield for aerospace applications [24] or the impact of floating ice floes and marine structures [25]. Moreover, in [26,27], Zhou et al. used such a model to study the bird strike problem on composite laminates. In all these cases, SPH was used to simulate the fragmentation of the solid bodies. On the other hand, similar models were used to simulate fluid–solid interaction between tsunami waves, debris, and coastal structures [28].

This paper aims to explore the capabilities of an ANSYS LS-DYNA hybrid numerical model combining SPH and FEM to investigate the impact of external objects on floating laminates compared to other boundary conditions. In fact, the literature review showed that this method has never been applied to impact loading on composite laminates to simulate the effect of the presence of water.

2. Materials and Methods

This section describes the software utilized and the hybrid FEM-SPH numerical model. Ansys LS-DYNA is one of the most widely utilized software packages for explicit simulations globally. It can simulate how materials react when subjected to intense, short-term loads. With a wide range of elements, contact formulations, material models, and other configurable parameters, it enables the simulation of intricate models while maintaining control over all aspects of the problem.

SPH is a particle-based method for modeling impacts, explosions, and fluid-structure interaction (FSI) scenarios. It was developed as an alternative to the finite element method to overcome issues like tangled meshes in severe deformation problems. SPH allows for modeling complex behaviors involving free surfaces, material interfaces, and solid fragmentation. Ansys LS-DYNA offers an SPH solver that represents the system’s state using particles, each with its material properties. These particles interact with one another within a specified range determined by a weight function. Unlike traditional methods, SPH does not rely on a grid structure; instead, the particles form the computational framework to solve the governing equations. Ansys LS-DYNA integrates the SPH method with finite element and discrete element methods, expanding its capabilities to handle a wide range of multiphysics problems, including incompressible flows, heat conduction, high explosives, and high-velocity impacts.

2.1. Geometry Model and Discretization

The numerical method described in this article is built upon a foundation of previous studies and developments in the experimental field of fluid-structure interaction [18]. Different test conditions were applied during the experimental investigation, while clamped air-backed, water-backed, and floating conditions (see Figure 1) were investigated by numerical method. The experiment consisted of a series of free falls of an instrumented cylindrical impactor of $19.8$ mm in diameter in the water. The total mass of $3.6$ kg of the impactor combined with the drop heights allowed the authors to obtain the selected impact energy of 10 J or 20 J. The impact velocity was 4 m/s. To verify the influence of the water and the effect of the water–structure interaction as well as the boundary conditions and the load distribution, the water was only on the back side of the laminate, and the impactor was on the other in the air. The traditional air-backed impact condition was used to understand differences in deformation and absorbed energy.

The investigation was performed on carbon fiber laminate specimens $100\times 150$ mm, obtained by overlapping ten T700 carbon fabric plies (0/90), resulting in 3 mm nominal thicknesses. The tank full of water used for the experimental tests is 300 mm in diameter.

In this study, Ansys Design Modeler was used for the CAD model of the experimental setup. Air influence was neglected and not modeled. As for boundary conditions (Figure 2a), all tank nodes were fixed to avoid the water outlet, while four lateral surfaces of the specimen were fixed to simulate the experimental conditions. Gravity was applied to the structure, while velocity was imposed on the impactor’s rigid body. Fixed support refers to constraining the movement of nodes within a model. It is commonly used for boundary conditions that impede specific degrees of freedom.

Three different meshing methods were retained to create a discretized smooth particle hydrodynamics–finite element method (SPH-FEM) model, Figure 2b. Beginning with the most straightforward approach, the impactor and tank were discretized using explicit finite element tetrahedral elements with sizes ranging from 5 mm to 12 mm. Both the impactor and container were modeled as rigid bodies. To set up the mesh for the specimen, mechanical four-node quadratic finite elements with a size of 5 mm were utilized, while Ansys Composite Preprocessing (ACP Pre) was essential for the realization of laminated carbon fiber specimens of 3 mm. Finally, the SPH method, with a size of 6 mm, was used to model fluid and to consider its interaction with the surfaces of structures discretized by shell finite elements [29]. This value was chosen as a compromise between calculation time and accuracy, as shown in the Section 3.

2.2. The SPH and FEM Methods

In smoothed-particle hydrodynamics (SPH), the continuum is discretized into pseudo-particles, comparable to finite elements, but without shared nodes among them. These pseudo-particles can have different neighboring elements throughout the analysis, allowing for better handling of large deformations [30,31]. However, mesh-free Lagrangian methods like SPH require frequent searches for neighboring particles, which can be computationally demanding. To ensure the conservation laws of continuum mechanics in SPH, two principles are utilized: kernel and particle approximations. These principles help accurately represent the material’s behavior and interactions within the SPH simulation [32]. The SPH method relies on the interpolation technique that expresses a function in the form of its values in a set of disordered points [33,34,35].

The reader interested in the details of the SPH method may refer to Appendix A.

In Figure 3 is presented the finite element (FE) mesh composed of TET4 (4 nodes tetrahedral) and HEX 8 (8 nodes hexagonal), both first-order 3D element types. Since the behavior of the tank and impactor under the low-velocity impact is not interesting for this investigation, it was modeled with four-node tetrahedral elements to reduce calculation time and discretization effort. The specimen is modeled by 8-node hexagonal elements obtaining uniform mesh for impact simulation. The quality of a mesh cell can be quantified in different ways. ANSYS performs several geometrical checks on mesh elements to determine their quality. One of them is represented as mesh metrics. The element quality is based on the square root of the cube of the sum of the square of the edge lengths for 3D elements. The specimen is modeled using ACP, and there are two possibilities, exporting it as a shell or solid model. During this study, we noticed that only a solid model could give us all the necessary details after applying the described loads and boundaries. Unluckily, shell elements (that could reduce the calculation time) created in ACP and then exported in Ls dyna could not adequately represent the specimen behavior. Perhaps this inaccuracy will be improved in the new version of the software.

2.3. Contacts in a Coupled FEM-SPH

LS-DYNA’s contact modeling is straightforward for most users, as the predefined contact definitions are usually sufficient. However, expert users can take advantage of the extensive options and over 70 contact types offered by LS-DYNA to enhance their applications. This flexibility requires a deep understanding and expertise on the user’s part to utilize these advanced features effectively. Good contact condition is essential in load determination and subsequent stress calculation. In this study of numerical models of drop tests, shell–particle contact was adopted. The shell–particle contact is implemented as $\ast CONTACT\_AUTOMATIC\_NODES\_TO\_SURFACE$ [29], defining the contact-target penalty algorithm in which SPH particles are contact while the shell surface is a target. The default contact thickness between SPH and Lagrangian parts is zero when discussing SPH-FEM coupling. In some situations, defining a contact thickness is recommended to detect and avoid penetration between an SPH and Lagrangian parts. The same contact with the same principle was set up between SPH particles and the tank geometry. Other contact conditions implemented in this numerical method were $\ast CONTACT\_AUTOMATIC\_SURFACE\_TO\_SURFACE$ between the impactor and the specimen. All the contacts mentioned above are frictionless.

The expression $\ast CONTACT\_$ determines the contact type and underlying contact algorithm. There are differences in:

• One-way vs. two-way contact that we used. In one-way contacts, only the slave nodes are checked for penetration, while in two-way contacts masters and slaves are reversed, and contact is checked both ways.
• Automatic vs. non-automatic contact, where automatic in the name indicates an improved algorithm invoked to compute contact. Automatic contacts are non-oriented and can detect penetration coming from either side of the shell element. Always take into account shell thickness; it is more robust.
• Contact algorithm, penalty or constraint. In the case of the penalty algorithm, a slave node penetrates the master segment. This penetration is detected, and the penetration depth (DP) is calculated. This “penalty” force is applied on the slave node $Force=Stiffness\ast DP$. The effect is to project the node back up to the surface of the master segment. A reaction force is applied to the master segment nodes, so the total force on the master nodes equals the applied slave node force. Both a normal $FN$ and a friction force $Ff=\mu \ast FN$ penalty force are applied to the slave node.

The SOFT parameter is on an optional card of $\ast CONTACT$ that should be set in LS DYNA and controls the contact formulation for each surface (SOFT = 0, penalty algorithm). The contact stiffness scale factor is ramped from 0 to 1.0 over a finite time to avoid a sudden, large application of contact force (in our case, 0.1 with depth 5).

Figure 4 represents the contact between the impactor and the specimen as Automatic_surface_to_surface contact, while the contact pair between the specimen and the water and water with the tank is defined as Automatic_node_to_surface contact.

2.4. Material Mathematical Models

The use of lamina-level failure criteria in predicting damage initiation within laminated structures has been common, despite recognized limitations. These criteria are employed in conjunction with a degradation scheme that reduces material properties once failure is initiated. Due to computational constraints, capturing all the mechanisms of failure in a simulation is not currently feasible. Explicit finite element codes, such as LS-DYNA, solve equations of motion through direct integration using explicit methods. LS-DYNA incorporates various material models (MAT cards) for modeling progressive failure (PFM) and continuum damage mechanics (CDM). Some material models used include MAT22, MAT54/55, MAT58, and MAT162 [36,37,38,39,40]. Each material model utilizes a distinct modeling strategy, including different failure criteria, degradation schemes, material properties, and a set of parameters required for computation that may not have an immediate physical interpretation. These material models provide flexibility in representing and simulating the behavior of materials under dynamic loading conditions. The impacted specimen was laminated using ACP (Pre), while for material model *MAT054 was applied. The LS-DYNA material model *MAT054, known as enhanced composite damage, is commonly employed in dynamic simulations to simulate damage progression. It offers the advantage of requiring fewer experimental input parameters than damage mechanics-based material models. *MAT054 has widespread use in crash and impact simulations involving fabric composite materials, even though its design caters explicitly to orthotropic materials like unidirectional (UD) tape composite laminates. It is essential to evaluate the fundamental elastic behavior, failure initiation, and post-failure characteristics of *MAT054 and its suitability for modeling fabric composite material systems. Different articles deal with this material, but unfortunately, none gives precise information that should be applied to all the factors that exist in the settings of this mathematical model; usually, the information is sparse and even unfitting. Unfortunately, the Ansys LS-DYNA manual does not contain detailed definitions of the input parameters used. Each ply is defined by one integration point through the thickness, with the prescribed orientation. The material stress–strain curves in the fiber (1) and matrix (2) directions are input parameters for the *MAT054 material card. These curves are generated using the values of

Table 1. Material properties of T700 uni-directional tape.

Property	Symbol	LS-DYNA Parameter	Experimental Value
Density	$\rho$	RO	$1.52$ kg/m${}^2$
Young Modulus (1)	$E_1$	EA	127 GPa
Young Modulus (2)	$E_2$	EB/EC	$8.41$ GPa
Shear Modulus	$G_{12}$	GAB	$4.21$ GPa
Major Poisson’s ratio	${\nu }_{12}$	/	$0.309$
Minor Poisson’s ratio	${\nu }_{21}$	PRBA	$0.02049$
Strength in tension (1)	$F_1^{tu}$	XT	$2.2$ GPa
Strength in compression (1)	$F_1^{cu}$	XC	$1.47$ GPa
Strength in tension (2)	$F_2^{tu}$	YT	$48.9$ MPa
Strength in compression (2)	$F_2^{cu}$	YC	199 MPa
Shear strength	$F_{12}^{su}$	SC	154 MPa

, taken from [41,42]. This study shows how the material model *MAT054 can be used to generate a model that closely approximates the experiment and captures all its significant features. As far as sensitivity of this mathematical model and variation of all those parameters, this article is based on previous research presented in [43] and on the studies of other authors [37,44].

Appendix B lists the parameters of the Element *MAT054 used in this paper.

The cylinder impacting on water simplified only with the half sphere was modeled of steel with the nominal Young modulus E of 200 GPa, a Poisson’s ratio $\nu$ of $0.28$, and a density $\rho$ of 8500 kg/m${}^3$. Concerning the water, a density of 1000 kg/m${}^3$ was considered, while the viscosity coefficient $\eta$ was defined as ${10}^{-6}$ m${}^2$ s${}^{-1}$. A linear equation of state (EOS) was required to simulate fluid behavior. An EOS determines the hydrostatic behavior of the material by calculating pressure as a function of density and energy. Properties of EOS coefficient: C1 is assumed to be 1647 ms${}^{-1}$, while S1 is $1.921$ and S2 is zero, as suggested directly in Ansys LS-DYNA material database [45,46].

3. Results

The numerical simulation results are compared with the experimental results obtained by Lopresto et al. [18] to validate the model’s accuracy regarding force–displacement, force–time, and energy–time curves. Diagrams of air-backed and water-backed clamped and floating (water-backed free) conditions were compared with experimental ones.

The simulations required a fine-tuning phase of all the model parameters to achieve convergence which represents a general term that refers to the ability of a numerical method or algorithm to approach a stable and accurate solution when increasing the resolution or refining the model. This is particularly important in explicit simulations, especially impact analysis or highly dynamic phenomena, and it can require a balance between model refinement and computational efficiency to achieve effective and efficient simulation. In this context, convergence can be evaluated in several ways, and one is shown in Figure 5, where the effect or the SPH model resolution is evaluated in comparison with the experimental results.

Time convergence is the ability to obtain stable and accurate results as the integration time interval (time step) decreases. In the case of impact on a floating object, different time step integration was explored in the way to arrive at the final one of 5 ms. This time definition is long enough to capture the system’s dynamic behavior [47].

Spatial convergence, on the other hand, is the ability to obtain accurate results as the computational grid density is increased or the model is refined. In our case, different mesh sizes and particle sizes were tested in the way to create a model that can capture the impact effect and, on the other side, is not computationally demanding.

Results convergence indicates that the results obtained from the numerical simulation are stable and approach a definitive solution as the model’s resolution increases or the time step decreases. Finally, energy convergence is related to the conservation of energy to obtain physically correct results.

In Figure 6, Figure 7 and Figure 8, the relation between the displacement of the impacted node and the contact force is shown in the case of a 20 J impact. From a qualitative point of view, the shapes of the numerical and experimental curves closely match, particularly during the rebound phase in the force–displacement diagram. The displacement obtained from the simulation increases with the contact force and reaches the peak value at the same time as in the experiment. It rebounds to zero after reaching the peak value, consistent with the experimental results. The curvatures of the numerical and experimental curves agree correctly. In both experimental and numerical curves, it is possible to evaluate the effect of the presence of the water on the area of the cycle, which represents the energy absorbed by the system formed by the laminate and the water body.

Figure 9 compares the simulated force–time diagrams for the three different boundary conditions. It can be noticed that Air and WB clamped curves show a similar increase until they reach the peak force before declining during the rebounding stage. The peak contact force, the time at the instant of reaching the peak force, and the total impact period predicted by the numerical model all match well with the measurements obtained in the experiments. As expected, the floating conditions, even in the idealized case of central vertical and symmetric loading, show a more complex behavior with respect to the clamped case, with or without water-backing.

In

Table 2. Comparison between the experimental and numerical results.

Variable	Units	Experimental	Numerical	Difference
Air-backed clamped
Max Displacement	mm	$6.2$	$5.7$	$-8\%$
Max Force	kN	$7.0$	$7.5$	$7\%$
Absorbed Energy	J	$7.0$	$9.5$	$+35\%$
Water-backed clamped
Max Displacement	mm	$5.8$	$5.6$	$-8\%$
Max Force	kN	$6.8$	$7.2$	$5.7\%$
Absorbed Energy	J	$10.0$	$13.0$	$30\%$
Floating
Max Displacement	mm	$12.0$	$12.5$	$-4\%$
Max Force	kN	$4.2$	$4.5$	$7\%$
Absorbed Energy	J	$11.0$	$15.0$	$+36\%$

, the comparison between the main 20 J impact dynamic results is presented. While there is a good agreement as regards maximum force and maximum displacement, the absorbed energy, calculated from the area of the load–displacement curve, is systematically overestimated by the simulation. This fact is not dependent upon the SPH model’s presence because it also happens in the air-backed-clamped case. Consequently, it may be attributed to the progressive damage model of the *MAT054 element.

Figure 10 shows a timelapse of a 20 J impact on a floating laminate. In the first two frames, the laminate is bending, while in the third frame, at maximum displacement, the laminate is almost flat again, and the energy has been transferred to the water.

Figure 11a,b show the shear stress distribution and the maximum principal stress distribution in the case of an air-backed 10 J impact, while in Figure 11c,d, the same quantities are shown for a 20 J impact. It can be seen that the value and the distribution of the shear stresses are compatible with the delaminated area measured in the experimental tests, and the value of the maximum principal stress shows that there is no penetration for the 10 J impact. At the same time, it is foreseen in the case of 20 J impact. The other cases show similar patterns with slightly lower values of stresses for the WB condition and much lower values for the free condition, as expected from the maximum computed load.

4. Discussion

The proposed hybrid SPH-FEM model can effectively predict the impact dynamics in its qualitative evolution. To obtain this result, fine-tuning several numerical parameters through an extensive simulation campaign was necessary. The relevant values are listed in Appendix B. The computer run time is around three hours on an average up-to-date workstation making this model suitable for more complex impact events on floating laminate structures. As a drawback, determining damaged areas, particularly regarding delamination, has to be performed in a post-processing phase, which is not automatic. Moreover, delamination is not considered in the *MAT054 progressive damage model, and this could be an issue for repeated low-energy, barely visible impacts, for example, at the waterline level of composite hulls.

5. Conclusions

In the present paper, a hybrid SPH-FEM model has been used to simulate low-velocity impact tests in different conditions to understand the critical issues during the navigation of a ship made of composite laminates. The results show that:

• The fluid–structure interaction significantly influences the dynamic response of a composite panel in good agreement with the literature’s experimental data;
• The absorbed energy is higher in the floating and water-backed clamped conditions due to the quota transferred to the water, thus resulting in lower damage to the laminates with respect to the air-backed clamped condition;
• Damage calculation by using *MAT054 element is not straightforward and does not include delamination.

In the future, other types of finite elements will be considered, and this model will be applied to more complex cases or geometries for which experimental data are not currently available or difficult to obtain, such as non-symmetric or angled impacts and simply supported or mixed boundary conditions.