Nomenclature

FE

finite element

SPH

smooth particle hydrodynamic

CFC

SAE channel filter class

EOS

equation of state

ACC

accelerometer

P

pressure transducer

${E_i}$

linear in internal energy

${C_i}$

material constant

r

position

m

mass

p

pressure

v

velocity

ρ

density

$\mu $

volumetric strain

2.0 Numerical model

SPH is a mesh-free method originally proposed by Lucy, Gingold, and Monaghan in 1977 [Reference Lucy29, Reference Gingold and Monaghan30]. It accurately and stably solves a set of partial differential equations without the need for any grid connecting the particles. The SPH method employs a discrete set of smooth particles with position, velocity and mass attributes to approximate the values and derivatives of continuous variables. Subsequently, it calculates the weighted averages of all neighbouring particles. A simplified overview of the SPH method is depicted in Fig. 1. In SPH, fluid is treated as a collection of particles, each possessing associated physical properties such as position r(t), mass m(t), density ρ(t), pressure p(t), and velocity v(t). The computational domain is then discretised into spatial bins of size 2 h, and particles within the particle’s own space and adjacent spaces are identified. For a given particle, only interactions with neighbouring particles within a distance less than 2 h are considered. The interactions are computed for the given particle and the selected neighbouring particles by solving the following NS equations:

(1) \begin{align}\frac{{d{v_i}}}{{dt}} = - \sum\limits_j {{m_j}} \left( {\frac{{{p_i}}}{{\rho _i^2}} + \frac{{{p_j}}}{{\rho _j^2}} + {\prod\nolimits_{ij}}} \right){\nabla _i}{W_{ij}} + {F_i},\end{align}

(2) \begin{align}\frac{{d{\rho _i}}}{{dt}} = \sum\limits_j {{m_j}} \left( {{v_i} - {v_j}} \right){\nabla _i}{W_{ij}}.\end{align}

Figure 1.

Brief overview of SPH method.

The position and velocity of the particle for the next time step are then updated based on the results of the following calculation. Figure 2 illustrates the detailed calculation cycle.

(3) \begin{align}{r_i}\!\left( {t + \Delta t} \right) = {r_i}\!\left( t \right) + \Delta t{v_i}\!\left( t \right),\end{align}

(4) \begin{align}{v_i}\!\left( {t + \Delta t} \right) = {v_i}\!\left( t \right) + \Delta t{a_i}\!\left( t \right).\end{align}

Figure 2.

The calculation cycle of SPH method.

It is widely recognised that structural loading for water impact is significantly different to hard surfaces [Reference Authority31]. As shown in Fig. 3, the aircraft’s kinetic energy is absorbed by the plastic deformation of the metal structure during impact on a hard surface, reducing the impact load and acceleration to survivable levels. However, in the case of water impact, the landing gear and conventional flat under floor structures perform poorly in absorbing impact energy [Reference Vignjevic and Meo32, Reference Ren and Xiang33]. The impact load is distributed along the skin, and energy absorption depends primarily on bending behaviour of the skin as well as the limited plastic collapse of supporting frames [Reference Taber34, Reference Ren and Xiang35]. In severe cases, aircraft skin rupture or inability to transmit impact loads could significantly reduce the efficiency of sub-floor structural energy absorption. Internal structures may also suffer secondary damage due to water ingress [Reference Thuis and Wiggenraad36]. Therefore, it is imperative to investigate the effects of water impact on different structures to enhance the survivability of passengers during aircraft emergency water landing.

Figure 3.

Comparison between impact on a rigid floor and on water.

Numerous accident cases have demonstrated that the aircraft skin rupture during ditching is the main reason affecting passenger survival rates [Reference Brooks, MacDonald, Donati and Taber37, Reference De Florio38]. To investigate the aircraft under floors impacting water surfaces and analyse the load on aircraft skin during water impact, the Italian Center of Aerospace Research (CIRA) conducted a water entry impact experiment at the Laboratory for Impact Tests on Aerospace Structures [Reference Borrelli, Ignarra and Mercurio39]. The test specimen of steel structure is dropped into the water at different velocities to generate an experimental database for validating numerical models. The test specimen is a semi-cylindrical aircraft skin reinforced with three frames positioned symmetrically on both sides and at the centre of the skin. To secure the structure on the drop tower trolley, L-shaped and U-shaped steels are mounted above the structure. Both the skin and frames are constructed using 2 mm thick hot-rolled steel (DD11). The material is defined by bilinear elastic plastic law and takes into account the strain rate effect by adding the stress-plastic strain curve. As the experimental impact velocities are not sufficient to cause damage to the material, the material model does not need to consider material failure. The structure had a length of 1,200 mm, a cross-sectional radius of 475 mm, a central circular hole with a diameter of 150 mm in the frames, a specimen weight of 66.4 kg and the trolley weight of 174 kg.

The structural model employed shell elements with a total of 7,592 elements. All shell elements are implemented in the Belytschko-Tsay format, and the selected thicknesses are consistent with the experimental setup. As there is no structural failure observed in the experiments, interconnection between models is established through merged nodes. The material model for the drop tower trolley utilised the rigid model, with the addition of mass points at the centre to simulate the actual trolley weight. To align with experimental constraints, only the vertical degrees of freedom for the trolley are retained.

The water model is structured using the SPH algorithm. The water domain consists of a cubic shape measuring 2,000 mm × 1,600 mm × 800 mm, composed of 40,000 particles. The Particle Approximation Theory employs an enhanced fluid formulation. The equation of state (EOS) is linear in internal energy ( ${E_i}$ ) and it defines the pressure-volume relationship as:

(5) \begin{align}\rho = {C_0} + {C_1}\mu + {C_2}{\mu ^2} + {C_3}{\mu ^3} + \left( {{C_4} + {C_5}{\mu ^2} + {C_6}{\mu ^3}} \right)\!{E_i},\end{align}

where in this case $\mu $ is the volumetric strain defined as:

(6) \begin{align}\mu = \frac{\rho }{{{\rho _0}}} - 1,\end{align}

$\rho $ is the current density, ${\rho _0}$ is the reference water density and ${C_i}$ (i = 1, 2, …, 6) are material constants. The following values are used: ${C_1}$ = 2.723 GPa, ${C_2}$ = 7.727 GPa, ${C_3}$ = 14.66 GPa, ${C_0} = {C_4} = {C_5} = {C_6} = 0$ [Reference Toso40].

The penalty methods are implemented with the contact interaction between the structural model and the water model, and the soft constraint mode is activated. The soft constraint option’s force scaling factor is set to 0.5. In the soft constraint option, the interface stiffness is based on the nodal mass and the global time step size. This method of calculating interface stiffness will typically give a much higher stiffness value than would be obtained using the bulk modulus. Therefore, this method is the preferred approach when water interacts with metals. Moreover, self-contact is enabled for the skin to prevent penetration. The structural model is initialised with impact velocities identical to those in the experiments, and gravitational effects are taken into account. Acceleration data are collected at the top plate of the model, representing the average acceleration at various points on the plate. Pressure values obtained by dividing the collected force by the corresponding area. The experiments involve drop tests with velocities of 3 m/s, 8 m/s and 10 m/s, and the overall simulation model are depicted in Fig. 4.

3.0 Numerical simulation results

3.1 Skin deformation

To analyse the impact behaviours of the aircraft skin structure by simulation, the experimental and numerical results obtained from Borrelli and Grimaldi are used to validate the devised numerical model [Reference Borrelli, Mercurio and Alguadich41, Reference Grimaldi, Benson, Marulo and Guida42]. The measured deformation sections are illustrated in Fig. 5. In Fig. 5(a), for an entry velocity of 8 m/s, the maximum relative displacement at Skin Section I is 31 mm. In Fig. 5(b), with an entry velocity of 10 m/s, the maximum relative displacement at Skin Section II is 45 mm. As shown in Fig. 5, there is a strong correlation between the experimental results and the present. Both values exhibit an error within 10% when compared to the experimental measurements obtained from the specimen post-impact. In the present results, the side of the skin is deformed outwards by the bending moment generated by the bottom deformation. This trend of the predicted displacements is closer to the experimental results and is not found in previous results. The symmetry of the deformation on both sides of the skin at different velocities aligns closely with the experimental findings, providing further evidence of the accuracy of the simulations.

Figure 4.

Steel structure finite element model and water particle model.

Figure 5.

Shape of the skin and measurement of deformation of section [Reference Grimaldi, Benson, Marulo and Guida42].

3.2 Acceleration and pressure

The positions for the accelerometer (ACC1, ACC2, ACC3 and ACC4) and pressure transducer (P1, P2, P4, P6 and P7) are illustrated in Fig. 6. Due to symmetry considerations, three locations on the skin are designated as pressure sensor locations in the simulation.

Figure 6.

Transducer deployment location [Reference Borrelli, Mercurio and Alguadich41].

The acceleration-time histories and pressure-time histories at different velocities are depicted in Figs. 7, 8 and 9. Sim1 represents the simulation results from Borrelli, Sim2 from Grimaldi and Sim1 employs filtering with SAE Channel Filter Class (CFC) 60, and Sim2 uses CFC180 for filtering. The water model of Sim1 utilises a mixed SPH-FE approach and Sim2 used the SPH method to model. In present simulation, both experimental and simulated data are subjected to filtering with the SAE Channel Filter Class (CFC) 60 Hz [Reference Cobb, Tyson and Rowson43]. In order to facilitate comparison, all the numerical curves are adjusted on the time scale to fit the experimental curves.

Figure 7.

At 3 m/s, time histories of acceleration and pressures [Reference Borrelli, Mercurio and Alguadich41, Reference Grimaldi, Benson, Marulo and Guida42].

Figure 8.

At 8 m/s, time histories of acceleration and pressure [Reference Borrelli, Mercurio and Alguadich41, Reference Grimaldi, Benson, Marulo and Guida42].

Figure 9.

At 10 m/s, time histories of acceleration and pressure [Reference Borrelli, Mercurio and Alguadich41, Reference Grimaldi, Benson, Marulo and Guida42].

As depicted in Figs. 8(a) and 9(a), Sim1 and Sim2 demonstrate a strong correlation in peak accelerations at impact velocities of 8 m/s and 10 m/s. However, as depicted in Fig. 7(a), at an impact velocity of 3 m/s, the water models of Sim1 and Sim2 exhibit excessive rigidity resulting in an overestimation of the acceleration peak. The results show that the water model fails due to excessive rigidity when the impact velocity is insufficient to induce plastic deformation in the specimen. With regards to pressure, accurate predictions are observed in soft areas (bays), while overestimations are evident in stiff areas (under the frames).

Building upon the findings of previous research, this study further refined the model parameters. The entire water pool model is represented using a particle-based approach. The particle spacing of the present simulation is chosen taking into account calculation accuracy and computational cost. To enhance stability in contact calculations, merged node connections are established between various models. The soft constraint formulation is implemented, and the scale factor for constraint forces of the soft constraint option is set to 0.5. The correlation between numerical and experimental peak accelerations at different velocities is significantly improved, with errors kept within 5%. Pressure values is the most challenging to measure. Beneath the stiff areas, the optimisation of the parameters demonstrated noticeable improvement, with errors controlled within 10%. In the soft areas, the predicted values aligned closely with experimental results when no plastic deformation occurred in the material. While there are still some errors in the predicted values at the lateral positions beneath the stiff areas, notable progress is observed compared to previous results. All errors of simulation are within 20%, with six peak data errors below 5%, showcasing strong correlations. Particularly, peak accelerations exhibit excellent agreement with experimental values across different impact velocities.

4.0 Simulation of different sub-floor shapes

Three different shapes of helicopter sub-floor structures are established, as shown in Fig. 10. The sub-floor shapes include the flat, the cylindrical surface with transverse curvature and the ellipsoid with double curvature. A coupled FE-SPH method is employed, which is not only capable of predicting non-linear structural collapse, but also represents the underlying kinematics of water and pressure transmissibility. As the weight of the floor structure is primarily concentrated on the trolley, weight variations among the models are less than 1% of the overall weight. Therefore, the impact of weight variation on the results can be disregarded. The impact velocities are 3 m/s, 8 m/s and 10 m/s, respectively.

Figure 10.

Under floor models.

Acceleration data are collected from the top plate, and the values represent the average acceleration at various points on the plate. As depicted in Fig. 11, the peak acceleration upon water impact exhibits a gradient distribution among the flat, cylinder and ellipsoid shapes. The flat surface experiences the highest impact acceleration, while the ellipsoid surface encounters the lowest. Experimental results demonstrate that appropriately shaped surfaces can effectively mitigate the impact during structural water entry. With increasing velocity, the deceleration effect persists, although its efficacy diminishes gradually. In short, at different speeds, the ellipsoid is the best of the three shapes.

Figure 11.

Time history of acceleration.

The P4 pressure sensor is positioned at the central location of the under floor. And it is the initial contact point with water after the descent of all structures. The pressure-time history of P4 is illustrated in Fig. 12. As the structure enters the water at 3 m/s, the peak pressures at different locations on various under floor structures are nearly identical. With an increase in velocity, the pressure levels on the flat and cylinder surfaces remain consistent, while the advantage of the ellipsoid shape becomes apparent, displaying significantly lower pressures on its skin surface compared to the other two structures.

Figure 12.

Time history of pressure P4.

The pressure-time curve for P1 is depicted in Fig. 13. As the structure enters the water at 3 m/s, the pressure at P1 on the flat surface is significantly higher than on the curved surface. With an increase in velocity, the difference between the two gradually diminishes. The pressures at P1 for the cylinder and ellipsoid surfaces are nearly identical, and at a velocity of 10 m/s, the slightly smaller peak pressure on the ellipsoidal surface compared to the cylinder.

Figure 13.

Time history of pressure P1.

From Fig. 14, it can be observed that, due to the curvature of the ellipsoid surface deflecting the water flow radially, the surface deformation changes relatively little, thereby minimising the loads at these locations. At 3 m/s, the pressure at P2 on the ellipsoid surface is comparable to that on the flat surface, as the curvature of the ellipsoidal surface is relatively small. With an increase in velocity, the superior structural resistance of the cylinder surface over the flat surface gradually becomes apparent. In general, the ellipsoidal sub-floor reduces pressures at the side observably due to the bi-directional curvature.

Figure 14.

Time history of pressure P2.

Figure 15 illustrates the plastic strain nephogram at various time for different models when the velocity is 8 m/s. The distribution of plastic strain on the ellipsoid surface is more uniform, with severe plastic deformations not extensively concentrated. The ellipsoid under floor has a significant inward deflection on either side of the frame. It is able to better transfer the load from the skin to the frame, so the main energy absorbing components are utilised in an efficient manner.

Figure 15.

Time histories of effective plastic strain at 8 m/s.

The deformation of the skin during ditching is investigated. Considering that the mechanical properties of the aircraft floor frame are much better than the skin, it is worthwhile to set the frame as a rigid body to compare the changes in skin deflection. In addition, a longitudinal stiffener was installed at the longitudinal centreline of the new set of models, aiming to analyse the influence of longitudinal stiffeners on skin deflection during water impact. Figure 16 shows the vertical displacement contour map when the skin undergoes maximum deformation. When there is no longitudinal stiffener, the maximum deformation positions of the skin are all on the centreline. On this line, the midpoint is supported by the frame and no deformation occurs. The location of maximum deformation of the skin is not right in the middle of the two frames, but more towards the edge frames. The area of deformation of the flat floor is the largest, and the deflection value is also the largest, at 21.2 mm. The maximum deflection of the cylindrical floor is 18.3 mm. The deformation area of the ellipsoidal floor is the smallest, and the maximum deflection value is also the smallest, at 16.6 mm. After installing a longitudinal stiffener, the original maximum deformation positions are supported by the stiffener, and the deflection of all skins are greatly reduced, as shown in Fig. 17. The maximum deflection of the flat, elliptical and ellipsoidal floors decreases by 35.8%, 41.5% and 72.9%, respectively. The maximum deflection of the ellipsoidal floor with a longitudinal stiffener is only 4.5 mm. The results indicate that under water impact conditions, the longitudinal stiffener can significantly reduce skin deformation, and the elliptical floor has better mechanical performance than flat and cylindrical floors.

Figure 16.

The vertical displacement of the skin at 3 m/s.

Figure 17.

Maximum deflection value of the skin at 3 m/s.