# Composite blade damaging under impact

Thomas Menouillard, Julien Réthoré, Harridh Bung, Arnaud Suffis

► To cite this version:

Thomas Menouillard, Julien Réthoré, Harridh Bung, Arnaud Suffis. Composite blade damaging under impact. Journal de Physique IV Proceedings, 2006, 134, pp.409-415. 10.1051/jp4:2006134063. hal-00938636

HAL Id: hal-00938636

https://hal.science/hal-00938636

Submitted on 4 Jun 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Composite blade damaging under impact

T. Menouillard<fnref n="1" />,<fnref n="2" />, J. Réthoré<fnref n="3" />, H. Bung<fnref n="2" /> and A. Suffis<fnref n="4" />

**Abstract.** Composites materials are now being used in primary aircraft structures, and other domains because of numerous advantages. A part of a continuous in-flight operating costs, gas turbine engine manufacturers are always looking for ways to decrease engine weight. This is the case of compressor blades which have to satisfy, for example, the standard bird strike or debris in order to measure the crashworthiness. Bird strike impacts are actually among the most challenging loads that composite blades must accommodate. Thus for the further development of composite structures, it becomes important to have available predictive tools for simulating the response of composite structures under crash or impact loads, which will allow to evaluate damage state in the structure in function of time. A composites damage model, without mesh dependency, is presented, and allows to obtain agreement with impact experiment. Examples of finite element simulations for the impact response of blade based on this materials model are developed. These numerical results correspond to a bird strike on an equivalent composites blade, and insists on damage evolution in structure.

## 1. INTRODUCTION

This paper will focus on dynamic composite damage modelling. First two different damage models will be presented: a linear strain damage model, and the Onera Damage Model (denoted ODM) [1]; these two models are developed for orthotropic materials, hence several damage parameters are useful for description. The Onera Damage Mechanics model is richer to describe composite material behaviour. Second a dynamic numerical problem for damage calculus is the mesh dependancy, which has been shown by O. Allix in [2]. Third numerical simulations will improve stability of models and accurate estimate damaged zone after impact.

## 2. COMPOSITE MATERIAL: DAMAGE MODELS

The composite material is modeled at the macroscopic scale, as a continuum. Let $S_0$ and $C_0$ be the compliance and stiffness tensors (6×6 matrix), respectively, of the initial material; $S_{eff}$ and $C_{eff}$ the elasticity tensors of the damaged material with all damage modes. The strain is denoted $\epsilon$, elastic strain $\epsilon_e$, stress $\sigma$, and the damage parameter $d$. The orthotropic characteristics make having several damage parameters, denoted $d_i$ with i being 1 to number of damage parameters.

### 2.1 Linear strain damage model

This part presents a linear strain damage model used for a composite. The damaged stiffness is expressed as: $E = (1 - d)E_0$.

And the damage law is function of strain as:

$$d = \begin{cases} d_c \cdot \left\langle \frac{\epsilon - \epsilon_0}{\epsilon_c - \epsilon_0} \right\rangle_+ & \text{if : } \epsilon \le \epsilon_c \\ d_c & \text{otherwise} \end{cases}$$

<footnote n="1">¹ LAMCOS, INSA Lyon, 18 rue des Sciences, 69621 Villeurbanne, France</footnote>

<footnote n="2">² CEA Saclay, DEN/DM2S/SEMT/DYN, 91191 Gif-sur-Yvette, France</footnote>

<footnote n="3">³ Faculty of Aerospace Engineering, Delft University of Technology, 2600 GB Delft, The Netherlands</footnote>

<footnote n="4">⁴ SNECMA, W/YKM3, Rond point René Ravaud, 77550 Moissy-Cramayel Cedex, France</footnote>
<figure><img src="image_2.png" /><figcaption>Figure 1. Evolution of linear damage as a function of strain.</figcaption></figure>

$$ \text{where } \langle X \rangle_{+} = \begin{cases} X & \text{if: } X \ge 0 \\ 0 & \text{else} \end{cases} $$

On figure 1 one plots the damage as a function of the strain; three parts are underlined: zero damage, a linear evolution from zero to critical value, and critical damage.

This model has been coded in the fast dynamics code EUROPLEXUS with 6 (one per strain tensor component) damage parameters ($d_1, d_2, d_3, d_{12}, d_{13}$ and $d_{23}$) using this type of law. The only three parameters of the model ($d_c, \varepsilon_0$ and $\varepsilon_c$) do not seem to allow to describe as well as expected composite material behaviour.

## 2.2 Onera Damage Model

The Onera Damage Mechanics model has been developed by ONERA [3, 4], and more particularly by J.-F. Maire [1, 5, 6] and N. Carrere. This model is supposed to described Ceramic matrix Composites (CMC) materials [7], but can easily be used for other composite materials. It actually takes into account:

* The initial anisotropy (i.e. the effect of damage on the elastic behaviour, with $H_0^i$);

* Defect growth law: damage evolution in function of fracture mode (with the variable $d_i$);

* Desactivation mechanisms take place when a component of strain vanishes and changes sign. The closing criterion is described by the variables $\eta_i$;

* Residual strain, which can be due to partial closing of defect even if external loading is suppressed, denoted by the variables $\varepsilon_s$ and $\varepsilon_r$.

Figure 2 describes the two last phenomena. The model we have been using is the strain ODM scalar model [8]. This means that fracture modes are linked to the anisotropic directions of the material (3 scalar parameters $d_1, d_2$ and $d_3$). Moreover, an essential feature of this model is its explicit character that makes easier its implementation in an explicit code.

<figure><img src="image_3.png" /><figcaption>Figure 2. Stress strain relationship for an ODM model.</figcaption></figure>
### 2.2.1 An active/passive damage model

This model takes into account the damage desactivation. This information of active/passive damage is described in the variable $\eta_i$ (for each damage). Figure 2, which presents a tension-compression test, shows that after a damaging load, we observe a complete damage desactivation ($\eta_i = 0$) characterised by a compressive Young's modulus equal to the initial value. The damage-effect tensors are modeled in $\mathbf{H}_i^0$; indeed figure 2 shows that compression in direction 1 makes growing transverse damages. The damaged elasticity tensor $\mathbf{C}_{\text{eff}}$ is function of the initial one $\mathbf{C}_0$ and the damage-effect tensors $\mathbf{H}_i^0$:

$\mathbf{C}_{\text{eff}} = \mathbf{S}_{\text{eff}}^{-1} = (\mathbf{S}_0 + \sum_i \eta_i \cdot \mathbf{d}_i \cdot \mathbf{H}_i^0)^{-1} = (\mathbf{C}_0^{-1} + \sum_i \eta_i \cdot \mathbf{d}_i \cdot \mathbf{H}_i^0)^{-1}$ where the matrix $6 \times 6 \mathbf{S}_0$ and $\mathbf{H}_i^0$ are defined as:

**Compliance** $\mathbf{S}_0$: $\mathbf{S}_0(1, 1) = 1/E_1^0$, $\mathbf{S}_0(2, 2) = 1/E_2^0$, $\mathbf{S}_0(3, 3) = 1/E_3^0$; $\mathbf{S}_0(1, 2) = -\nu_{12}/E_1^0$, $\mathbf{S}_0(1, 3) = -\nu_{13}/E_1^0$, $\mathbf{S}_0(2, 3) = -\nu_{23}/E_2^0$; $\mathbf{S}_0(4, 4) = 1/G_{12}^0$, $\mathbf{S}_0(5, 5) = 1/G_{23}^0$, $\mathbf{S}_0(6, 6) = 1/G_{13}^0$, and 0 otherwise.

**Damage effect tensors** $\mathbf{H}_0^1, \mathbf{H}_0^2, \mathbf{H}_0^3$:

$\mathbf{H}_0^1$: $\mathbf{H}_0^1(1, 1) = h_1^n/E_1^0$, $\mathbf{H}_0^1(4, 4) = h_1^p/G_{12}^0$, $\mathbf{H}_0^1(6, 6) = h_1^{pn}/G_{13}^0$, and 0 otherwise.

$\mathbf{H}_0^2$: $\mathbf{H}_0^2(2, 2) = h_2^n/E_2^0$, $\mathbf{H}_0^2(4, 4) = h_2^p/G_{12}^0$, $\mathbf{H}_0^2(5, 5) = h_2^{pn}/G_{23}^0$, and 0 otherwise.

$\mathbf{H}_0^3$: $\mathbf{H}_0^3(3, 3) = h_3^n/E_3^0$, $\mathbf{H}_0^3(5, 5) = h_3^p/G_{23}^0$, $\mathbf{H}_0^3(6, 6) = h_3^{pn}/G_{13}^0$, and 0 otherwise.

### 2.2.2 Computational model

This part focuses on the algorithm of the ODM scalar model [8]. From knowing strain tensor $\underline{\varepsilon}$ at time $n + 1$, it aims at evaluating stress $\sigma$ at time $n + 1$. First, positive strain tensor has to be extracted; so

$\underline{\varepsilon}^+ = \underline{P} \cdot \langle \underline{\varepsilon}^{\text{diag}} \rangle_+ \underline{P}^T = \underline{P} \cdot \begin{bmatrix} \langle \varepsilon_I \rangle_+ & 0 & 0 \\ 0 & \langle \varepsilon_{II} \rangle_+ & 0 \\ 0 & 0 & \langle \varepsilon_{III} \rangle_+ \end{bmatrix} \underline{P}^T$ where $\varepsilon_I, \varepsilon_{II}$ and $\varepsilon_{III}$ are the eigen values of strain tensor.

Then, the thermodynamic forces $y_i^n, y_i^t$ are computed in function of positive strain as following: $y_i^n = 1/2 \mathbf{C}_{ii}^0 \cdot \varepsilon_i^+ \cdot \varepsilon_i^+$ for $i \in \{1, 2, 3\}$ and $y_i^t = 1/2 (\mathbf{b}_1 \cdot \mathbf{C}_{66}^0 \cdot \varepsilon_{13}^+ \cdot \varepsilon_{13}^+ + \mathbf{b}_2 \cdot \mathbf{C}_{44}^0 \cdot \varepsilon_{12}^+ \cdot \varepsilon_{12}^+)$, $y_i^t = 1/2 (\mathbf{b}_3 \cdot \mathbf{C}_{55}^0 \cdot \varepsilon_{23}^+ \cdot \varepsilon_{23}^+ + \mathbf{b}_4 \cdot \mathbf{C}_{44}^0 \cdot \varepsilon_{12}^+ \cdot \varepsilon_{12}^+)$, $y_i^t = 1/2 (\mathbf{b}_5 \cdot \mathbf{C}_{66}^0 \cdot \varepsilon_{13}^+ \cdot \varepsilon_{13}^+ + \mathbf{b}_6 \cdot \mathbf{C}_{55}^0 \cdot \varepsilon_{23}^+ \cdot \varepsilon_{23}^+)$.

Afterwards, the damage law of the scalar model is: $d_i = g_i^n(y_i^n) + g_i^t(y_i^t)$ where

$$g_i^n(y_i^n) = d_{c(i)}^n \left( 1 - e^{-\left(\frac{\sqrt{\varepsilon_i^t} - \sqrt{\varepsilon_{6(i)}^t}}{\sqrt{\varepsilon_{c(i)}^t}}\right)^{\mu}} \right)$$

and $g_i^t$ is the same function too.

Then, desactivation parameters are $\eta_i$:

$$\eta_i = \begin{cases} \frac{1}{2} \left( 1 - \cos \left( \frac{\pi \bar{\varepsilon}_i + \Delta \varepsilon_i^f}{2 \Delta \varepsilon_i^f} \right) \right) & \text{if: } \Delta \varepsilon_i^f \le \bar{\varepsilon}_i \\ 0 & \text{if: } -\Delta \varepsilon_i^f \le \bar{\varepsilon}_i \le \Delta \varepsilon_i^f \\ 0 & \text{if: } \bar{\varepsilon}_i \le -\Delta \varepsilon_i^f \end{cases}$$

where $\Delta \varepsilon_i^f = (1 + a_i^f d_i) \Delta \varepsilon_i^0 (d_c^n(i), p^n(i), y_o^n(i), y_c^n(i), d_c^t(i), p^t(i), y_o^t(i), y_c^t(i), a_i^f, \Delta \varepsilon_i^0$ are parameters).
<figure><img src="image_4.png" /><figcaption>Figure 3. Damage evolution as a function of strain for the ODM model.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Figure 4. Evolution of desactivation parameter in function of principal strain.</figcaption></figure>

Then the incremental residual strains are:

$\Delta\varepsilon_s = \mathbf{S}_0 : (\sum_i \Delta\eta_i \cdot d_i \cdot \mathbf{C}_{eff} : \mathbf{H}_i^0 : \mathbf{C}_{eff}) : \bar{e}$ and $\Delta\varepsilon_r = \mathbf{S}_0 : (\sum_i \zeta_i \cdot \eta_i \cdot \Delta d_i \cdot \mathbf{C}_{eff} : \mathbf{H}_i^0 : \mathbf{C}_{eff}) : \bar{e}$ where $\bar{e} = \varepsilon - \varepsilon_0$ ($\varepsilon_0$ allows to take into account residual strain of material forming, which can depend on temperature). New residual strains can also be computed: $\varepsilon_s^{n+1} = \varepsilon_s^n + \Delta\varepsilon_s$ and $\varepsilon_r^{n+1} = \varepsilon_r^n + \Delta\varepsilon_r$, and eventually the new stress becomes: $\sigma^{n+1} = \mathbf{C}_{eff} \cdot \varepsilon_e = \mathbf{C}_{eff} \cdot (\varepsilon_s^{n+1} - \varepsilon_r^{n+1})$, so $\Delta\sigma = \sigma^{n+1} - \sigma^n$.

### 3. MESH DEPENDENCY: THE DELAY EFFECT

The idea of a limitation of the damage rate was the basis of the delayed damage model proposed by the LMT Cachan within the framework of the studies on the rupture of composite laminates [2, 9]. In this model, two parameters are introduced. One is the characteristic time $\tau_c$ representing the inverse of the maximum damage rate. The equation of the damage rate evolution reads:

$$ \left\{ \begin{array}{ll} \frac{\partial d}{\partial t} = \frac{1}{\tau_c} (1 - e^{-a(D_{nc} - d)_+}) & \text{if : } d < d_c \\ d = d_c & \text{otherwise} \end{array} \right. $$

where a is the second parameter of the model, and $D_{nc}$ the classical definition of damage without delay effect. However the determination of these two parameters ($a, \tau_c$) has to be done for each damage parameter (i.e. 3 for ODM model). For an isotropic material, determination of characteristic length is described in reference [10].

### 4. NUMERICAL SIMULATIONS

The different numerical simulations are static tension-compression, static compression-tension, dynamic tension/compression and impact on a blade. The first and second cases allow to validate the numerical results (stress, strain, damages) in static with the damage law. The third case allows to validate results in dynamic, by taking into account mesh dependency. And to finish, the last case is the application of the damage calculus for an impact on a blade, which could be compared to a bird strike for example. The orthotropic material axes (denoted 1, 2, 3) are the geometrical axes (x, y, z) in this whole case.
## 4.1 Compressive/tensile cases

### 4.1.1 Static case

Figure 5 presents the experiment and the boundary conditions. Figure 6 shows the two loadings used: loading 1 is tension then compression, and loading 2 is compression then tension. Figures 7 and 8 shows damage as a function of time for both loading cases. As a result it underlines the transverse effect of damage.

### 4.1.2 Dynamic case

Without delay effect, damaged zone is restricted to the elements directly linked with loads [10], whereas using delay effect allows the damage response not to depend on the mesh. Thus the theoretical damaged zone in this test is the fitted zone. It improves the use of delay effect.

<figure><img src="image_6.png" /><figcaption>Figure 5. Geometry of tensile/compressive case.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Figure 6. Quasi-static loading.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>Figure 7. Damage evolution in function of time for loading 1: curve 1 (respectively 2, 3) represents damage 1 (respectively 2 and 3).</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Figure 8. Damage evolution in function of time for loading 2: curve 1 (respectively 2, 3) represents damage 1 (respectively 2 and 3).</figcaption></figure>

## 4.2 Blade impact

This experiment can be compared to a bird strike on fan blade. The blade is blocked on its down-basis, and the dynamic loading is the same than in figure 9 on a small area of the blade. The dimensions are the following: 600 mm × 150 mm × 30 mm. The result is the opposed edge to the load is actually broken in real experiment (bird strike at SNECMA in 2003).

### 4.2.1 Linear strain damage model

The model used here is the linear strain damage model with delay effect. The initial speed of the small area is 50m/s. The figure 10 presents the results of simulation: the damage $d_{13}$ on deformed mesh at three
<figure><img src="image_10.png" /><figcaption>**Figure 9.** dynamic loading.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>**Figure 10.** Damage 13 on deformed mesh for three times (15, 55 and 80µs) after impact.</figcaption></figure>

different times (15, 55 and 80µs) after impact. We can observe that not only the loaded zone presents a significant damage, but the up left edge (opposed to the load) too. The real experiment underlines the same behaviour: after impact this opposed edge is broken.

### 4.2.2 ODM model

The model used is the Onera Damage Mechanism one. Table 1 presents the values for the different parameters of the model.

Figure 11 shows the damage 3 on deformed mesh for four times (15, 55, 80 and 110µs). Figure 12 shows the displacement of the edge opposed to the load as a function of time. It is underlines that the edge displacement in direction 1 is relatively important.

Similarly this second simulation reveals that damaged zones are the loaded one and the opposed edge too.

<table><caption>**Table 1.** Parameters of ODM model.</caption><thead><tr><th>Parameters</th><th>y<sub>o</sub><sup>n</sup>(i), y<sub>o</sub><sup>t</sup>(i)</th><th>y<sub>c</sub><sup>n</sup>(i), y<sub>c</sub><sup>t</sup>(i)</th><th>p<sup>n</sup>(i), p<sup>t</sup>(i)</th><th>d<sub>c</sub><sup>n</sup>(i), d<sub>c</sub><sup>t</sup>(i)</th><th>ξ<sub>i</sub></th><th>a<sub>(i)</sub><sup>f</sup></th><th>Δε<sub>f</sub>(i)</th><th>b<sub>i</sub></th><th>h<sup>p</sup>(i), h<sup>n</sup>(i), h<sup>pn</sup>(i)</th></tr></thead><tbody><tr><td>Values</td><td>5,5 Pa</td><td>3410 Pa</td><td>0,90</td><td>2,85</td><td>0,15</td><td>0.</td><td>5,51e<sup>-3</sup></td><td>1.</td><td>1.</td></tr></tbody></table>

<figure><img src="image_12.png" /><figcaption>**Figure 11.** Damage 3 on deformed mesh for four times (15, 55, 80 and 110µs) after impact.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>**Figure 12.** Displacement of the opposed edge to the load as a function of time.</figcaption></figure>
## 5. CONCLUSIONS

In brief, modelling composite materials makes appear damage-induced anisotropy, such as microcracks opening, closure, shearing, residual strains. These phenomena are taken into account in the Onera Damage Mechanism model, which describes elastic behaviour coupled with both damage desactivation effects, unilateral effects, which could have been verified. The numerical implementation of the ODM model in the fast dynamics code EUROPLEXUS has been achieved for treating dynamic damage response. As seen above, to estimate the damage in a composite structure after impact, not only damage description but mesh dependency have to be solved. Afterwards simulations can estimate the local damaged zones with good accuracy. As this model can be applied to wide variety of materials, it may be used for composite fan blade. Nevertheless the identification of the several parameters is as difficult as rich is the model.

## Acknowledgements

The authors acknowledge Prof. Combescure for useful discussions about numerical simulations. They wish to thank ONERA, and more particularly J.F. Maire and N. Carrere, for their authorization to use their Damage Model. They are grateful to SNECMA (SAFRAN Group) for its financial support and cooperation.

## References

[1] Maire J.F. and Chaboche J.L., A New Formulation of Continuum Damage Mechanics (ODM) for Composite Materials, Aerospace Science and Technology 4 (1997) 247-257.

[2] Allix O. and P. Ladevèze, Damage analysis of interlaminar fracture specimens, Composite Structures 31 (1995) 61-74.

[3] Maire J.F. and Carrere N., Onera Damage Mechanics.

[4] Chaboche J.L., Continuum Damage Mechanics: Part I and II, ASME J. Of Applied Mechanics 55 (1988) 59-72.

[5] Maire J.-F. and Lesne P.-M., An Explicit Damage Model for the Design of Composites Structures, Composites Science and Technology 58 (1998) 773-778.

[6] Maire J.-F. and Carrere N., Modèles ODM, ONERA (2005). in French.

[7] Chaboche J.L. and Maire J.-F., New progress in micromechanics-based CDM models and their application to CMCs, Composites Science and Technology 61 (2001) 2239-2246.

[8] Quilici S., Industrialisation des modèles ODM de l'ONERA, Lot n°2 : Formulation en déformations, (2003). in French.

[9] Allix O. and Deü J.F., Delay-damage modelling for fracture prediction of laminated composites under dynamic loading; Eng. Trans. 45 (1997) 29-46.

[10] Suffis A., Lubrecht T. and Combescure A., Damage model with delay effect: Analytical and numerical studies of the evolution of the characteristic damage length, International Journal of Solids and Structures 40 (2003) 3463-3476.