1 Introduction

Ever since they began to share the sky with the birds a century ago, aircrafts have been perpetually suffering from bird strikes. In fact, about 90% of all foreign object damage (FOD) can be traced to avian origins. FAA statistical data reveals that the gas turbine engines are the most vulnerable to bird strikes. In reality, almost all severe accidents are closely related to the failure of the engines, which are the sole suppliers of thrust and power to the entire aircraft. Thus, it is obvious that the problem of a bird striking a fan blade system requires careful investigation. In the 1970s, predominantly experimental techniques were used for evaluating the mechanical integrity and fail safe resistance of aircraft engine parts subjected to bird strike and ingestion (French 1974). Wilbeck (1977) found that the stresses generated during a high-speed impact greatly exceed the tissue strength of the bird, and that the bird actually behaved like a fluid with negligible viscosity. Barber et al. (1978) revealed that the loads produced by the bird strike were adequately duplicated by representing the bird as a circular cylinder with the same mass, density and compressibility as the bird. A few years later, artificial bird models began to substitute real birds (Wilbeck and Rand 1981). Today, artificial birds, usually made from gelatin, are widely adopted by many aerospace companies (Richard 2000).

With the development of advanced numerical techniques and computing technology, Aiello et al. (1982) numerically studied the mechanical performances of hollow blades with composite in-lays. They reported that these configurations would improve the crashworthiness of shroudless hollow fan blades. The centrifugal stiffening effect was also found to be beneficial in reducing blade deformation under bird strike through the preliminary dynamic analyses of Shioya and Stronge (1985) and Schuette (1990). Miyachi et al. (1991) further established that the centrifugal force effect on local deformation is smaller than that on the global scale. Furthermore, the mechanical properties of the blade were determined to be critical to assessing the strike force (Alexander 1981). Subsequently, Schuette (1990) and Martin (1990) proved that the bird properties at high-speed impact scenarios were key parameters as well. The explicit numerical studies of Schuette (1990) and Martin (1990) concerning a bird striking fan blades revealed that the bird actually behaved like compressible nonlinear fluid, and this was further verified by high-speed photography studies (Gao and Li 1990; Teichman and Tadros 1991).

It should be noted, however, that it was only in the late 1980s that the explicit time integration scheme was developed and applied to carry out large deformation simulations, such as bird strike problem. One of these works emerged from Rolls-Royce (Lawson and Turley 1987), where the large wide chord fan blade was tested for its aerodynamic performance. With the help of commercial explicit integration software incorporated with contact-impact algorithm, the bird’s slicing effect by the blades was simulated, and achieved good agreements with their respective experimental tests (Frischbier 1997; Letellier et al. 1997; Audic et al. 2000; Guan et al. 2004). A bird strike is characterized by loads with high intensity and short duration. The exposed target material undergoes high strain rates, large deformations and inelastic strains. In addition, a greater interaction exists between the impact loads and the response of the structure. With the development of explicit finite element codes in the late 1980s, it became possible to numerically analyze this kind of events with a certain degree of accuracy. The unique difficulties of bird strike simulations include: (i) the method adopted in modeling the bird, (ii) the constitutive behavior of the bird and target materials at high impinging speeds, and (iii) the geometric intricacy of the bird and the target.

It has long been noted that the finite element method (FEM) is very powerful in the numerical simulation of the crashworthiness and failure analyses. Canonical Lagrangian formulations, when applied to the bird model, will result in two difficulties; namely, reducing the time step and resulting in adversely distorted elements. When solving dynamic transient non-linear problems, the time step Δt is closely related to the physical length scale of the smallest element l min within the model, i.e.,

$$ \Updelta t = \frac{{l_{{\min }} }} {c}, $$
(1)

where c is the sonic speed in the local medium. For soft-bodied birds, large strain distortions will develop during high-speed impingement. The resulting excessively compressed meshes might decrease the time step to impractical lower values such that the analysis may take unacceptable excessive time. Another difficulty is that elemental volume may become negative thus producing negative stiffness matrix leading to failure of the analysis. To address these problems, a numerical elimination procedure for elements reaching either the maximum failure strain or the minimum reference value of the time step has been developed and employed by many researchers such as Stoll and Brockman (1997), Anghileri and Bisagni (2000), Airoldi and Cacchione (2005). A typical problem arising from element deletion is the artificial oscillations in the contact forces due to the discretised nature of the simulated contact algorithms, especially for coarse meshes. Once the frontal elements are deleted, the contact force will decrease dramatically until the impactor comes into contact again with the target, and this introduces artificial noise into the contact forces. To reduce these noises, a ‘regularization’ strategy has been implemented by Stoll and Brockman (1997) to gradually impose the contact constraint as a node descends through a buffer zone above the contact surface, and solid tetrahedral elements to model the bird. As a result, a refined mesh is highly desirable to reduce this type of artificial oscillations, although this inevitably increases the CPU time. To circumvent element deletion, Langrand et al. (2002) attempted another technique; namely, the mass scaling of the bird elements by increasing the local mass density ρ deliberately to keep the time step constant. However, this inadvertently causes the mass of the bird to become several hundred times greater than the original bird at the final time step.

In view of the recent rapid advancement in high performance computing (HPC), the strategy of employing dense meshes can be used to alleviate the problems associated with Lagrangian formulations. Modern HPC systems can tolerate very small time steps of 10−9 s or even 10−10 s, where it was used to be 10−6 s to 10−7 s just a few years ago. Thus, by refining the mesh of the Lagrangian bird, the artificial oscillations associated with the impact force are expected to be reduced to negligible levels (Mao et al. 2007).

In addition to the canonical Lagrangian formulation, there are other methods such as the Eulerian and arbitrary Lagrangian Eulerian (ALE) approaches. A very refined mesh is usually required for the Eulerian method to capture the material response, making it very time consuming. The ALE approach, on the other hand, suffers significant numerical dissipation and the solid–fluid interface is also rather difficult to describe in cases involving rotating blades. As reported by Anghileri et al. (2005), the meshes used in discretizing the bird were extremely distorted that the accuracy of the ALE approach become questionable, even when they were refined at the contact interface. Other researchers, such as Langrand et al. (2002), also compared the Lagrangian with the ALE approach for bird strike simulations and reached the same conclusion; specifically, the ALE is less accurate.

Besides the finite element method, the discrete element methods (DEMs) have also been developed and implemented for large deformation simulations, such as Feng et al. (2004) on the planetary ball milling process. Two notable DEMs, namely the nodal masses (NM) model and smooth particle hydrodynamics (SPH), have been found suitable for bird strike simulations. However, the NM method neglects the energy absorbing mechanism due to the deformation of the bird (Anghileri and Bisagni 2000). The meshless SPH method is relatively new, whose drawback is that the connectivity between particles has to be calculated repeatedly (Audic et al. 2000).

In this paper, we employ explicit 3-D finite element analysis (LS-DYNA) to carry out comprehensive studies on the transient nonlinear response of a bird striking a flexible fan blade. In connection with this, the complex and intricate geometrical configurations of different bird species have perpetually posed a problem for developing a sufficiently simplified and consistent bird model. On the other hand, it has been highly recommended by the International Birdstrike Research Group that the bird model, once standardized, should become the norm for all bird impact testing thereafter (Bowman 2004). Many researchers have simplified the bird torso as a hemispherical-ended cylinder, such as Frischbier (1997), Langrand et al. (2002), McCarthy et al. (2004), and Airoldi and Cacchione (2006). The ellipsoid geometry is also a well accepted choice which has been suggested by the International Birdstrike Research Group (Richard 2000), and has been used by Guan et al. (2004). Besides these two configurations, the straight-ended cylinder has also been adopted by Brockman and Held (1991), but its application remains somewhat infrequent. These three configurations, which are typical of artificial bird geometries, are shown schematically in Fig. 1. Nevertheless, the differences among the impact behavior associated with these different configurations have not been reported in the open literature. In our present study, the configuration effect will be examined using a canonical 4 lb bird, which strikes either a rigid panel or a typical engine fan blade. Furthermore, the effect of the length-to-diameter aspect ratio, which represents the biometric property of different bird species, is also examined in our present investigation.

Fig. 1
figure 1

The three configurations considered for the artificial bird model. (a) Straight-ended cylinder. (b) Hemispherical-ended cylinder. (c) Ellipsoid

Full size image

The literature indicates that most analyses assume normal impact. In real situations, the trajectory of the bird relative to the aero-engine fan blade is generally random, and the bird may impinge the fan blade at different angles of incidence, as shown in Fig. 2(a). In fact, even for the same engine, the incidence angle is subject to change due to different flight paths and missions. It is important to note that when the incidence angle is low, the impact scenario is more akin to that of the fan blades slicing the bird, and this is known as the slicing effect (Letellier et al. 1997; Frischbier 1997; Guan et al. 2004). In the present study, we will examine the effect of the incidence angle when a canonical 4 lb bird strikes a typical engine fan blade as shown in Fig. 2(a).

Fig. 2
figure 2

Details of finite element model of an artificial bird striking a fan blade. (a) Incidence angle of a bird model striking a fan blade. (b) Isometric view of impact system. (c) Top view of the studied blade. (d) Bottom view of the studied blade

Full size image

2 Finite element modeling

2.1 Bird properties

The mechanical properties of typical avian tissues at low speeds are neither uniform nor homogeneous. However, at progressively higher speeds, this nonuniformity and inhomogeneity become increasingly negligible, and the bird can safely be considered as a homogeneous jet of fluid impinging a structure (Wilbeck 1977). Thus, the constitutive material law of homogenized fluidic materials can be used; viz,

$$ \sigma _{{ij}} = - P\delta _{{ij}} + 2\rho \gamma \ifmmode\expandafter\dot\else\expandafter\.\fi{e}_{{ij}} . $$
(2)

There are different hydrodynamic models which have been successfully used to describe the material properties of the bird in compression and among which the most popularly used is the polynomial fitted pressure equation

$$ P = C_{0} + C_{1} \mu + C_{2} \mu ^{2} + C_{3} \mu ^{3} , $$
(3)

in which μ is the non-dimensional mass density changing ratio of the bird tissue. In this work, we model the loading from the bird as a pressure pulse on the structure. This is a feasible approach and thus we can also avoid modeling the bird disintegrating into numerous parts. The Brockman compressible modulus (Brockman and Held 1991) are employed in the present simulations, whereby

$$ \left\{ \begin{aligned}{} & C_{0} = 0 \\ & C_{1} = 2,323\,{\text{MPa}} \\ & C_{2} = 5,026\,{\text{MPa}} \\ & C_{3} = 15,180\,{\text{MPa}} \\ \end{aligned} \right. $$
(4)

The bird model is meshed with three dimensional 8-node fully-integrated solid elements (Solid-164 in LS-DYNA) with a characterized ratio between bird diameter to mesh size of 32, as shown in Fig. 1, which is found to provide efficient simulation runs without compromising on accuracy. In addition, the hourglass energy is found to be well controlled in the present simulation.

2.2 The flexible fan blade

A sector of the fan disk, composed of a single blade and a hub section, as shown in Fig. 2b and c, is used in the present simulation. The hub is assumed to be fixed and rigid in comparison with the blade. The bird impacts the fan blade at a location 85% of the radius away from the center of rotation, and the angle of attack (AOA) of the inclined airfoil at this location is 60°. The deformable blade is discretized using the two-dimensional shell elements (Shell-163 in LS-DYNA), with 31 nodes in the axial direction and 61 nodes in the radial direction found to provide sufficiently converged results. In addition, the Hughes-Liu formulation was used to eliminate the hourglass modes. The hourglass coefficient was set to 0.1, quadratic bulk viscosity to 1.5, and linear bulk viscosity to 0.06.

Fan blades of modern aeroengines are typically made of Titanium alloy Ti-6Al-4V. Due to the high strain rates associated with this problem, the constitutive law used is of the viscoplastic type originally devised by Perzyna (1968)

$$ \sigma _{y} {\left( {\varepsilon ^{P}_{{eff}} ,\ifmmode\expandafter\dot\else\expandafter\.\fi{\varepsilon }^{P}_{{eff}} } \right)} = \sigma _{y} {\left( {\varepsilon ^{P}_{{eff}} } \right)}{\left( {1 + {\left( {\frac{{\ifmmode\expandafter\dot\else\expandafter\.\fi{\varepsilon }^{P}_{{eff}} }} {C}} \right)}^{{\frac{1} {P}}} } \right)}, $$
(5)

where \( \ifmmode\expandafter\dot\else\expandafter\.\fi{\varepsilon }^{P}_{{eff}} \) is the effective plastic strain rate, and \( \sigma _{y} (\varepsilon ^{P}_{{eff}} ) \) is the initial quasi-static yield stress of the blade material. C and P are strain rate sensitive parameters determined by experiments. The magnitudes of the relevant parameters are listed as follow:

$$ \left\{ \begin{aligned}{} & E = 1.14 \times 10^{{11}} \,{\text{Pa}} \\ & \nu = 0.33 \\ & \rho _{{blade}} = 4.429 \times 10^{3} \,{\text{kg}} \cdot {\text{m}}^{{{\text{ - 3}}}} \\ & \sigma _{y} (\varepsilon ^{P}_{{eff}} ) = 1.14 \times 10^{8} \,{\text{Pa}} \\ & C = {\text{40}}{\text{.0}}\,{\text{s}}^{{{\text{ - 1}}}} \\ & P = {\text{5}}{\text{.0}} \\ \end{aligned} \right.. $$
(6)

3 Results and discussions

3.1 Rigid validation test

The aim of this test is to assess the accuracy and validity of the finite element simulations. In this test section, we selected a 4 lb bird at an initial mass density of 934 kg·m−3 striking a stationary large rigid target at normal incidence. The initial velocity of the bird was set to 225 m s−1. The bird was modeled as a cylinder with hemispherical ends, with length-to-diameter ratio of 2. Convergence tests pertaining to mesh density, integration time steps and contact were carried out, and it was found that 39,424 hexahedral elements were required to model the bird. In the convergence test carried out, when the number of elements for the bird model was refined from 26,412 to 39,424, the result in terms of global impact force changed by less than 4%.

Figure 3 shows a typical example of the velocity vectors for the hemispherical-ended cylindrical bird model striking the rigid target at a normalized time \( T = t/T_{0} = 0.05 \). The simulations began at the instant the traveling bird impinges the target. Both side and isometric views are shown in that figure. They show that the bird has been deformed, with its right tip contacting the target (Fig. 3a). The maximum velocity of the fluidic bird elements, which appears near the frontal surfaces of the bird close to the target, was found to be approximately 123 m s−1 at this instant of impact.

Fig. 3
figure 3

Velocity vectors of the bird at T = 0.05 when impacting a rigid target. (a) Side view. (b) Isometric view

Full size image

The time histories of the different phases of a bird impacting the rigid panel are depicted in Fig. 4. The hydrodynamic fluid-like behavior of the bird can be clearly observed from T = 0.5 onwards. During the impact process, the momentum of the bird will be progressively absorbed by the target; whereas the kinetic energy of the bird will be dissipated in terms of heat, as well as manifested by the elimination of some of the bird elements. Finally, the bird elements, with the umbrella-like configuration, do not possess additional forward momentum; and the impact process is completed at the normalized time T = 1. Nevertheless, following this impact event, the bird elements continue to deform and expand outwards.

Fig. 4
figure 4

Deformation history of a bird impacting a rigid target. (a) T = 0, (b) T = 0.25, (c) T = 0.50, (d) T = 0.75, (e) T = 1

Full size image

Figure 5 shows the variation of the normalized impact pressure over time. It is compared with the experimental data from the GARTEUR Bird Strike Group (Willows and Driffill 1999) as well as the numerical results of Langrand et al. (2002). It is observed that the present numerical results correspond well with the experimental data. When compared with the simulation results reported in the literature, the present result exhibits fewer oscillations than those prevalent in earlier studies, and is clearly more stable. This can be attributed to the fine mesh density used in discretizing the bird geometry. Specifically, we used bird-to-element length ratio of 32, compared with 8 used by Stoll and Brockman (1997), and Airoldi and Cacchione (2006). The peak value of the pressure from the present simulation is about 20% lower than that of the experimental data, and this difference is probably due to the use of water compressible modules for the bird’s material model, which underestimates the behavior of a real bird.

Fig. 5
figure 5

Variation of normalized impact pressure versus normalized time for a hemispherical-ended cylindrical bird impinging a rigid panel

Full size image

3.2 Effects of bird aspect-ratio

The initial mass and velocity of the hemispherical-ended bird are the same as those of the rigid impact case being 4 lb and 225 m s−1 respectively. However, besides the aspect-ratio 2.0:1, two other ratios of 1.5:1 and 2.5:1 are also examined and compared. Table 1 lists all diameters of artificial birds at these three aspect ratios. For the rigid impact cases, the normalized impact pressures are plotted in Fig. 6a. Their peak values are much the same, being 2.75 for aspect ratio 1.5:1, 2.62 for aspect ratio 2.0:1, and 2.61 for aspect ratio 2.5:1. The impulses for the three aspect ratios are shown in Fig. 6b. It is found that the impulse, after normalization, is 0.29 for the aspect ratio 1.5:1, 0.23 for the aspect ratio 2.0:1, and 0.26 for the aspect ratio 2.5:1.

Table 1 Details of the different bird aspect ratios considered
Full size table
Fig. 6
figure 6

Effect of aspect-ratio of artificial bird in rigid impact case. (a) Normalized impact pressures. (b) Impact impluse

Full size image

3.2.1 Flexible impact at aspect-ratio 2.0:1

For the flexible fan blade impact case, firstly the case of aspect ratio 2.0:1 has been examined. Sample time histories of the deformation of the blade and the bird for this case are shown in Fig. 7. The figure shows the evolution of the bird deformation during impact. The severity of the localized contact at the initial stages of impact will lead to plastic deformation in the considered blade. The displacements, velocities and accelerations in X-, Y- and Z-directions of one blade point in the near field of the impact region are shown in Fig. 8. As expected, the displacements, velocities and accelerations in X- and Y-directions are much lower than those in the Z-direction. The maximum velocity of this near field point in the Z-direction is as high as 32% of the initial bird velocity. The maximum instantaneous acceleration is found to be 1.6 × 106 m s−2 occurring just after impact and decreases quite rapidly.

Fig. 7
figure 7

Deformation history for a hemispherical-ended cylindrical bird impinging a flexible fan blade. (a) T = 0, (b) T = 0.25, (c) T = 0.50, (d) T = 0.75, (e) T = 1, (f) T = 1.34

Full size image
Fig. 8
figure 8

Time histories of displacement, velocity and acceleration of near field point. (a) Displacements in X-, Y-, Z-directions. (b) Velocities in X-, Y-, Z-directions. (c) Accelerations in X-, Y-, Z-directions

Full size image

Three points on the blade, corresponding to three representative distances to the point of impact, have been designated for a more detailed examination. These points are the near field, the middle field and the far field impact points, and their coordinates are provided in Fig. 2b. Figure 9 shows the time histories of the normal, shear and von Mises stresses of the near field point. Due to the bending of the blade, the normal stresses in the Y-direction are of the highest magnitude; whereas those in the X- and Z-directions are considerably lower. The maximum normal stress in the Y-direction is found to be about 295 MPa. The shear stresses are almost one order lower than the normal stresses at that point. Consequently, the von Mises stresses reflect the influence of the dominant σ y stress. Figure 10 shows plots of the von Mises plastic strains of the three representative points. The von Mises plastic strain of the near field point is several times larger than that of the middle field point, which is in turn larger than that of the far field point. The maximum von Mises plastic strain is about 0.028.

Fig. 9
figure 9

Normal, shear and von Mises stresses of near field point. (a) Normal stresses in X-, Y-, Z-directions. (b) Shear stresses in XY, YZ, XZ planes. (c) von Mises stress

Full size image
Fig. 10
figure 10

von Mises plastic strain histories of near, middle and far field points

Full size image

The impact force and impulse variations in the X-, Y- and Z-directions are shown in Fig. 11. As in the case of the respective velocities and accelerations, the impact forces in the X- and Y-directions are much lower than the impact force in the Z-direction, which is almost of the same magnitude as the total impacting force. As such, in the subsequent results presented, only the impact force in the Z-direction is considered. The maximum normalized impact force in that direction is found to be 1.48 which is equivalent to about 360 kN. From the impulse plot, we also observe that the momentum is transmitted from the bird to the blade throughout the impact duration which is measured in milliseconds. Due to the fixed boundary conditions at the root of the blade, losses result and only 37.5, 8.0 and 1.7% of the initial momentum of the bird are transmitted to the blade in the Z-, X- and Y-directions, respectively. The present simulation is most successful in suppressing the stairway profiles during the process of the impulse accumulation, which was observed in the simulation of Stoll and Brockman (1997) due to the presence of undesired spiky contact force histories in their results.

Fig. 11
figure 11

Normalized impact pressure (a) and impulse (b) variations in X-, Y-, Z-directions in fan blade impact case

Full size image

The kinetic, deformation and hourglass energies of the bird, blade and the system that is composed of both the bird and blade, have been plotted in Fig. 12. The kinetic energy is transmitted from the bird to the flexible blade. Meanwhile, the total kinetic energy decreases, as it is transformed into deformation energy of the flexible blade, as can be seen from Fig. 13. The deformation energy of the blade is found to be much larger than that of the bird, due to the inelastic fluid-like properties of the bird elements. The hourglass energy of the bird is found to be significantly higher than that of the blade. In the present analysis, the hourglass energy of the whole system is only 6% (or less) of the total energy.

Fig. 12
figure 12

Energy distributions for the bird, blade and bird-blade system. (a) Kinetic energy. (b) Deformation energy. (c) Hourglass energy

Full size image
Fig. 13
figure 13

Total, kinetic, deformation and hourglass energies of the bird-blade system in fan blade impact case

Full size image

3.2.2 Flexible impact at different aspect-ratios

Figure 14a illustrates the evolution of the normalized impact pressure between the fan blade and hemispherical-ended cylindrical birds with the three different aspect ratios. Since their configurations are similar, their evolution histories are also found to be moderately comparable. For all of the cases considered, there are three significant peaks within their impact profiles, and all of them reach their respective maximum values at the second peak, at magnitudes of 1.66 for aspect ratio 1.5:1, 1.78 for aspect ratio 2.0:1, and 1.70 for aspect ratio 2.5:1. The relative differences are less than 7%. Next the normalized impulse variations for the three aspect ratios are shown in Fig. 14b. It can be seen that the impulse evolutions are quite similar as well. Within the normalized time T = 0 to 0.6, the impulse increases in an almost linear manner. The rate of increase in the impulse plateaus levels off after T = 0.6. When the impact force reaches zero, we note that 34, 37 and 32% of the initial momentum of the bird has been transmitted to the fan blade, for the three respective aspect ratios. Thus, we can draw from this study that the aspect ratio does not have significant effects on the results of normalized impact pressure and impulse, both for the impact case of the rigid panel and flexible fan blade.

Fig. 14
figure 14

Effect of aspect-ratio of artificial bird in fan blade impact case. (a) Normalized impact pressure. (b) Normalized impluse

Full size image

3.3 Effects of bird configuration

The initial mass and velocity of the bird are the same as those of the rigid impact case being 4 lb and 225 m s−1 respectively. However, besides the hemispherical-ended cylinder, two other configurations including the straight-ended cylinder and ellipsoid are also examined and compared. Table 2 lists the detailed parameters of the three configurations at aspect ratio of 2.0:1.

Table 2 Details of the different bird configurations considered
Full size table

3.3.1 Rigid impact effects

Figure 15a shows the rigid impact force variations for different configurations of the bird model. It is found that the straight-ended cylinder reaches its maximum impact force of 799 kN at about 0.049 ms after initial contact, and there is only one dominant peak in its contact force profile. However, there are two quite distinct peaks corresponding to the hemispherical-ended cylinder and ellipsoid models. The hemispherical-ended cylinder reaches its maximum force of 527 kN at 0.054 ms, and its second peak of 333 kN at 0.134 ms, although this second peak is not as significant as the first one. The ellipsoid model, on the other hand, reaches its first peak of 360 kN at 0.111 ms; and its second peak of 378 kN at 0.28 ms. Interestingly, this second peak is found to be the maximum magnitude in its force variation time history. Among the three different bird configurations, the maximum impact force for the straight-ended cylindrical bird is the highest. This is because the instant the three types of bird configurations impinging the rigid target, the straight-ended cylindrical bird has comparatively the largest instantaneous contact area. This also explains why there is only one significant peak in its force variation, and it is because all the elements of this bird model will decelerate immediately and dramatically once it impacts the target. However, this situation is quite different for the hemispherical-ended cylinder and ellipsoidal birds. Since the bird in the high-speed impact scenarios behaves as a liquid, the shear stress between the bird’s elements flowing on neighboring streamlines, compared with the local high pressure, is comparatively negligible. As pointed out earlier, this has been deduced by Wilbeck (1977), Barber et al. (1978), and Gao and Li (1990) in their individual experiments. So our present investigation treats the bird hydro-dynamically, and neglects the fluidic viscous effect of the bird tissue. Therefore, for the hemispherical-ended cylinder and ellipsoidal birds, only the elements directly downstream of those in contact with the target will experience the dramatic deceleration; whereas the other elements will not decelerate significantly until their fronts come into contact with the target. Thus, when impact progresses for the hemispherical-ended cylinder or ellipsoidal birds, the contact forces decrease after the initial peaks but will increase again, leading to their respective second peaks. This second peak is quite distinct in the ellipsoidal case, but rather less pronounced in the hemispherical-ended cylindrical bird case. In the ellipsoid model, the contact area reaches its maximum value somewhat late into the impact event. It occurs at approximately the midpoint of the entire impact duration, which explains the presence of its second peak.

Fig. 15
figure 15

Effect of configuration of artificial bird in rigid impact case. (a) Impact forces. (b) Normalized impact pressures

Full size image

The normalized impact pressure evolutions for the three bird configurations are plotted in Fig. 15b. It should be noted that because the cross-sectional areas vary along the geometrical axes of the hemispherical-ended cylinder and the ellipsoid models, the normalization is based on their respective averaged cross-sectional areas. The maximum normalized impact pressure is found to be 3.73 at T = 0.055 for the straight-ended cylinder, 2.61 at T = 0.053 for the hemispherical-ended cylinder, and 2.00 at T = 0.27 for the ellipsoid model. Thus the maximum impact pressure for the straight-ended cylinder impact case is about 43% higher than the hemispherical-ended cylinder case; and the latter is about 30% higher than the ellipsoid case.

3.3.2 Flexible impact effects

For the impact with flexible fan blade, the normalized impact pressures between the birds and fan blade are shown in Fig. 16a for all the three bird configurations. The dominant peaks are 1.75 at T = 0.094, 1.35 at T = 0.189, and 0.77 at T = 0.32 for the straight-ended cylindrical bird; and 1.62 at T = 0.06, 1.78 at T = 0.16, and 1.46 at T = 0.45 for the hemispherical-ended cylindrical bird. For the ellipsoidal bird, the peaks are 0.79 at T = 0.06, 1.53 at T = 0.16 and 1.34 at T = 0.39. From this figure it can be seen that there are multiple peaks in the impact force evolution, which is mainly due to the coupling between the bird and the flexible blade. The blade, upon being impacted by the bird accelerates in the same direction as the bird’s initial velocity, i.e., the Z-direction. When this happens, the frontal elements of the bird which are in close contact with the blade would expand and their high pressures would thus decrease. Therefore, the total impact force would also decrease dramatically. Consequently, with the lowered impact force from the bird, the bending of the blade would decelerate, and as a result the contact between the bird and the blade becomes intense once more. The impact force between the bird and blade would increase again. The whole procedure will repeat several times before the impact force eventually vanishes. It is also interesting to note that all the three types of bird configurations have three significant peaks within their variations of the impact forces. This is mainly because the bird mass and density are similar for the three configurations. By varying the mechanical properties of either the bird or the blade, including their masses, densities and/or other properties, it can be expected that the number of significant peaks for the impact force may change.

Fig. 16
figure 16

Effect of configuration of artificial bird in fan blade impact case. (a) Normalized impact pressure. (b) Normalized impluse

Full size image

On the other hand, we can observe from Fig. 16a that the normalized impact pressure profiles for the three bird configurations are quite different from the preceding rigid impact cases. This is especially so for the straight-ended cylindrical bird, which no longer possesses the highest peak impact pressure value among the three bird configurations. As there is an attack angle of 60° between the bird trajectory and the blade (Fig. 2), the initial contact area between the straight-ended cylindrical bird and the blade is not as large as the rigid panel case. Therefore, the initial contact force between the bird and target is also significantly reduced. More importantly, the blade is not rigid but flexible, and the maximum impact force between the straight-ended cylindrical bird and blade is found to be only 47% of the rigid panel impact case. Similar conclusions have been reached by Shimamura et al. (2004) that the target’s flexibility strongly influences the impact force. The impact force resulting from the hemispherical-ended cylindrical bird impinging the blade is about 68% of the former rigid impact case; and correspondingly 76% for the ellipsoidal bird model. Comparing the three configurations, it is noted that the ellipsoidal bird model has a maximum impact force which is lower than the other two configurations by about 15%.

The normalized impulse between the bird and blade is shown in Fig. 16b. The momentum of the bird is transmitted to the fan blade over the duration of the impact process. However, only 24.4, 37.5 and 40.2% of the initial momentum of the bird are finally transmitted to the blade in the Z-direction, for the straight-ended cylinder, hemispherical-ended cylinder and ellipsoidal bird, respectively.

3.4 Effects of angle of incidence

In this section, the modeling of the artificial bird and its initial velocity of 225 m/s are the same as those of the rigid impact case. For the 3-dimentionally contoured fan blade, its tangential surface at the impact point is found to be approximately perpendicular to the Z-direction. The incidence angle θ is varied from 60°, 30°, 0° (normal incidence), to −30°.

Figure 17 shows the variations of the impact force for the four incidence angles considered in this investigation. The force variations in the case where θ = 60° are quite small relative to other angles. Thus, the slicing force of the blade at this high incidence angle would be much smaller than the other cases. For the case of θ = 30°, the impact process has the relative longest duration of about 2.2 ms. The normal incidence case of θ = 0° has the maximum impact force values among the four cases considered. The force variation profile of the case of θ = −30° is very similar to the normal incidence case except that its impact force magnitudes are generally lower by about 20%.

Fig. 17
figure 17

Evolution of impact forces for various incidence angles

Full size image

The maximum impact force and maximum plastic strain of the blade are plotted in Fig. 18 for the four cases. For the case of θ = 60°, there is no plastic strain sustained by as the impact forces are too low to generate any permanent deformations. For the case of θ = 30°, the maximum plastic strain of the blade remains small even though the impact forces for this case increased significantly compared to the case of θ = 60°. This is because the blade yield stress (σ y ) increases with the increase of the impact forces due to the higher strain rate. Considering the four cases, both the maximum impacting force and plastic strain reach their respective peak values at θ = 0°. The results for the case of θ = −30° are very similar to those of θ = 60°, with both displaying considerable impact force values, but much smaller plastic strains than that of the normal incidence case.

Fig. 18
figure 18

Maximum impact force and plastic strain values of the blade for different incidence angles

Full size image

Figure 19 shows the final displacements of nodes 3, 800 and 1831 for the four cases considered. As there is no plastic strain for the case of θ = 60°, the final displacements for all the nodes are zero. The displacements are found to be of maximum values for the normal incidence case of θ = 0°. Figure 20 shows the final deformation modes of the blade for the four cases. Again, there is no permanent deformation for the case of θ = 60°. It is quite clear that the final deformation of the normal impact case is the most significant among the four different cases. It should be noted that the shearing forces have been neglected in the present investigation as the fluidic bird model is assumed to be non-viscous. This assumption is quite reasonable when considering the bird impinging the blade at small incidence angles. However, for the slicing effect when the incidence angle is large, the shear force will play a significant or even dominant role. Thus, refinements to the numerical model and further experimental data will be necessary.

Fig. 19
figure 19

Final displacements in the Z-direction at typical nodal positions

Full size image
Fig. 20
figure 20

Final deformation modes for different incidence angles

Full size image

3.5 Effects of bird impact velocity

The initial impact velocity of the bird was varied from 83 m s−1 to 318 m s−1, with the properties and modeling of the bird being the same as those of the rigid impact case. The velocities considered for the bird cover the normal velocity range of commercial airplanes in flight.

The normalized impact forces are shown in Fig. 21 for the respective velocities 83, 225 and 318 m s−1. The maximum normalized forces for these three cases are 1.59, 1.48 and 1.36, respectively. For the case of \( \ifmmode\expandafter\dot\else\expandafter\.\fi{w}_{0} \) = 83 m s−1, the impact force drops to relatively low magnitudes at about normalized time T = 0.7, and vanishes at T = 1.08. It is also noted that there are more oscillations in the impact force compared to the other two higher velocity cases considered. As for the case of \( \ifmmode\expandafter\dot\else\expandafter\.\fi{w}_{0} \) = 318 m s−1, the force reaches its maximum value and then gradually decreases while slightly oscillating until it vanishes at T = 1.54. The impulses resulting from these three velocity cases are plotted in Fig. 22. The impulse variations for the three cases considered exhibit the same trend. Since the mass of the bird remains constant, normalization of the impulse equations effectively meant the division by \( \ifmmode\expandafter\dot\else\expandafter\.\fi{w}_{0} \) . Accordingly, when the impulse values are normalized by their corresponding initial impact velocity \( \ifmmode\expandafter\dot\else\expandafter\.\fi{w}_{0} \), the higher velocities result in overall lower normalized impulse curve.

Fig. 21
figure 21

Normalized impact forces for different impact velocities

Full size image
Fig. 22
figure 22

Normalized impulse values for different impact velocities

Full size image

It is also found that the maximum impact force increases in a quadratic manner with velocity. Their relationship can be approximated as

$$ F_{{\max }} = 5.439\ifmmode\expandafter\dot\else\expandafter\.\fi{w}_{0} ^{2} + 431.8\ifmmode\expandafter\dot\else\expandafter\.\fi{w}_{0} - 24830,\;\;{\left( {318{\text{ m}}\,{\text{s}}^{{ - 1}} \ge \;\ifmmode\expandafter\dot\else\expandafter\.\fi{w}_{0} \; \ge \;83{\text{ m}}\,{\text{s}}^{{ - 1}} } \right)}. $$
(7)

It should be noted that this relationship is not general but specific to the presently described bird strike case. In addition, with the increase of the impact velocity, the maximum plastic strain of the blade increases as well due to higher striking forces from the bird, and follows the following relationship for the presently considered case

$$ \varepsilon _{{p,\max }} = (0.07\ifmmode\expandafter\dot\else\expandafter\.\fi{w}_{0} ^{2} + 7\ifmmode\expandafter\dot\else\expandafter\.\fi{w}_{0} - 270) \cdot 10^{{ - 5}} ,\;\;{\left( {318{\text{ m}}\,{\text{s}}^{{ - 1}} \ge\; \ifmmode\expandafter\dot\else\expandafter\.\fi{w}_{0}\; \ge \; 83{\text{ m}}\,{\text{s}}^{{ - 1}} } \right)}. $$
(8)

3.6 Effects of bird size

The total mass of the bird is 4 lb (or 1.82 kg), which is presently used as the upper limit of the bird mass criteria in the bird-aeroengine strike scenarios (Paul 2000). Nevertheless, this standard inevitably leaves out many large flocking bird species. To investigate the effects of heavier bird strikes, 6 lb (2.72 kg) and 8 lb (3.63 kg) birds are considered here. The mass density of bird models, with feathers removed, has been enumerated by the databases of the International Birdstrike Research Group (Seamans et al. 1995) as

$$ \rho _{0} = 1148 - 63 \times \log _{{10}} (1000 \times m) = 959 - 63 \times \log _{{10}} (m), $$
(9)

where m and ρ 0 are the respective mass and initial density of the bird in SI units. The diameter of the bird model is given by

$$ D = 2{\left( {\frac{m} {{2\pi \rho _{0} }}} \right)}^{{\frac{1} {3}}} {\left( {\frac{L} {D} - \frac{1} {3}} \right)}^{{ - \frac{1} {3}}} , $$
(10)

where L/D is the length-to-diameter aspect ratio. Based on the relationships of Eqs. 9 and 10, the lengths and diameters of 4 lb (1.82 kg), 6 lb (2.72 kg) and 8 lb (3.63 kg) birds can be obtained and these data are listed in Table 3.

Table 3 Large bird properties
Full size table

Hemispherical-ended cylinder models with aspect ratio 2.0:1 are employed to simulate the birds. The initial velocity of the birds is taken to be 225 m s−1 and impinging at normal incidence. As expected, the maximum impact force increases with the mass of the bird (Fig. 23). However, it is also noted that the increase rate decreases slightly due to the lower density effect dictated by Eq. (9). In addition, the maximum plastic strain of the blade struck by a 6 lb (2.72 kg) bird is found to be even less than that experienced by the 4 lb (1.82 kg) bird. This again can be attributed to the reduced density of the heavier bird, which decreases logarithmically with respect to the mass of the bird.

Fig. 23
figure 23

Maximum impact force and plastic strain for birds with different masses

Full size image

4 Conclusions

With the advent of large fan engines with increased frontal areas and no inlet guide vanes, the understanding of the bird strike problem and its proper mitigation now take on even more significance. While manufacturers undertake a great deal of in-house testing prior to certification, it is also highly desirable to have predictive numerical models to assess the bird impact resistance of new generation engine fan blades. To validate the present finite element formulation, comparison with reported experimental data for a bird striking a rigid target was carried out. The Lagrangian bird, which was modeled with a hydrodynamic constitution law, is found to be appropriate for the bird strike simulation.

To address the diversities of bird geometry modeling, three different frequently used configurations and aspect ratios, representing diverse biometric bird species, have been examined. It is found that the initial contact area between the bird and target in the early phase of the impact event would have a significant effect on the peak impact force value. For the bird model impacting a rigid target, the maximum impact force for the case of a straight-ended cylindrical bird is about 43% higher than that of a hemispherical-ended cylindrical bird, which in turn is 30% higher than that of an ellipsoidal bird. The impact force profile is also found to be highly dependent on the deformation of the flexible fan blade. On the other hand, the length-to-diameter aspect ratio of the bird striking both rigid panel and flexible fan blade is found to have little influence on the impact results.

As to the influence of the incidence angle, it is found that the maximum impact force exerted by the bird on the fan blade, as well as the accumulated plastic strain, would reach peak values at normal incidence. Due to the increase in the yield stress of the blade material at higher strain rates, decrease in the incidence angle will result in more pronounced increases in the impact force values, as compared to less pronounced increases in plastic strains. The maximum impact force increases with larger birds, but this trend is not necessarily the case for the maximum plastic strain of the blade due to the reduced density of the larger bird.