1 Introduction

Over the last decades, several industrial sectors including automotive [1], aerospace [2] and energy [3, 4] have focused their attention on composite materials, which improve the general performance of structures due to their high in-plane mechanical properties and low weight [5]. However, these materials show weak resistance against dynamic loading along the through the thickness direction due to their intrinsic stratified nature. In particular, under Low-Velocity Impact (LVI) conditions such as tool drop, bird strikes and debris impact, the material is prone to generate Barely Visible Impact Damage (BVID) that can locally weaken the mechanical properties of the structure and cause its eventual and catastrophic failure [6].

In order to overcome this limitation, several solutions have been proposed in literature. One is modification of the single component that can be carried out by reinforcing fibre-matrix interface bonding [7, 8] or modifying the toughness of the matrix material [9] in order to improve the mechanical performances as failure strength and stiffness of the composite. Another solution is the modification of the inner structure of traditional composites, for example introducing a discontinuity pattern within the composite’s body, to increase absorbed energy and toughness due to their brick and mortar configuration [10].

Hybridisation is also an excellent solution to increase the impact properties of composite materials and it can be carried out by the introduction of one or more additional reinforcing materials within the polymeric matrix. In this context, Safri et al.. in their extensive review investigated several hybrid reinforcement configurations including synthetic-synthetic (i.e. glass-carbon), synthetic-natural (i.e. carbon/flax) and natural-natural fibres (i.e. flax-hemp) assessing that synthetic-synthetic hybrid reinforcement is the one showing the best mechanical performances, but also that the synthetic-natural has an excellent potential for its application in advanced sectors. As an example of synthetic-synthetic hybridisation, Hosur et al. investigated the effect on impact resistance using carbon and glass twill fabrics as reinforcement for epoxy resin systems under low velocity impact conditions. Results showed that hybrid composites are able to tolerate a higher contact force in comparison with traditional composites reporting also an increase of the structure’s stiffness. Following the same concept, also natural fibres can be used in the hybridisation process to improve mechanical properties and increase its bio-friendless. An example of synthetic-natural hybridisation is illustrated by Sarasini et al.. [11] that used a carbon/flax reinforcement to fabricate a hybrid composite that was tests under low impact velocity conditions. Results showed a positive effect on the structure performance with better impact absorption, a reduced BVID generation and higher post-impact mechanical properties when compared to traditional tradition carbon fibre composites.

Other research works instead focussed their attention in reinforcing the material by inserting a secondary reinforcement (hybridisation) such as copper [12] or Shape-Memory Alloy (SMA) wires [13, 14] within the composite lamination sequence, that increases impact resistance and, simultaneously, enables additional non-structural features such as Structural Health Monitoring (SHM) via strain monitoring, or damage evaluation via in-situ embedded thermography [15, 16].

In this context, another approach to enhance the impact resistance of composite materials is the use of non-Newtonian fluids as they allow for additional energy absorption mechanisms that are activated before the brittle fracture of the fibrous reinforcement. In contrast to Newtonian Fluids which show a constant viscosity in function of the load applied, non-Newtonian fluids are characterised by a viscosity which is dependent on the solicitation applied to the fluid itself. Shear Thickening Fluids (STFs) are a specific category of Non-Newtonian fluids, whose viscosity depends on the strain rate applied to the material: the higher is the magnitude of the applied load, the higher is the fluid’s viscosity[17]. These fluids are generally colloidal solutions in which solid insoluble particles are dispersed and suspended in a solvent. Several theories have been developed in order to explain the shear-thickening phenomenon. Firstly, Hoffman [18] observed that when a shear rate is applied on the STF colloidal solution, an order-disorder transition occurs and hydroclusters (solid particles agglomerates) are formed The presence of these agglomerates increases the solution’s hydrodynamic resistance against motion, leading to a subsequent increase of viscosity and, consequently, more energy is absorbed during the displacement of the fluid. Based on this work, other authors [19] investigated on the use of STF and its rheology evaluating a similar behaviour. Afterwards, the same Hoffman [20] improved his rheological theory introducing the contact force model and explaining that shear thinking phenomenon is not driven by order-disorder transition, but it is initiated by a hydrodynamic instability and amplified by contact forces generated during their particles interaction. In particular, when a low shear rate is applied on the colloidal solution, the particles organise themselves in ordered layers and the fluid is able to easily flow between them. Afterwards, increasing the strain rate, the shear-thickening is triggered, and the hydrodynamic instability moves the particles away from these ordered layers. When this happens, the particles interact via direct contact and/or clustering, jamming the space between the layers and forming obstacles for the fluid flow (hydroclusters). Increasing the magnitude of the applied shear rate, the hydroclusters start to be directly in contact despite the repulsive presence of the fluid forming contact force networks and increasing the size of the obstacles. This leads to observe a considerable increase of the viscosity of the solution since the fluid has to flow around these massive obstacles and thus a higher energy is required to obtain the same fluid displacement at lower viscosities. The validity of this theory has been confirmed in several works [21, 22] in which the phenomena of shear-thickening and jamming are analysed, discusses and applied in protective applications.

Hence, by embedding a layer of STF within the laminate’s body it is possible to exploit this unique property to dissipate more energy during an impact event and reduce the extent of the internal damage. Several researchers have studied STF coupled with different structures demonstrating its ability to improve vibrational and impact properties. For example, Fischer et al.. [23] tuned the stiffness and damping properties by embedding STF as core in a polyvinyl chloride (PVC) sandwich and controlling the volume fraction of particles dispersed within the fluid. Testing the structure under dynamic loading, relevant results were found as the presence of STF was able to use part of the applied stress on the structure to activate the shear-thickening effect increasing the viscosity and, consequently, the damping ability of structure. Another example of impact resistance improvement is reported by Wang et al., [24] who studied the impact response of a sandwich structure in which shear-thickening fluid was used as core (contained within a frame) for an aluminium sandwich. High absorbed impact energy and a reduced contact force were observed for STF-aluminium sandwich in comparison to other typologies of cores due to the synergic action of STF in damping and distributing stress on a wider area for energies between 10 and 40 J. Gürgen et al.. [25] instead investigated several typologies of structures using STF as bulk or impregnating element in order to demonstrate its effectiveness in protection under stab conditions. Results showed a significant improvement of the impact response in particular in case when different ply orientations are used. Majumdar et al.. [26] instead designed an optimised process to create STF-impregnated Kevlar fabrics reporting an significative increase of impact energy when compared to neat Kevlar fabrics under ballistic conditions. An advanced use of STF as reinforcing component is shown by Gürgen et al..[27] that used multiphase STF to reinforce a Kevlar fabric. The multiphase STF was prepared using two typologies of particles and introduced into the Kevlar fabrics. Ballistic tests were then carried out on this composite textile demonstrating superior performance in comparison of single-phase STF and neat fabric, but an increase of weight is also observed.

A different approach was introduced by Pinto et al.. [28] who used the STF directly as a reinforcement layer in the lamination sequence of layered composite materials with a negligible weight increase. STFs at different concentrations were successfully integrated into the laminates and results showed that for high impact energy (40J), the hybrid CFRP laminate was able to absorb 50% more of impact energy than traditional laminate reducing the damaged area by 30%.

Although the benefits of STF have been widely investigated from experimental aspects of these works, there is a need for implementing a numerical analysis of the STF/CFRP system that may allow the prediction of its mechanical response and the optimisation of the STF characteristics. This lack of work is due to the complexity of simulating the STF non-linear behaviour coupled with a structural analysis. Examples of the STF simulation for impact resistance improvement can be found in literature [29, 30] for impregnated dry fabric, and simulated the presence of the fluid by modifying the contact friction coefficient. Following this approach, Lee et al.. [31] developed a numerical model to describe the STF-impregnated Kevlar fabric behaviour under ballistic impact analysing the fundamental mechanisms of STF-Kevlar interaction and its damage suppression nature. Another example in simulating STF/fabric systems using the frictional approach is illustrated by Gürgen [32] that compared the effect of single phase and double phase STF in reinforcing textile subjected to high velocity impact conditions reporting good correlation between numerical and experimental data especially when two typologies of particles are used to create the STF. Rabb et al.. [30] also followed the same STF modelling concept but they introduced a hybrid particles-element numerical approach to simulate the beneficial effect of STF on impact resistance impregnating Kevlar fabrics under ballistic conditions. This numerical model shows good agreement with experimental results, but it fails in predicting ballistic limit configuration.

Although these numerical works are efficient to predict the STF effect in impregnating dry fabrics, their numerical approach is not suitable when the fluid is used as reinforcing layer in composite structures. Indeed, in the case of dry fabrics, the STF contribute can be correctly modelled via the modification of the algorithm contact adjusting the friction parameters between the warp and weft of the fabric. For composite structures, instead, this approach cannot be used since the fluid is embedded within the structure and actively bears loads during the dynamic event. Thus, an accurate and efficient fluid-structure interaction analysis is required in order to evaluate the correct stress distribution in the entire structure.

Based on this premise, the main focus of this work is the development of an innovative numerical FEA model able to predict the mechanical behaviour of STF/CFRP hybrid laminate under LVI conditions using a fluid-structure interactions algorithm. An explicit 3D FEA code was developed to evaluate the impact response of the different mesh configurations and optimise the damage suppression ability. 3D structured mesh configurations were developed by locating STF at three different depths along the thickness of the laminate. Furthermore, impact output of STF/CFRP hybrid laminate were compared to traditional CFRP results evaluating damaged area reduction and mechanical properties enhancement. An experimental impact campaign was carried out in order to collect impact data for model validation. In this paper, the description of the shear thickening behaviour of STF is reported in Sect. 3 while the manufacturing procedure and the experimental setup used during the impact campaign are reported in Sect. 4. Section 5 illustrates the FE modelling approach used to develop the simulation and optimise the position of STF layer. Validation of the model and numerical results are reported and discussed in Sect. 6 whilst conclusions are illustrated in Sect. 7.

2 Shear-Thickening Fluid (STF) Properties

Viscosity is mainly defined as the fluid resistance to gradual tensile or shear deformations. The microscopic causes of viscosity can be found in the friction among fluid molecules that needs to be overcome to start the flow. Using viscosity coefficient \( \mu \), it is possible to correlate shear stress (τ) with shear rate deriving shear displacement (γ) in time (t) using the Eq. 1:

$$ \tau =\mu \frac{\partial \gamma }{\partial t} $$
(1)

Using the definition of Manley, Newtonian Fluids are those fluids that show an absolute viscosity independent from the velocity gradient in a straight and parallel flow. The representation of the characteristic of these fluids is a straight line with \( \mu \) as slope (Fig. 1). If the viscosity is not constant, depends from several process parameters including absolute temperature (T), shear rate \( \left(\dot{\gamma}=\partial \gamma /\partial t\right) \) and on the intrinsic properties of the n compounds that compose the fluid \( \left({M}_k={\sum}_1^n{M}_k\right) \). Therefore, the shear rate τ can be written as:

$$ \tau =\mu \left(T,\dot{\gamma},{M}_k\right)\frac{\partial \gamma }{\partial t} $$
(2)
Fig. 1
figure 1

Schematisation of shear stress-shear strain plot for Newtonian and Non-Newtonian fluids

Non-Newtonian fluids that report an increase of viscosity with the increase of the shear rate are called dilatant fluids (shear-thickening), while fluids whose viscosity decreases with an increase of shear rate are called pseudo-plastic (shear-thinning).

In Fig. 1, it is reported a schematisation of the different rheological behaviours.

Considering a non-Newtonian fluid suspension made of particles (subscript p) and solvent (subscript s), the constitutive equation of suspension viscosity (η) is modified in comparison of the previous definition (Eq. 2) as described by Stickel et al. [33] taking into consideration several parameters such as hydrodynamic radius (rp), density (ρp), and concentration (Cp) of particles, density (ρs), viscosity of solvent (ηs), and Boltzmann constant (kB):

$$ \eta =f\left({r}_p,{r}_s,{\rho}_p,{\rho}_s,{C}_p,{\eta}_s,{k}_BT,\dot{\gamma},t\right). $$
(3)

In order to reduce the number of variables, an alternative constitutive equation for these fluids is reported in literature [28] where dimensionless quantities are considered, including volume of particles in suspension (\( \zeta \)), Reynolds number function \( \left({\mathit{\operatorname{Re}}}_{\dot{\gamma}}\right) \) and Peclét number \( \left({Pe}_{\dot{\gamma}}\right) \), functions of shear rate:

$$ \eta =f\left(\zeta, {\mathit{\operatorname{Re}}}_{\dot{\gamma}},{Pe}_{\dot{\gamma}}\right) $$
(4)
$$ {\mathit{\operatorname{Re}}}_{\dot{\gamma}}=\frac{\rho_s{r}_p^2\dot{\gamma}}{\eta_s} $$
(5)
$$ {Pe}_{\dot{\gamma}}=\frac{\dot{\gamma}{r}_p^36\pi {\eta}_s}{k_BT} $$
(6)

Using the Peclét number, it is possible to correlate the rheological properties to the mechanical conditions on the fluid connecting the microscopic phenomenology to the macroscopic effects. When low shear rate is applied on the fluid \( \left({Pe}_{\dot{\gamma}}<1\right) \), the Brownian motion can reset the fluid structure in an equilibrium state that leads to a constant viscosity (Newtonian fluid). On the contrary, if a high strain rate is applied on the fluid, the deformation operates on a shorter time than Brownian motion and the solution is moved out of its equilibrium state, leading to a macroscopic viscosity increase (thickening). In this work, a suspension of silica particles (SiO2) in Polyethylene Glycol (PEG) is considered as secondary reinforcement for hybrid CFRP structures. The synthesis process of this STF and its compatibility with CFRP system is demonstrated in several works [19, 28, 34, 35] and it is not the focus of this work. On a microstructural level, when the shear rate is applied, an STF is characterised by a double structure formed by the arrangement of the silica particles and the solvent. Indeed, at low shear rate, the particles are able to organise themselves in a layered-flow structure due to the repulsive forces (Brownian forces) that ease their movement across the solvent molecules, leading to an initial reduction of viscosity. Increasing the shear rate above a critical threshold (shear-thickening onset) [20, 21], an hydrodynamic instability is progressively introduced into the system leading the particles to start losing their ordered state and interact via contact and/or clustering creating agglomerates that reduce the flow ability of the solution. Increasing the shear rate, a higher number of particles are able to interact, and these agglomerates start to create networks via contact forces that further decrease the mobility of the solution. Thus, since more energy is required for the fluid displacement, the viscosity of the entire system increases. Once the shear rate is removed, the jamming is removed and the stability of solution (initial viscosity) is restored.

3 Samples Manufacturing and LVI Experimental Setup

STF characterisation was carried out using rheological analysis obtaining the instantaneous viscosity-strain rate curve (reported in Fig. 2 used to calibrate the numerical model.

Fig. 2
figure 2

Viscosity- shear rate curve used in the numerical model to consider non-Newtonian contribute of STF

Hybrid STF-CFRP samples used during the experimental campaign were obtained following the procedure illustrated by Pinto et al. [28]. To manufacture the samples 8 plies of unidirectional carbon fibres prepreg were cut in rectangular shape 150 mm x 100 mm x ~ 0.25 mm and laminated following a quasi-isotropic lamination sequence [0/90/+45/–45]S. STF was fabricated using microsized particles of silica (SiO2) (25% wt.) dispersed in PolyEthylene Glycol (PEG). An ultrasonic bath was used during the mixing process in order to guarantee a good distribution of silica particles in the medium. The silica was added gradually to the PEG and manually mixed. The mixture was immersed in the ultrasonic bath after reaching the gel state. The final product was then stored under vacuum overnight to eliminate eventual air bubbles. This STF was then embedded and uniformly distributed into the CFRP structure as part of the stacking sequence using a roller to obtain a layer thinner than the carbon ply. This led to obtain a minimal weight increase in comparison with traditional laminates as also demonstrated in previous works [28]. Analysing the cost of the STF, silica particles and PEG can be easily obtained on the market for a moderate price that is negligible if compared to the cost of carbon fibre prepreg used to manufacture the hybrid and traditional structures. Final dimensions of impact samples were 150 mm x 100 mm with an average thickness of 2 mm. Traditional CFRP samples were used as reference. Images of the manufactured samples are reported in Fig. 3.

Fig. 3
figure 3

STF/CFRP plates impacted during the experimental campaign

In order to generate a BVID (fibre breakage, matrix failure and delamination) within the laminate structure and evaluate the impact response enhancement of the hybrid STF/CFRP structure, a drop-tower impact machine (Fig. 4) was used equipped with a weight of 12.867 kg and a hardened-steel hemispherical tip of 20 mm according to ISO 6603-2:2000 standard. All the plates were placed into the machine support and properly constrained to avoid vibrations and plates movements. Three different impact energies (6J, 20 and 40 J) were used during the experimental campaign with an impact velocity of 0.99 m/s, 1.80 m/s and 2.55 m/s respectively. In order to evaluate the damaged area within the impacted samples, an UltrasonicSciences Ltd. C-scan system employing a high-resolution 35 MHz was used.

Fig. 4
figure 4

Drop tower impact machine used during the experimental campaign

4 Numerical Implementation

In order to simulate the dynamic structural behaviour of hybrid CFRP material, an explicit 3D Finite Element (FE) model was developed using a Ls-Dyna code. The composite plate was modelled using 3D brick (solid) elements with a constant stress formulation (0.25 mm thick) in order to predict the through-the-thickness stress component that is crucial in composite failure evaluation. Rectangular (150 mm x100mm x 2 mm) CFRP plates were simulated using a quasi-isotropic stacking sequence [0/90/+45/–45]S. FE mesh size was fixed with an element size of 0.5 mm x 0.5 mm after a trade-off procedure between results accuracy and computational time. The impactor was modelled as a hemispherical elastic solid body of 12.864 kg and 20 mm in diameter. A penalty-based contact was defined to guarantee the interaction between the CFRP and the impacting body. Impact velocity used for each case of the experimental campaign was defined for all the nodes of the impacting body. Figure 5 illustrates the FE mesh of the numerical model:

Fig. 5
figure 5

FE mesh of CFRP reporting composite plate impact and impact direction. In detail, the STF layer and plies layup is reported

The STF layer was located at three different positions (as reported in Fig. 6) in order to evaluate its effect on the impact resistance of different composite laminates. STF_T (Fig. 6a) is the configuration with the STF placed close to the impact surface, STF_B (Fig. 6b) is the one with the STF far from the impact area and STF_C (Fig. 6c) has the STF layer on the midplane of the laminate.

Fig. 6
figure 6

Location of STF layer position along CFRP thickness: a)STF-T: STF placed close to the impact surface; STF-B:STF placed far from impact surface; and STF-C: STF placed on the midplane of the laminate

In order to simulate the boundary conditions used during the experimental setup, the bottom and top edges of the plate were constrained in all the translational degrees of freedom. All the nodes relative to the impact tip were constrained along the x and y direction, ensuring motion only along through-the-thickness direction (z-axis).

5 Arbitrary Lagrangian Element (ALE)

In order to model the STF in Ls-Dyna, the Arbitrary Lagrangian Element (ALE) formulation was used [36]. This element type is used in applications where Fluid-structure Interaction (FSI) must be taken into account. In this formulation, the mesh elements are not fixed in space like for the Eulerian element formulation or attached to the material as for Lagrangian element formulation, but the single domain is able to adapt its shape arbitrarily to optimise the elements shape, leading to a more accurate and computationally efficient multi-material simulation. The use of ALE element formulation is well-known in fluid-dynamics [37], non-linear structural [38] and fluid-structural coupled [39] analyses in which large mesh deformation is considered. The integral conservation equations are obtained from the well-known Eulerian forms in which the material velocity v (bold for vector) is replaced with the convective velocity c which represents the relative velocity between the material and the mesh (c=v-v). The spatial coordinates \( \chi \) are expressed using a particle reference system which is mapped using a combination of Lagrangian and Eulerian reference systems [40]:

$$ \mathrm{Mass}\kern4.25em {\left.\frac{d\rho}{d t}=\frac{\partial \rho }{\partial t}\right|}_{\chi }+\boldsymbol{c}\cdotp \mathbf{\nabla}\rho =-\rho \mathbf{\nabla}\cdotp \mathbf{v} $$
(7)
$$ \mathrm{Momentum}\kern4.25em \rho \left({\left.\frac{\partial \boldsymbol{v}}{\partial t}\right|}_{\chi }+\left(\boldsymbol{c}\cdotp \mathbf{\nabla}\right)\boldsymbol{v}\right)=\mathbf{\nabla}\cdotp \boldsymbol{\sigma} +\rho \boldsymbol{b} $$
(8)
$$ \mathrm{Total}\ \mathrm{Energy}\kern4.25em \rho \left({\left.\frac{\partial E}{\partial t}\right|}_{\chi }+\boldsymbol{c}\cdotp \mathbf{\nabla}\mathbf{E}\right)=\mathbf{\nabla}\cdotp \left(\boldsymbol{\sigma} \cdotp \boldsymbol{v}\right)+\boldsymbol{v}\cdotp \rho \boldsymbol{b} $$
(9)
$$ \mathrm{Internal}\ \mathrm{Energy}\kern4.25em \rho \left({\left.\frac{\partial e}{\partial t}\right|}_{\chi }+\boldsymbol{c}\cdotp \mathbf{\nabla}\mathrm{e}\right)=\boldsymbol{\sigma} :{\mathbf{\nabla}}^{\boldsymbol{S}}\boldsymbol{v}. $$
(10)

Where \( \rho \) is the density (kg/m3), v is the material velocity vector (m/s), σ is the Cauchy stress tensor (p.a.), b is specific body force vector (N/m3), E and e are the specific total and internal energy (J/m3) where only mechanical contributes are included, and \( {\mathbf{\nabla}}^S\mathbf{v} \) indicates the strain rate tensor which is the symmetrical part of the velocity gradient \( \boldsymbol{\nabla}\boldsymbol{v}\left({\boldsymbol{\nabla}}^{\mathrm{S}}\boldsymbol{v}=\frac{\mathbf{1}}{\mathbf{2}}\boldsymbol{\nabla}\boldsymbol{v}+{\boldsymbol{\nabla}}^T\boldsymbol{v}\right) \). The arbitrary component of the constitutive equations is reflected by the left hand side of the equations, while the Eulerian form is on the right-hand side. The Ls-Dyna (R10) used to perform the numerical analysis solves these equations for each time step and for each node, in order to accurately simulate the mechanical behaviour for the STF material. Another important parameter to define when ALE element formulation is involved, is the boundary conditions applied on these elements. Two main conditions for material surface need to be satisfied: (i) no particles can cross the material boundary and (ii) the stress has to be continuous across the surface. Considering these conditions, a fluid-structure interaction is considered for this work to model the ALE element boundary conditions. This interaction is formulated considering the particles velocity coupled with the nodes of rigid or deformable structures. Thus, the condition (ii) has to be coupled with Eq. 7, 8, 9 and 10 in order to avoid that the viscous fluid overlaps or detaches the structure during the interaction. These conditions are called continuity of displacement and continuity of velocity and are expressed with the equation below (Eqs. 12 and 13 respectively):

$$ \boldsymbol{n}\cdotp \boldsymbol{w}=0\ or\kern0.5em \boldsymbol{n}\cdotp \boldsymbol{v}=\boldsymbol{n}\cdotp \boldsymbol{v}^{\prime } $$
(11)
$$ \boldsymbol{u}={\boldsymbol{u}}_{\boldsymbol{s}} $$
(12)
$$ \boldsymbol{v}={\boldsymbol{v}}_{\boldsymbol{s}} $$
(13)

where u represents the displacement. S subscript indicates the structure contribute while no subscript indicates the fluid one. Following this approach, fluid mesh movements are independent from the fluid itself due to fluid ALE formulation, and two nodes (fluid and structure) are placed in the same interface point, where they are constrained to be contiguous and permanently aligned during sliding. Dynamic conditions at the interface are verified along the structure boundaries but velocity and stress compatibility are required since the stress in the fluid must be equal to the stress in the structure. Based on this approach, it is possible to develop a STF/CFRP FE model that integrates an FSI method using an ALE formulation and predicts the mechanical response of the entire CFRP/STF system under dynamic loading.

ALE elements were simulated using an element formulation 11 (1 point ALE multi-material element) and a user-defined (UD) material card ALE_VISCOUS. This user-defined card was developed considering the strain rate \( {\dot{\varepsilon}}_{jj} \) applied on the element. Both hydrostatic εjj and deviatoric εji strains are considered [41]:

$$ {\dot{\varepsilon}}_{jj}=\frac{d{\varepsilon}_{jj}}{dt} $$
(14)
$$ {\dot{\varepsilon}}_{ji}=\frac{1}{2}\frac{d{\varepsilon}_{ji}}{dt} $$
(15)

The viscosity tensor terms \( {\eta}_{ji} \) were computed using a Carreau model where the dependency from deviatoric strain rate is considered in Eq. 16 for each step of Ls-Dyna calculations:

$$ {\eta}_{ji}\left({\dot{\varepsilon}}_{ji}\right)={\eta}_0+\left({\eta}_{\infty }-{\eta}_0\right){\left(1+{\left(K{\dot{\varepsilon}}_{ji}\right)}^2\right)}^{\frac{n-1}{2}} $$
(16)

Where n is the power law index that is > 1 for shear-thickening fluids, K is characteristic time, \( {\eta}_0 \) is the value of viscosity at zero shear rate and \( {\eta}_{\infty } \) is the value for the viscosity at maximum shear rate. These parameters are evaluated and implemented using an experimental viscosity/deviatoric strain rate curve (Fig. 2). The shear stress tensor can be computed from the calculated viscosity tensor using the formula:

$$ {\boldsymbol{\tau}}_{\boldsymbol{ji}}=2{\eta}_{ji}{\dot{\boldsymbol{\varepsilon}}}_{ji} $$
(17)

The total stress tensor is calculated as:

$$ {\boldsymbol{\sigma}}_{\boldsymbol{ji}}={\boldsymbol{\tau}}_{\boldsymbol{ji}}-\left(-\frac{\sigma_{jj}}{3}-{E}_B{\varepsilon}_{jj}\right){\delta}_{ji} $$
(18)

Where EB is the bulk modulus (GPa) and \( {\delta}_{ji} \) is the Kronecker delta.

In ALE material card, material density was set to 1583(kg/m3) while cut-off pressure was set to -1·1030(Pa) since no cavitation is assumed in the real impact event. An instantaneous viscosity-shear rate curve (Fig. 2) is, then, assigned to consider the Non-Newtonian behaviour during the impact event. A Linear Polynomial Equation of State (EOS) was used to take into consideration the fluid nature and a bulk modulus (3 GPa) was defined to evaluate the density variation during fluid movements. In order to constrain the STF layer within the specified CFRP section portion, an elastic shell box (Lagrangian elements) was created surrounding the fluid mesh to avoid fluid leaking. The box thickness was considered 100 times thinner than a single layer to avoid eventual contributes in the stress calculation. Density and elastic modulus of the box were set up similar to STF material. Boundary conditions were defined using the dedicated boundary conditions card (LAGRANGIAN_IN_SOLID) to define the interface between structural and fluid material. All default settings were unchanged except for the parameter CTYPE and DIRE where penalty coupling and compression transfer only were chosen since the fluid is unable to bear tension.

6 Composite Material Model

In order to guarantee an optimal and appropriate prediction of the structure failure, a material model LAMINATED_FRACTURE_DAIMLER_PINHO was used to model the CFRP material using MAT_261 in Ls-Dyna library. This material is characterised by an orthotropic elastic behaviour until failure and a continuum fracture-based damage model after the damage initiation, which is able to predict with a certain accuracy fibre breakage and matrix failure modes. In particular, compressive and tensile damage criteria for fibre breakage are considered to analyse the fibre contribute in the structural failure while intra-laminar tensile and shear matrix failure criteria are used to analyse the matrix behaviour during failure. Once a typology of damage is initiated on a specific element of the ply, an accumulation of damage is considered until the fracture conditions are satisfied. After this, the element needs to be deleted (eroded) as they are unable to bear any load.

The element erosion procedure was carried out when one of the multiple fracture conditions was met. Table 1 shows the mechanical properties used in the model:

Table 1 -MAT_261 orthotropic material card paramters: density (RO), modulus of elasticity (E11, E22, E33), poisson’s ratios (PR12 ,PR31, PR32), shear modulus (G12, G23, G31), normal and transverse strenght under traction and compression (X1t, X1c, X2t, X2c), shear strenght (S12, S23, S31), compresisve fibre failure energy (ENKINK), tensile fibre failure energy (ENA), Intralaminar matrix tensile energy failure (ENB), Intralaminar matrix trasnverse shear energy failure (ENT) and Intralaminar matrix longitudinal shear energy failure (ENL)

Toughness values for fibre fracture are reported in [42]. Damage propagation between laminae (delamination) was implemented using a TIEBREAK_CONTACT [43, 44], defined between two adjacent plies. The two surfaces nodes were considered tied until the interface failure criteria are satisfied. Afterwards, the nodes were untied and a typical penalty-contact was defined on both the surfaces. Guidelines on the typical values for interlaminar failure properties to define this contact were obtained from literature [45].

7 Results and Discussion

In order to allow the optimisation of the STF layer along the CFRP thickness, a validation of the numerical model is required. The results for reference and STF_T samples obtained during the experimental campaign are reported in Fig. 7 and compared with the relative numerical results to demonstrate the ability of the simulations to predict the impact response of the STF.

Fig. 7
figure 7

Experimental and numerical Force-Displacement curves at 20 J used for numerical model calibration

As shown in the graphs (Fig. 8a), a good correlation between experimental and numerical results is identified allowing a good calibration of the numerical model and a good reliability for its further optimisation. Similar force peaks, maximum displacements and absorbed energies are evaluated at 20J with an error between simulated and experimental data of 0.7%, 1% and 5% respectively for reference curves, and 8%, 1% and 4% the STF curves for STF (STF_T) samples. C-scan images are then reported in Fig. 8 in which damage extension was evaluated for traditional CFRP (Fig. 8b) and STF/CFRP (Fig. 8c) at 20 J and compared to the numerical one (Fig. 8d and Fig. 8e respectively) using the eroded element area.

Fig. 8
figure 8

Comparison between experimental and simulated damage: (a) experimental reference at 20J, (b) FEM reference at 20J, (c) experimental STF/CFRP at 20J, (d) simulated STF/CFRP at 20J and (e) bar plot of damage extension

Analysing the damaged area extension, it is possible to observe that the numerical model can simulate the damage generation mechanisms of STF/CFRP and traditional CFRP structures with an error between simulated and experimental results of 9% and 12% respectively.

Considering these results, it is then possible to confirm that the numerical model is able to predict the impact behaviour of the STF/CFRP system allowing the optimisation of the STF position to maximise its effect on the impact properties.

Numerical results of the optimisation analysis of impact properties are reported in Fig. 9 where the output data of different positions for the STF layer are shown at 6J, 20 and 40 J. In Fig. 10 output impact data are shown illustrating values of maximum contact force, maximum displacement of the plate, absorbed energy, and internal energy at force peak.

Fig. 9
figure 9

Experimental and numerical Force-Displacement curves at 6J (no experimental curves) (a), 20 J (b) and 40J (c)

Fig. 10
figure 10

Impact results from experimental and numerical tests: (a) Maximum contact force, (b) Maximum displacement, (c) absorbed energy, (d) internal energy and (e) Damaged area extension

As it is possible to see from the curves, the introduction of the STF layer modifies the impact response of the CFRP in function of the impact energy applied on the sample.

Analysing the curves at 6J and the relative impact results, it is possible to observe that the thickening effect of the embedded fluid leads to a variation of maximum displacement (+ 10%, + 11% and + 15%) while it does not seem to affect the maximum force (+ 6%, -2% and − 4%) for all the analysed STF positions (STF_T, STF_C and STF_B respectively) with a small stiffness reduction given by the slope of the Force-Displacement curve. The reason of this behaviour is that the impact response of the system is altered by the presence of STF that offers a reinforcing due to its thickening ability. In other words, the higher is the strain rate applied during the impact event, higher is the viscosity of the fluid (thickening) and, consequently the resistance that the fluid is able to oppose against motion [34]. During the impact event at 6J, the applied strain rate values are moderate and they have only a partial reinforcing effect on the material (STF_T case) therefore this is unable to compensate the structural detriment given by the presence of the additional interface (STF), resulting in lower mechanical properties than traditional CFRP [46]. Considering the 20J impact results, it is possible to observe similar values in both contact force and maximum displacement (+ 5% and − 3% for STF_T, + 4% and + 1% for STF_C and + 5% and − 1% for STF_B respectively) with no appreciable stiffness variations in comparison with the reference. It is important to highlight that all the STF-reinforced samples show no signs of severe structural failures as opposed to the traditional CFRP that displays a severe load drop in the force-displacement plot (see Fig. 9). This is due to the higher strain rate applied on the sample that improves the thickening effect of the STF, reinforcing the structure and mitigating damage. STF reinforcing effect is clearly shown in the curves relative to 40J impact (Fig. 10) where an improvement of impact properties is reported (maximum contact force) for all the STF configurations with + 36% for STF_T configuration, + 25% for STF_C and + 22% for STF_B in comparison with the reference sample. In particular, the latter shows complete penetration and a residual energy of 7J is evaluated on the projectile.

In order to investigate more in depth, the STF effect on the mechanical response of the structure and understand the role of STF in the damage suppression phenomenon, absorbed and internal energy are evaluated and reported in Fig. 10c and d respectively. Analysing how the impact energy is transferred to the laminate during the impact event, it is possible to consider different contributes as reported by the Eq. 19 where Einternal represents the energy stored into the material as internal energy while Edissipation is the energy dissipated by the system:

$$ {E}_{impact}={E}_{internal}+{E}_{dissipation} $$
(19)

Using Ls-Dyna post-processing tools, the internal energy accumulated within the global system at a certain time step can be estimated considering the internal energy as a sum over six directions of the product of stress, strain and volume of the single elements. The internal energy contribute includes both the elastic energy of the material and STF energy used during particles interaction:

$$ {E}_{internal}\left({t}_n\right)={E}_{elastic}+{E}_{STF}={E}_{internal}\left({t}_{n-1}\right)+\sum \limits_i\left({\sigma}_i{\varepsilon}_iV\right). $$
(20)

Where \( {E}_{elastic}\left({t}_{n-1}\right) \) is the internal energy contribute for the previous time step, the index i represents the i-th component of stress and strain and V is the volume of the element.

The dissipated energy (Eq. 21) is defined as the sum of the energy contribute Edamage required to generate damage in the structure and the energy Enon−linear that considers the energy dissipated in non-linear phenomena including damping, friction and vibration phenomena. All the energy dissipated by the system is considered in the absorbed energy output data. No numerical energy loss are considered during this analysis

$$ {E}_{dissipation}={E}_{damage}+{E}_{non- linear}={E}_{absorbed} $$
(21)

When the STF is introduced into the CFRP structure, the presence of the additional reinforcement changes the impact response of the laminate and the distribution of impact energy within the STF/CFRP structure. The thickening ability of STF increases the internal energy threshold of the system since more energy is required to start the fluid’s flow, allowing the system to tolerate a higher amount of stored energy before starting any structural failure (i.e. cracks or delamination). Hence, since less energy is used by the system to dissipate impact energy via damage generation, the extent of the damaged area is reduced. Moreover, due to the dependency of STF viscosity from the applied strain rate (Fig. 2), it is possible to evaluate an enhancement of STF ability in storing impact energy when greater solicitations (higher impact energy) are applied on the structure.

Based on these considerations and analysing Fig. 11a, 6J impact results for STF/CFRP samples show a consistent variation in absorbed energy (-31% for STF_T, -21% for STF_C and − 5% for STF_B) and internal energy for STF_T, STF_C and STF_B configurations (+ 16%, + 6% and − 4% respectively) in comparison with traditional CFRP. This is due to the aforementioned STF ability to transfer a part of energy from the dissipated energy contribution to the internal energy one allowing an improvement in the impact response. Indeed, this behaviour is confirmed by analysing the damaged area extension of the three configurations (STF_T, STF_C and STF_B) where no signs of damage are reported. Analysing the 20 J case (Fig. 11b), it is possible to observe a lower absorbed energy in comparison with traditional CFRP for all STF configurations (STF_T:-42%, STF_C:-33%, STF_B:-21%) and a higher internal energy (STF_T: +76%, STF_C: +53% and STF_B: +28%). This is due to the improved STF ability at higher impact energy levels in storing more energy as internal one with a consequent reduction in energy used to generate damage. This is confirmed evaluating the variation of damaged area between traditional CFRP and STF/CFRP structures (STF_T, STF_C and STF_B) that showed a reduction of 48%, 21% and 23%. This enhancement is further illustrated by the 40J (Fig. 11c) impact results where STF effect is evident in protecting the STF/CFRP structure from the impact penetration. Indeed, all the impact energy that can cause critical damage to the structure is stored into STF layer as internal energy via particles interaction, minimising the dissipated energy contribute and the damage generation with a reduction of absorbed energy of -27% for STF_T, -20% for STF_C and − 13% for STF_B in comparison with traditional CFRP where the total perforation is identified. A similar trend is also shown by evaluation of the damaged area for which reductions of 44%, 33%, 31% between STF_T, STF_C and STF_B respectively and traditional laminates are reported. On the other hand, the internal energy is strongly affected by the STF contribute with an increment of + 161% for STF_T, + 136% for STF_C and + 106% for STF_B.

Fig. 11
figure 11

Percentage variation in comparison with traditional CFRP (reference) of absorbed energy, internal energy and damaged area at 6J,20 and 40 J for the three STF/CFRP configurations: (a) STF_T, (b) STF_C and (c) STF_B

Based on these results, it is clear that STF_T configuration (STF closer to the impact surface) is the one that shows the highest impact resistance improvement among the tested configurations. A possible explanation for the higher impact performance of STF_T is found analysing STF position and its relative distance from the area where the impact occurs since the strain rate distribution is more intense in the surroundings of the impact location. STF_C configuration is intermediate in terms of impact performances since the presence of STF along the middle plane of the structures intensifies shear stress localisation due to the mechanical properties difference between STF and CFRP. Hence, STF has only a partially enhancing effect on the impact properties. The lowest impact performances improvement is shown by the STF_B configuration due to the presence of TPU in the farthest area from the impact location that generates a consistent stress intensification effect in the tensile portion of the laminate due to the additional interfaces inserted within the structure. However, an impact improvement is still observed when a high strain rate is applied on the structure (20 and 40 J).

8 Conclusion

This work proposes and optimises the use an STF layer within a CFRP laminate to improve low-velocity impact performances. A Ls-Dyna code was developed and implemented using an ALE-based fluid-structure interaction analysis to simulate the STF behaviour and its interaction with the CFRP. STF/CFRP hybrid samples and traditional CFRP plates (used as reference) with different levels of impact energy (20 and 40 J) were experimentally investigated and results used to validate the numerical model in order to be effective in predicting impact behaviour of STF/CFRP material. The optimisation of STF layer along the CFRP thickness was numerically carried out comparing the impact response of STF/CFRP plates with three different STF configurations along the laminate thickness (close to the impact surface, on the middle plane and far from impact surface). Results obtained for high energy impacts showed a considerable increase in contact force (~ 36%), internal energy, and a reduction in absorbed energy (~ 42%) and damaged area (~ 35%) due to the presence of STF that was able to store an additional fraction of impact energy as internal energy via the interaction of the particles. Perforation shown by traditional CFRP samples impacted at 40 J was totally prevented by the insertion of STF at any position within the structure, but its effect was maximised when STF was embedded as close as possible to the impact surface. Based on these considerations, hence, STF/CFRP structures can be successfully simulated and optimised for their applications in several advanced sectors including aerospace, automotive and energy (wind blades) representing an important step-up in the production of high-resistant hybrid composite structures.