# Modélisation numérique de la propagation instable et limitée de fissure due à la présence de zones locales fragiles

Haïthem Adouani, Bing Tie, Clotilde Berdin, Denis Aubry

► To cite this version:

Haïthem Adouani, Bing Tie, Clotilde Berdin, Denis Aubry. Modélisation numérique de la propagation instable et limitée de fissure due à la présence de zones locales fragiles. 8e Colloque national en calcul des structures, CSMA, May 2007, Giens, France. hal-01495832

HAL Id: hal-01495832

https://hal.science/hal-01495832v1

Submitted on 26 Mar 2017

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 CC0 1.0 - Universal - International License
# Modélisation numérique de la propagation instable et limitée de fissure due à la présence de zones locales fragiles

Haithem Adouani*, Bing Tie*, Clotilde Berdin*, Denis Aubry*

*Laboratoire de Mécanique des Sols,

Structures et Matériaux (CNRS UMR 8579)

École Centrale Paris, Grande Voie des Vignes

92295 Châtenay-Malabry Cedex, France

haithem.adouani@ecp.fr, bing.tie@ecp.fr, clotilde.berdin@ecp.fr, denis.aubry@ecp.fr

**RÉSUMÉ.** Assurer l'intégrité des structures et composants industriels passe par la considération de la présence de défauts et par la détermination des conditions d'amorçage, de propagation et d'arrêt des fissures dues à la présence de ces défauts. Notre étude est motivée par la prédiction de la propagation dynamique et de l'arrêt de fissures dans les composants industriels de grande taille. Il s'agit d'étudier l'extension instable de ces fissures lorsqu'elles s'amorcent dans des zones locales fragiles et de prédire leur arrêt dû en particulier à des zones plus tenaces. La propagation instable et limitée de fissure est associée au phénomène dit de « pop-in » [4]. Par ailleurs, ces fissures se propageant à des vitesses de l'ordre de 10% de la vitesse des ondes de pression dans le matériau, il est donc nécessaire de déterminer le rôle de la dynamique du phénomène dans la propagation et l'arrêt de ces fissures.

**ABSTRACT.** The assessment of some industrial security components often needs to predict the initiation of brittle crack. A second level of security can be achieved in predicting crack arrest, especially when crack initiates in local brittle zone and extends towards tougher material zones. In this case, one observed the so-called “pop-in” phenomenon [4], which comes from limited brittle crack propagation. Since brittle crack propagates at speed about 10% of dilatational wave speed, elastic waves can be triggered, so it is necessary to investigate if they may influence the crack advance and arrest.

**MOTS-CLES.** Propagation instable de fissure, pop-in, solveur non linéaire, zone cohesive.

**KEYWORDS.** Instable crack propagation, pop-in, non linear solver, cohesive zone.
# 1. Introduction

Au cours de cette étude, nous proposons d'analyser le phénomène du pop-in induit par la propagation instable de fissure due à la présence d'une zone locale fragile en développant une approche numérique appropriée, qui permet une modélisation fiable de la propagation limitée de fissure. La fissuration est modélisée par la théorie de la zone cohésive. Ce modèle combine à la fois un critère énergétique et un critère en contrainte locale afin de quantifier les mécanismes de dissipation d'énergie dans la zone du processus de fissuration. Le phénomène de dissipation est traité par l'introduction d'une énergie de surface définie au niveau de la surface cohésive, ce qui constitue la principale différence par rapport aux approches locale et globale de la rupture. Par conséquent, il est nécessaire de créer explicitement les surfaces cohésives de fissuration dans les modèles de discrétisation numérique. Dans notre étude, nous avons choisi d'utiliser des éléments finis d'interface pour représenter les fissures. Ces éléments ont été développés dans un code éléments finis orienté objets OOFE (Object Oriented Finite Element program) pour des structures 2D et 3D. Le comportement mécanique de ces éléments d'interface est régi par la loi de comportement cohésive adoucissante développée pour évaluer l'endommagement et la rupture des matériaux étudiés [2, 3]. Le solveur dynamique développé pour l'intégration temporelle de la propagation de fissure est discuté. Comme la géométrie et par conséquent le maillage de la structure étudiée évoluent au cours du temps d'une manière continue et rapide durant l'analyse transitoire, il est primordial d'assurer la stabilité numérique de l'intégration en temps. L'étude théorique de la stabilité de notre solveur dynamique est présentée.

# 2. Modèles et solveurs numériques

Ce chapitre présente les modèles adoptés et les outils numériques développés dans cette étude.

## 2.1. Notations et problème mécanique

Nous considérons un domaine $\Omega$ borné par sa frontière $\partial\Omega$, qui est la réunion des trois parties $\Gamma_u$, $\Gamma_\sigma$ et $\Gamma_c$. $\Gamma_\sigma$ est la partie où les forces surfaciques externes $g_\sigma$ sont appliquées, $\Gamma_u$ est celle encastrée et $\Gamma_c$ est une fissure contenue dans le volume $\Omega$. $\Omega$ est également soumis à des chargements volumiques externes $f$. La fissure $\Gamma_c$ contient une zone de processus de fissuration, notée $\Gamma_{cohe}^{\pm}$, qui évolue au cours du temps et dont le comportement est régi par une loi cohésive. Selon la théorie de la zone cohésive, les tractions cohésives $T^+$et $T^-$ sont appliquées respectivement sur $\Gamma_{cohe}^+$ et $\Gamma_{cohe}^-$ (Figure 1-a).
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Domaine $\Omega$ fissuré par $\Gamma_c$ incluant $\Gamma_{cohe}^{\pm} = \Gamma_{cohe}^{+} \cup \Gamma_{cohe}^{-}$.</figcaption></figure>

La théorie cohésive définit une loi de comportement entre les forces cohésives $\mathbf{T} = \mathbf{T}^{-} = -\mathbf{T}^{+}$ et les déplacements d'ouverture des deux lèvres des surfaces cohésives, selon le modèle suivant :

$$ \mathbf{T}([\mathbf{u}]_x) = \mathbf{H}([\mathbf{u}]_x) \cdot [\mathbf{u}]_x = (1 - D([\mathbf{u}]_x)) \cdot \mathbf{K}^{+} \cdot [\mathbf{u}]_x \quad [1] $$

Où $[\mathbf{u}]_x = \mathbf{u}(x^+) - \mathbf{u}(x^-)$ est le saut de déplacement des surfaces cohésives, $\Gamma_{cohe}^{\pm}$, $\mathbf{H}([\mathbf{u}]_x)$ est le tenseur non linéaire de la loi cohésive qui dépend d'une variable d'endommagement D selon une loi de Crisfield [2] (Figure 1-b) et $\mathbf{K}^{+}$ est la raideur élastique initiale de la zone cohésive.

Les équations d'équilibre dynamique du domaine $\Omega$ s'écrivent comme suit :

$$
\begin{align*}
\mathbf{Div}\sigma(\mathbf{u}(x)) + \mathbf{f} &= \rho \ddot{\mathbf{u}} && \text{sur } \Omega \\
\sigma \cdot \mathbf{n} &= g_{\sigma} && \text{sur } \Gamma_{\sigma} \\
\sigma(\mathbf{u}^{\pm}) \cdot \mathbf{n}^{\pm} &= \mathbf{T}^{\pm} && \text{sur } \Gamma_{\text{cohe}}^{\pm} \quad [2] \\
\mathbf{T}^{+} + \mathbf{T}^{-} &= \mathbf{0} && \text{sur } \Gamma_{\text{cohe}}^{\pm} \\
\mathbf{u} &= \mathbf{0} && \text{sur } \Gamma_{\mathbf{u}} \\
[\mathbf{u}]_x \cdot \mathbf{n}^{-} &\ge 0 && \text{sur } \Gamma_{\text{cohe}}^{\pm}
\end{align*}
$$

Où $\sigma$ est le champ de contraintes, $\mathbf{u}$ est le champ de déplacement dans le domaine $\Omega$ et $\mathbf{n}$ et $\mathbf{n}^{\pm}$ sont respectivement les vecteurs normaux sortants aux surfaces $\Gamma_{\sigma}$ et $\Gamma_{\text{cohe}}^{\pm}$. Les surfaces cohésives sont maillées par le biais d'éléments finis d'interface. Ainsi, des éléments finis cohésifs sont développés dans le code OOOFE. Des calculs réalisés sur des éprouvettes CT 25 en 2D et en 3D ont validé le développement de ces éléments cohésifs [1].

## 2.2. Solveurs numériques

Des solveurs non linéaires du type Newton modifié sont implémentés. La discrétisation en temps se fait de manière incrémentale. Le calcul des variables du problème à l'instant $t^{k+1}$ à partir des variables déterminées à l'instant $t^k$ se fait de manière itérative. Nous avons développé deux solveurs. Pour le premier, la matrice de rigidité initiale est utilisée comme matrice d'itération.
Pour améliorer la convergence de la résolution non linéaire, le deuxième solveur, plus coûteux en temps, réactualise au début de chaque étape de temps la matrice de rigidité qui dépend de la rigidité $(1 - D^k)K^+$ des éléments joints à l'étape $t^k$. Des difficultés de convergence sont mises en évidence lorsque le pop-in a lieu. Pour cela, une stratégie de discrétisation adaptative en temps est adoptée et s'avère indispensable. Par cette stratégie, le pas de temps $\Delta t$ peut être diminué ou augmenté d'une étape de temps à une autre selon l'état de convergence de la résolution itérative.

## 2.3. Stabilité dynamique

Les résultats des calculs présentés dans la section suivante portent sur des problèmes quasi-statiques. Néanmoins, un solveur non-linéaire dynamique est en cours de développement et sa propriété de stabilité a été étudiée analytiquement. Ce solveur est basé sur la méthode de Galerkin espace-temps discontinue en temps [4, 5]. Cette méthode consiste à considérer un domaine espace-temps $S = \Omega \times [0, T[$, où $[0, T[$ représente l'intervalle du temps du problème dynamique. Le domaine espace-temps est subdivisé en plusieurs tranches espace-temps $S_n = \Omega \times [t_n, t_{n+1}[$. Les inconnues du problème, qui sont les champs de déplacement **u** et de vitesse **v**, sont discontinues en temps entre deux tranches espace-temps successives (Figure 2).

<figure><img src="image_3.png" /><figcaption>**Figure 2.** Illustration d'un domaine espace-temps</figcaption></figure>

Les sauts en temps des différents champs mécaniques sont donnés par:

$$ \mathbf{u}_n^{\pm} = \lim_{\delta \to 0} \mathbf{u}(t_n \pm \delta); [\mathbf{u}_n]_{t_n} = \mathbf{u}_n^{+} - \mathbf{u}_n^{-} \quad [3] $$

La méthode de Galerkin discontinue en temps, dans le cas d'un problème élasto-dynamique, garantit la stabilité inconditionnelle du solveur élasto-dynamique même si le maillage change au cours du temps [5]. Dans le cas où la structure contient une fissure qui se propage au cours du temps, ce qui est une situation nouvelle, nous montrons que cette stabilité inconditionnelle est conservée. En effet, l'analyse de stabilité est basée sur une approche énergétique [5]. Ainsi, en notant $E(t)$ l'énergie totale à l'instant $t$, nous montrons que :

$$ E(t_{n+1}^{-}) - E(t_n^{-}) = -E_J(t_n^{-}) - \int_{t_n^{-}}^{t_{n+1}^{-}} (T([\mathbf{u}]_x), [\mathbf{v}]_x)_{\Gamma_{\text{cohe}}} \le 0 \quad [4] $$

Où $E_J$ est l'énergie dissipée dans les discontinuités en **u** et **v**, qui est propre à la méthode de Galerkin discontinue en temps. La relation ci-dessus montre clairement que l'énergie cohésive dissipée renforce la stabilité inconditionnelle établie par la méthode de Galerkin discontinue en temps.
### 3. Modélisation du pop-in

Ce chapitre présente des résultats numériques obtenus sur une éprouvette CT25 présentant des zones locales fragiles en pointe de fissure (Figure 3). Les deux matériaux, la zone locale fragile (la zone rouge sur la figure 3) et le métal de base, ont un comportement élasto-plastique. La zone locale fragile présente une contrainte d'écoulement plus élevée que le matériau de base. La contrainte critique cohésive est plus basse dans la zone locale fragile alors que les énergies dissipatives critiques sont identiques. L'éprouvette est chargée en déplacement imposé.

<figure><img src="image_4.png" /><figcaption>**Figure 3.** Rupture d'une éprouvette CT25 présentant une zone locale fragile</figcaption></figure>

#### 3.1. Choix du solveur et paramètres numériques

Avec l'arrivée du pop-in, qui correspond au début de la chute de force sur la courbe force-déplacement (Figure 4-a), la convergence de la résolution non linéaire est très difficile. La stratégie de pas de temps adaptatif s'est avérée indispensable pour modéliser le pop-in. En effet, pendant le pop-in, il est nécessaire d'utiliser un pas de temps très petit, qui est de l'ordre $\Delta t = 10^{-11} U_{\text{total}}$. Pour le solveur utilisant la matrice de rigidité initiale comme matrice d'itération, une autre difficulté de la résolution non linéaire a lieu à la fin du pop-in où la force recommence à augmenter. À ce moment, si on garde le petit pas de temps pour le pop-in, on obtient des oscillations (Figure 4-a). La première stratégie adoptée pour éviter ces oscillations consiste à autoriser l'augmentation du pas de temps après le pop-in (Figure 4-b).

Une deuxième alternative consiste à utiliser le solveur où la matrice de rigidité est réactualisée en début de chaque étape, qui réduit notablement les oscillations numériques (Figure 4-c).

<figure><img src="image_5.png" /><figcaption>**Figure 4.** Courbes force-déplacement montrant la resolution des oscillations numériques</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Figure 5.** Courbe force-déplacement a) et avancée de la fissure b) en fonction de l'ouverture totale pour deux lois d'écrouissage différentes du métal de base</figcaption></figure>

**3.2. Etudes paramétriques du pop-in**

L'influence des différents paramètres mécaniques sur le phénomène du pop-in a été étudiée. Ici nous montrons comme exemple, l'influence de la plasticité du métal de base sur le pop-in. Une loi d'écrouissage isotrope de type $\sigma_{eq} = \sigma_{y0} + C(\varepsilon^p)^n$ est utilisée pour le métal de base avec différentes valeurs de C. Les résultats montrent que plus C est petit, plus l'amorçage du pop-in est retardé, ce qui donne une longueur de fissure légèrement plus grande (Figure 5).

**4. Conclusions**

Pour modéliser le pop-in, la théorie de la zone cohésive a été adoptée. Des difficultés de convergence de la résolution non linéaire dans le cas du pop-in sont mises en évidence et des stratégies sont proposées et testées pour surmonter ces difficultés. Par exemple, la stratégie de pas de temps adaptatif se révèle indispensable pour modéliser le pop-in avec un pilotage en déplacement imposé. Les outils numériques développés dans le code de calcul OOFE ont été validés sur une éprouvette CT25 ayant une zone locale fragile en pointe de fissure. Des études sur l'influence de paramètres mécaniques sur le pop-in ont été réalisées. Par ailleurs, la modélisation de la propagation dynamique de la fissure par la méthode de Galerkin discontinue en temps montre analytiquement que le schéma d'intégration en temps ainsi défini est inconditionnellement stable dans le cadre d'un modèle de zone cohésive.

**References**

[1] Adouani H. et al., « Adaptive numerical modeling of dynamic crack propagation », *Conférence Internationale DYMAT 2006*.

[2] de Borst R., "Numerical aspects of cohesive zone models", *Engineering Fracture Mechanics*, **70** (2003), 1743-1757.

[3] Crisfield M. A. et al. "Some aspects of the non-linear finite element method". *Finite Elements in Analysis and Design*, **27** (1997), 19-40.

[4] Dawes M. G. "Quantifying pop-in severity in fracture toughness tests". *Fatigue Fracture Engineering Materials and Structures*, **14** (1991), 1007-1014.

[5] Tie B. et al., « Adaptive computation for elastic wave propagation in plate/shell structures under moving loads ». *European Revue of FEM*, **12** (2003), 717-736.