r7.02.21 Loi de comportement cohésive#

r7.02.21 dépendante de la mixité du chargement : CZM_TURON#

Résumé:

La loi de comportement cohésive CZM_TURON permet de modéliser l’endommagement d’une interface fissurante telle que l’ouverture de la fissure génère des forces de cohésion entre les lèvres de celle-ci.

Elle relie les forces cohésives, homogènes à une contrainte, à la discontinuité de déplacement à travers l’interfaceselon une loi bi-linéaire, qui dépend du taux de mixité du chargement: une interface soumise à un chargement mixte s’endommagera plus tôt qu’en mode pur.

Ce modèle a été développé par l’équipe de Turon de l’Université de Girona (Espagne) [ 1 ]. Elle est par exemple utilisée dans le cadre de la modélisation de joints collés en matériau composite dans les pâles d’éoliennes.

La loi est supportée par des éléments finis de joint dont les degrés de liberté sont les sauts de déplacement au passage de l’interface (voir R3.06.13).

Nous présentons ici la formulation de la loi cohésive ainsi que les variables internes du modèle. L’existence d’instabilités lors de la propagation de fissures nécessite un pilotage par la prédiction élastique dont on présente les éléments spécifiques à ces lois.

Définition de la loi CZM_TURON#

Définitions#

Commençons par introduire les paramètres du modèle. On a:

  • les densités d’énergie de surface critique en mode normal pur et en mode tangentiel pur, notées \({G}_{c,N}\) (mot clé GC_N) et \({G}_{c,T}\) (mot-clé GC_T ) qui traduisent le coût énergétique de la fissuration dans chacun des deux modes;

  • les contraintes critiques à la rupture notées \({\sigma}_{N}^{c}\) (mot-clé SIGM_C_N) et \({\sigma}_{T}^{c}\) (mot-clé SIGM_C_T);

  • la raideur initiale de l’interface \(K\) (mot-clé K);

  • le paramètre

    ../../../../_images/Shape218.gif

du critère de Benzeggagh-Kenane (mot-clé ETA_BK);

  • le paramètre adimensionnel C_RUPT qui permet de définir une rigidité résiduelle «post-rupture» (lorsque l’élément est complètement endommagé et que la contrainte cohésive est nulle), comme une fraction de la raideur initiale. Il vaut 0,001 par défaut et ne peut excéder 0,1;

  • le paramètre CRIT_INIT qui permet de choisir le type de critère d’initiation de l’endommagement, parmi “YE” (critère elliptique) ou “TURON” (critère de type Benzeggagh-Kenane).

Les paramètres matériau à proprement parler, à savoir GC_N, GC_T, SIGM_C_N, SIGM_C_T, K et ETA_BK sont à identifier à partir d’essais quasi-statiques standards. Par exemple, dans le cadre de la simulation d’un joint collé, l’adhésif étant modélisé avec le modèle de zone cohésive, on identifiera les paramètres en réalisant les essais suivants sur des éprouvettes «sandwich»:

  • l’essai DCB («Double Cantilever Beam») qui sollicite le joint en mode normal pur permet d’identifier les paramètres GC_Net SIGM_C_N;

  • l’essai ENF («End-notched flexure»), correspondant à une flexion 3-points, sollicite le joint en mode tangent pur et permet d’identifier les paramètres GC_Tet SIGM_C_T;

  • l’essai MMB («Mixed-mode bending») est une combinaison des deux essais précédents et permet de soumettre le joint à un chargement mixte, le taux de mixité étant à la main de l’expérimentateur qui peut ajuster le dispositif. En complément des essais précédents, il est utilisé pour caler le paramètre ETA_BK.

En ce qui concerne le paramètre \(K\) , on constate généralement qu’il influence peu la réponse globale simulée et qu’il est donc difficile de l’identifier par ajustement des courbes numérique et expérimentale. En pratique, il est souvent supposé ad hoc .

Remarque :

On adopte les notations suivantes: pour tout vecteur \(x\) on note:

\(x={x}_{N}N+{x}_{{T}_{1}}{T}_{1}+{x}_{{T}_{2}}{T}_{2}\)

\({x}_{T}={x}_{{T}_{1}}{T}_{1}+{x}_{{T}_{2}}{T}_{2}\)

\(⟨x⟩=⟨{x}_{N}⟩N+{x}_{{T}_{1}}{T}_{1}+{x}_{{T}_{2}}{T}_{2}\) avec \(⟨{x}_{i}⟩=max({x}_{i},0)\)

\(\Vert x\Vert =\sqrt{x.x}\) et \({\Vert x\Vert }_{+}=\sqrt{⟨x⟩.⟨x⟩}\)

En particulier, on notera \(\lambda ={\Vert \Delta \Vert }_{+}\) le déplacement équivalent de l’élément cohésif et \({\lambda}_{T}=\Vert {\Delta}_{T}\Vert\) la norme du saut tangentiel.

On appelle enfin:

  • \({\Delta}^{0}\) la valeur que doit attendre le saut équivalent pour que l’endommagement s’initie. On l’appelle donc saut d’initiation de l’endommagement;

  • \({\Delta}^{f}\) la valeur que doit atteindre le saut équivalent pour que l’élément casse. On l’appelle donc saut de propagation de la fissure (mais cette appellation ne doit pas tromper le lecteur: les modélisations avec loi CZM ne présentent pas de fissure à proprement parler, par exemple directement maillée ou bien définie par des level sets).

La raideur initiale \(K\) permet de régulariser l’énergie libre en 0. La présente loi est donc une loi de type régularisée. Notons que dans le cas où l’interface est censée être infiniment rigide, par exemple lorsque l’interface modélise un trajet de fissuration dans un matériau fragile, alors la raideur initiale doit être suffisamment petite pour que l’énergie soit bien régularisée, mais cela peut engendrer en retour un assouplissement global de la structure qui n’est pas physique.

Dans le cadre de la modélisation d’une interface possédant une certaine «souplesse»qui est physique (par exemple une couche de colle entre deux adhérents qui va se déformer élastiquement avant de commencer à s’endommager), ce problème ne se pose pas.

Il convient de rappeler qu’il ne doit pas y avoir une trop grande différence de raideur entre les éléments cohésifs de type JOINTet les éléments directement voisins (généralement volumique), au risque de rencontrer des problèmes de conditionnement de la matrice de rigidité du système, propre à ce type d’éléments [ 5 ]. En pratique, on peut vérifier que le paramètre \(K\) et la quantité \(\frac{E}{{h}_{\mathit{elem}}}\) sont du même ordre de grandeur, \(E\) étant le module d’Young de l’adhérent et \({h}_{\mathit{elem}}\) étant la taille des mailles adjacentes aux éléments joint. On peut, dans une certain mesure, adapter la taille des mailles adjacentes. Mais si les problèmes de conditionnement ne peuvent pas être levés compte-tenu d’une trop grande différence entre les propriétés matériaude l’interface cohésive et de l’adhérent, on s’orientera plutôt vers des éléments de type INTERFACE[ 5 , 6 ], sur lesquels sont disponibles certaines lois CZM [ 7 ].

Équations de comportement#

Dans cette partie on présente l’énergie de surface :math:`Psi ` ainsi que le vecteur contrainte qui

../../../../_images/Shape316.gif

en dérive. La formulation proposée assure l’irréversibilité de l’endommagement.

Énergie libre de surface#

L’énergie libre par unité de surface \(\Psi\) , définie sur une discontinuité \(\Gamma ` , dépend du vecteur saut de déplacement entre les lèvres de la fissure :math:\)Delta ` (noté :math:`Delta ` pour alléger l’écriture) et d’une variable interne qui gère l’irréversibilité de la fissuration.

On définit tout d’abord la variable \(v\) quimémorise la plus grande norme du saut (déplacement équivalent)atteinte au cours de l’ouverture. Sa loi d’évolution entre deux incréments de chargement successifs - et + s’écrit:

\({v}^{+}=max({\lambda}_{+},{v}^{-})\)

On définit également la variable seuil notée \(r\) qui donne le saut minimal à atteindre pour (re)commencer à endommager:

\({r}^{+}=max({\Delta}^{0,}{v}^{-})\)

Définissons enfin la variable d’endommagemet \(d\) dont on donnera la loi d’évolution dans la suite. L’élément est sain lorsque \(d=0\) et complètement endommagé, cassé, lorsque \(d=1\) .

On note que du point de vue thermodynamique, on utilise ici la variable \(d\) est la variable interne du comportement. En revanche, du point de vue de l’implémentation numérique, on utilisera la variable \(v\) comme variable interne, les autres variables \(r\) et \(d\) étant calculées à partir de \(v\) .

On se donne l’énergie libre par unité de surface suivante:

(4561)#\[\Psi (\Delta ,d)=(1-d){\Psi}^{0}(\Delta )\]

avec \({\Psi}^{0}(\Delta )=\frac{1}{2}{\Delta}_{i}{D}_{ij}^{0}{\Delta}_{j}\) pour \(i,j=1,3\)

\({\Psi}^{0}\) étant une fonction convexe dans l’espace des sauts de déplacement.

Les valeurs négatives de \({\Delta}_{N}\) n’ont pas de sens car l’interpénétration des lèvres est rendue impossible par le contact. Pour empêcher l’interpénétration après dé-cohésion complète, on propose la modification de l’équation selon la forme suivante:

\(\Psi (\Delta ,d)=(1-d){\Psi}^{0}(\Delta )-d{\Psi}^{0}({\delta}_{Ni}⟨-{\Delta}_{N}⟩)\)

où \({\delta}_{ij}\) est le symbole de Kronecker.

Expression de la contrainte#

Les équations de comportement sont obtenues en différenciant l’énergie libre par rapport aux sauts de déplacement:

(4562)#\[{\tau}_{i}=\frac{\partial \Psi }{\partial {\Delta}_{i}}=(1-d){D}_{ij}^{0}{\Delta}_{j}-d{D}_{ij}^{0}{\delta}_{\mathit{Nj}}⟨-{\Delta}_{N}⟩\]

où \({D}_{ij}^{0}={\delta}_{ij}K\) est le tenseur de raideur initial, avec \(K\) la raideur de pénalisation.

En supposant que \(N\equiv 1\) et que \(({T}_{1,}{T}_{2})\equiv (2,3)\) , l’équation peut s’écrire de manière vectorielle:

\(\tau =\left\lbrace \begin{array}{c}{\tau}_{1}\\ {\tau}_{2}\\ {\tau}_{3}\end{array}\right\rbrace =(1-d)K\left\lbrace \begin{array}{c}{\Delta}_{1}\\ {\Delta}_{2}\\ {\Delta}_{3}\end{array}\right\rbrace -dK\left\lbrace \begin{array}{c}⟨-{\Delta}_{1}⟩\\ 0\\ 0\end{array}\right\rbrace\)

L’énergie dissipée pendant l’évolution de l’endommagement est notée :math:`Xi ` et s’obtient en différenciant l’énergie libre par la variable interne d’endommagement:

\(\Xi =-\frac{\partial \Psi }{\partial d}\dot{d}\)

Le modèle de comportement donné à l’équation est bien évalué si la variable d’endommagement \(d\) peut être calculée en chaque instant du processus de déformation, ce qui est possible en se donnant un critère d’initiation et une loi d’évolution pour l’endommagement.

Évolution de l’endommagement#

L’évolution de l’endommagement est formulé sous la forme d’une fonction seuil dans l’espace des sauts de déplacement. Elle s’écrit, en chaque instant \(t\) :

\(F({\lambda}^{t},{r}^{t}):={\lambda}^{t}-{r}^{t}\le 0\) \(\forall t\ge 0\)

où \({r}^{t}\) est le seuil d’endommagement à l’instant actuel.

Si on note \({r}^{0}\) est le seuil d’endommagement initial alors on doit avoir \({r}^{t}\ge {r}^{0}\) en chaque instant. Le seuil d’endommagement ne peut qu’augmenter. L’endommagement s’initie lorsque la norme du saut de déplacement \(\lambda ` dépasse le seuil initial :math:`{r}^{0}\) , qui se calcule notamment à partir des caractéristiques matériau.

Dans la présente formulation, on utilisera une expression équivalente plus facile à traiter sur le plan algorithmique:

\(\overline{F}({\lambda}^{t},{r}^{t}):=G({\lambda}^{t})-G({r}^{t})\le 0\) \(\forall t\ge 0\)

où \(G\) est une fonction scalaire monotone allant de 0 à 1. L’évolution du seuil et de l’endommagement sont donnés par leurs accroissements :

\(\dot{r}=\dot{\mu}\)

\(\dot{d}=\dot{\mu}\frac{\partial \overline{F}(\lambda ,r)}{\partial \lambda }=\dot{\mu}\frac{\partial G(\lambda )}{\partial \lambda }\)

la variable :math:`mu ` étant le multiplicateur de Lagrange qui permet de définir les conditions de charge – décharge suivant les relations de Kuhn et Tucker:

\(\dot{\mu}\ge 0\) (a); \(\overline{F}({\lambda}^{t},{r}^{t})\le 0\) (b); \(\dot{\mu}\overline{F}({\lambda}^{t},{r}^{t})=0\) (c)

Par soucis pédagogique, on fera explicitement l’intégration – très classique – des relations de Kuhn et Tucker.Pour cela, on sépare les deux «phases» :

Phase A (non-endommageante) :

Tant que \({\lambda}^{t}<{r}^{t}\) , on a \(G({\lambda}^{t})<G({r}^{t})\) par croissance de la fonction \(G\) , et donc \(\overline{F}({\lambda}^{t},{r}^{t})<0\) (cohérent avec (b)).

Par (c), on a nécessairement que \(\dot{\mu}=\dot{r}=0\) . Donc \(r\) reste constant, égal à sa valeur en début de phase, notamment à \({r}^{0}\) s’il s’agit de l’instant initial.

Phase B (endommageante) :

Si on regarde l’autre expression de (b), où \(\overline{F}({\lambda}^{t},{r}^{t})=0\) , alors \({\lambda}^{t}={r}^{t}\) (on ne peut avoir \({\lambda}^{t}>{r}^{t}\) sinon (b) n’est plus satisfait). Évacuons le cas où \(\dot{\mu}=0\) correspondant au cas où \(\lambda ` n'évolue pas. Alors :math:\)dot{mu}>0` d’après (a) et donc \(r\) croît, en même temps que \(\lambda ` . Par ailleurs, on peut écrire que :math:\)dot{d}=dot{r}frac{partial G(r)}{partial r}` d’où l’intégration directe de \(d\) .

Les deux phases peuvent s’alterner suivant le trajet de chargement, et on peut donc intégrer la variable seuil \(r\) et la variable d’endommagement \(d\) de la manière suivante:

\({r}^{t}=max({r}^{0},{max}_{s}\lambda )\) \(0\le s\le t\)

\({d}^{t}=G({r}^{t})\)

On notera dans la suite \({r}^{0}={\Delta}^{0}\) le seuil d’endommagement.

La fonction scalaire \(G\) définit l’évolution de la variable d’endommagement. Pour un taux de mixité donné, on propose une fonction de la forme suivante:

\(G(\lambda )=\frac{{\Delta}^{f}(\lambda -{\Delta}^{0})}{\lambda ({\Delta}^{f}-{\Delta}^{0})}\)

qui correspond à une équation de comportemententre la contrainte et le saut équivalents de forme bilinéaire, avec \({\Delta}^{0}\) est le saut de déplacement à l’initiation de l’endommagement, donc \({\Delta}^{0}={r}^{0}\) , et \({\Delta}^{f}\) est le saut de déplacement final lorsque l’endommagement atteint sa valeur maximale de 1.

On peut facilement écrire l’équation des droites dans les deux régimes en partant de l’équation écrite avec la contrainte et le saut équivalent :math:`tau =K(1-d)lambda =K(1-G(r))lambda ` :

  • en régime élastique (phase non-endommageante), avec \(\lambda <r\) , on a \(\tau =K\frac{{\Delta}^{0}({\Delta}^{f}-r)}{r({\Delta}^{f}-{\Delta}^{0})}\cdot \lambda ` , qui correspond bien à l’équation d’une droite puisque :math:`r\) n’évolue pas. En particulier, \(\tau =K\cdot \lambda ` lorsque :math:`r={\Delta}^{0}\) ;

  • en régime dissipatif (phase endommageante), avec \(r=\lambda ` , on a :math:\)tau =Kfrac{{Delta}^{0}{Delta}^{f}}{{Delta}^{f}-{Delta}^{0}}-Kfrac{{Delta}^{0}}{{Delta}^{f}-{Delta}^{0}}cdot lambda ` qui est bien l’équation d’une droite affine décroissante.

La loi bilinéaire de ce type peut être visualisée sur la Figure .

../../../../_images/1000020100000260000001C3A9563FFCBAA62FFC.png

Figure 2.2.1: Relation bi-linéaire entre les contrainte et saut de déplacement équivalents [ 1 ]

Les sauts de déplacement à l’initiation et à la propagation s’obtiennent à partir des surfaces seuil correspondantes, qu’il reste donc à définir pour compléter la formulation des équations de comportement.

Dans le présent modèle, ces critères prennent en compte l’interaction entre les modes de chargement normal et tangentiel. Ils dépendent donc du taux de mixité du chargement défini de la manière suivante:

\(\beta =\frac{{\Delta}_{T}}{⟨{\Delta}_{N}⟩+{\Delta}_{T}}\)

Critère de propagation de la fissure (surface seuil pour \(d=1\))#

Le critère de propagation de la «fissure» en mode mixte proposé est établi à partir du taux de restitution d’énergie critique. Pour les systèmes de type joint de composite, un critère courant est celui de Benzeggagh-Kenane, issu du domaine de la mécanique de la rupture:

\({G}_{c}={G}_{\mathit{Nc}}+({G}_{T,c}-{G}_{N,c}){\left(\frac{{G}_{T}}{{G}_{\mathit{tot}}}\right)}^{\eta}\)

où \({G}_{\mathit{tot}}={G}_{N}+{G}_{T}\) .

On précise que le critère, écrit dans le cadre de la mécanique de la rupture, fait intervenir les quantités \({G}_{I,c}\) , \({G}_{\mathit{II},c}\) et \({G}_{\mathit{III},c}\) qui sont les énergies critiques selon les trois modes de sollicitation classiques d’une fissure. Mais dans le cadre d’une interface modélisée par une loi d’interface CZM, on ne peut décrire que

\({G}_{\mathit{II},c}={G}_{\mathit{III},c}={G}_{T,c}\)

puisque la définition de ce qu’on appelle communément les modes II et III est faite à partir du front de fissure. Or ici, nous ne modélisons pas une fissure mais une interface endommageante, pour laquelle on peut au mieux définir une direction normale et une direction tangentielle.

On peut réécrire ce critère à partir des sauts de déplacement en notant que l’énergie critique, correspondant à l’aire sous la courbe bilinéaire de la Figure , s’écrit:

\({G}_{c}=\frac{1}{2}K{\Delta}^{0}{\Delta}^{f}\)

ce qui donne

\({\Delta}^{f}=\frac{{\Delta}_{N}^{0}{\Delta}_{N}^{f}+({\Delta}_{T}^{0}{\Delta}_{T}^{f}-{\Delta}_{N}^{0}{\Delta}_{N}^{f}){\left(\frac{{G}_{T}}{{G}_{\mathit{tot}}}\right)}^{\eta}}{{\Delta}^{0}}\)

avec \({\Delta}_{N}^{0}\) et \({\Delta}_{T}^{0}\) les critères d’initiation de l’endommagement, et \({\Delta}_{N}^{f}\) et \({\Delta}_{T}^{f}\) les critères de propagation en modes purs, calculés à partir des données matériau pour chaque mode compte-tenu de la forme bilinéaire de la loi (voir figure ):

\({\Delta}_{N/T}^{0}=\frac{{\sigma}_{c,N/T}}{K}\) et \({\Delta}_{N/T}^{f}=\frac{2{G}_{c,N/T}}{{\sigma}_{c,N/T}}\)

On notera que le taux de mixité peut être mesuré par une autre quantité \(B\) à partir des taux de restitutiond’énergie :

\(B=\frac{{G}_{T}}{{G}_{\mathit{tot}}}\)

Les deux définitions sont reliées par la formule suivante comme cela est démontré dans [ 1 ]:

\(B=\frac{{\beta}^{2}}{1+2{\beta}^{2}-2\beta }\)

Critère d’initiation de l’endommagement (surface seuil pour \(d=0\))#

Lorsque le chargement est en mode N ou T pur, l’endommagement est initié lorsque les forces cohésives atteignent leur valeur critique respective \({\sigma}_{c,N/T}\) (contraintes critiques). Sous un chargement mixte, il faut prendre en compte une interaction entre les modes. Il existe différents modèles définissant le critère l’initiation en fonction du taux de mixité. La présente loi de comportement propose le modèle de Benzeggagh–Kenane (qui est donc utilisé à la fois pour le critère d’initiation et de propagation – et ce avec la même puissance \(\eta ` ). On désigne ce critère par «critère de Turon», le fait d’adopter la même équation de couplage des modes pour l’initiation et la propagation étant un choix fait par les auteurs de [ :ref:`1 <RefNumPara__371857051>\) ]. Un autre critère purement elliptique, désigné par «critère de Ye» est également proposé.

  • Le critère elliptique de Ye

C’est un critère simple et plutôt répandu. Il est basé sur un interaction quadratique des contraintes critiques:

\({\left(\frac{⟨{\sigma}_{N}⟩}{{\sigma}_{c,N}}\right)}^{2}+{\left(\frac{{\sigma}_{{T}_{1}}}{{\sigma}_{c,T}}\right)}^{2}+{\left(\frac{{\sigma}_{{T}_{2}}}{{\sigma}_{c,T}}\right)}^{2}=1\)

ce qui s’écrit de manière équivalente en termes de saut de déplacement:

\({\left(\frac{⟨{\Delta}_{N}⟩}{{\Delta}_{N}^{0}}\right)}^{2}+{\left(\frac{{\Delta}_{{T}_{1}}}{{\Delta}_{T}^{0}}\right)}^{2}+{\left(\frac{{\Delta}_{{T}_{2}}}{{\Delta}_{T}^{0}}\right)}^{2}=1\)

On a donc bien un critère elliptique dans le cas général anisotrope où les contraintes critiques et les sauts de déplacement associés sont différents pour les modes normal et tangent.

Ce critère est accessible via la valeur du mot-clé CRIT_INIT = ‘YE’ dans DEFI_MATERIAU / RUPT_TURON.

Le critère proposé dans la loi CZM_LIN_REGest un cas particulier de ce critère, et correspond au cas où le joint est isotrope, c’est-à-dire lorsque les propriétés normale et tangente sont égales.

  • Le critère de Turon

C’est un critère de Benzeggagh – Kenane pour l’initiation. Bien qu’il est difficile d’avoir des données expérimentales validant le choix de l’un ou l’autre des deux critères évoqués, l’idée est ici de formuler les surfaces seuil d’initiation et de propagation de manière dépendante. Il s’écrit:

\({\sigma}_{c}^{2}={\sigma}_{c,N}^{2}+\left({\sigma}_{c,T}^{2}-{\sigma}_{c,N}^{2}\right){B}^{\eta}\)

et, de manière équivalente:

\({({\Delta}^{0})}^{2}={({\Delta}_{N}^{0})}^{2}+\left({({\Delta}_{T}^{0})}^{2}-{({\Delta}_{N}^{0})}^{2}\right){B}^{\eta}\)

Ce critère est accessible via la valeur du mot-clé CRIT_INIT = ‘TURON’ dans DEFI_MATERIAU.

Visualisations des critères#

Pour visualiser le couplage des modes par les critères de Ye ou de Turon, on peut regarder les surfaces seuil dans l’espace des déplacements, et voir notamment l’influence du paramètre :math:`eta ` sur le seuil.

On rappelle les propriétés normale et tangente issues de la référence [ 1 ], caractérisant l’endommagement intralaminaire d’un matériau composite (i.e. à l’une des interfaces entre les différents «plis» qui constituent le composite: c’est donc un endommagement interne au composite):

\(K=1.0\times {10}^{15}N/{m}^{3}\)

\({\sigma}_{c,N}=80\times {10}^{6}\mathit{Pa}\)

\({G}_{c,N}=969N/m\)

\({\sigma}_{c,T}=100\times {10}^{6}\mathit{Pa}\)

\({G}_{c,T}=1719N/m\)

ce qui donne les sauts de déplacements à l’initiation de l’endommagement suivant:

\({\Delta}_{N}^{0}=0.08\mu m\)

\({\Delta}_{T}^{0}=0.10\mu m\)

\({\Delta}_{N}^{f}=24\mu m\)

\({\Delta}_{T}^{f}=34\mu m\)

Avec ce jeu de paramètres, on voit que les propriétés normale et tangente ne sont pas extrêmement différentes (rapport 1.2 pour les propriétés d’initiation de l’endommagement, 1.4 pour les propriétés de propagation).

Sur les figures (a) et (b), on constate que les surfaces seuil du critère de Turon peuvent être au-dessus ou en-dessous du critère elliptique suivant les valeurs du paramètre puissance \(\eta ` et que les deux critères sont même très proches pour certaines valeurs de :math:\)eta ` . C’est notamment le cas pour la valeur issue de l’article qui a été calée à partir d’essais expérimentaux, à savoir \(\eta =2.84\) . L’effet du paramètre \(\eta ` est plus facile à appréhender sur les figures (c) et (d): plus le paramètre :math:\)eta ` est grand, plus on «tarde» à passer à la valeur du mode tangente à mesure que le taux de mixité augmente de 0 (mode normal pur) à 1 (mode tangent pur). On a également représenté, en pointillé, les critères «circulaires», c’est-à-dire si les propriétés normale et tangente sont égales, ce qui correspond au cas de la loi CZM_LIN_REG.

On peut visualiser ces mêmes courbes sur la figure avec des propriétés déterminées dans le cadre de la caractérisation expérimentale d’un joint interlaminaire (entre deux composites distincts) qu’on peut trouver typiquement dans une pâle d’éolienne. On ne donne pas ici explicitement les propriétés (les sauts d’initiation et d’endommagement peuvent néanmoins être lus sur les courbes): l’important est que le rapport entre les propriétés normale et tangente est plus important, surtout pour ce qui est de la propagation de la fissure. Sur la figure

Cette loi de comportement ayant été développée initialement pour modéliser l’endommagement de joints collés, on n’a pas considéré de cas où les propriétés normales sont supérieures aux propriétés tangentes, car ce n’est pas représentatif de la réalité: une sollicitation normale (pelage) est généralement plus critique pour un joint collé qu’une sollicitation en cisaillement.

../../../../_images/10000201000002BC000002BC211F6F88947E4A51.png ../../../../_images/1000020100000280000001E0AFDDB041C9413EB9.png
../../../../_images/10000201000002BC000002BC6C34F544215B564A.png ../../../../_images/1000020100000280000001E0342329C958C3E661.png

Figure 2.2.2: Visualisation des critères de Turon (pour différentes valeurs de η) et de Ye avec les propriétés matériau issues de [ 1 ]: (a) surfaces seuil d’initiation de l’endommagement; (b) évolution du saut d’initiation en fonction du taux de mixité; (c) surfaces seuil de propagation de la fissure; (d) évolution du saut de propagation en fonction du taux de mixité. Attention, les échelles selon X et Y sont différentes.

../../../../_images/10000201000002BC000002BCA9ED872048DFA853.png ../../../../_images/1000020100000280000001E0C685089825E76B87.png
../../../../_images/10000201000002BC000002BC2E82FCBD33B5F614.png ../../../../_images/1000020100000280000001E058409BE72DF7EA03.png

Figure 2.2.3: Visualisation des critères de Turon (pour différentes valeurs de η) et de Ye avec les propriétés matériau issues de d’une campagne de caractérisation de joint: (a) surfaces seuil d’initiation de l’endommagement; (b) évolution du saut d’initiation en fonction du taux de mixité; (c) surfaces seuil de propagation de la fissure; (d) évolution du saut de propagation en fonction du taux de mixité. Attention, les échelles selon X et Y sont différentes.

Cette formulation assure une transition continue entre entre les surfaces seuil d’endommagement, où la surface seuil associée à \(d=0\) est la surface d’initiation de l’endommagement, et celle associée à \(d=1\) la surface de propagation de la fissure.

On peut le visualiser sur la Figure, obtenu avec les valeurs de paramètres données ci-dessus et issues de [ 1 ] par le fait que les surfaces seuil (évoluant en fonction de l’endommagement \(d\) ) sont bien imbriquées les unes dans les autres à mesure que \(d\) augmente et ne se croisent pas. Ainsi, si la mixité du chargement change, le taux de mixité correspondant à l’angle dans le diagramme \(({\Delta}_{N},{\Delta}_{T})\) , on n’aura pas d’ «auto-cicatrisation» de l’élément cohésif (l’évolution de \(d\) est bien monotone croissante), et l’énergie de dissipation sera toujours positive. Pour un élément «auto-cicatrisant», l’endommagement \(d\) pourrait diminuer au changement de taux de mixité, comme cela est schématisé sur la Figure .

../../../../_images/10000201000001F20000014F8714A290D42D93F2.png

Figure 2.2.4: Surface seuil d’endommagement dans l’espace des déplacements (avec ici shear ≡ T et 3 ≡ N) [ 1 ]

../../../../_images/100002010000021E000000D0C46026FFD3725C3D.png

Figure 2.2.5: Exemple d“« auto-cicatrisation » d’un élément cohésif lorsque le taux de chargement évolue de la valeur A à la valeur B [ 1 ]

Matrice tangente#

La matrice tangente est nécessaire pour l’implémentation numérique de la loi de comportement dans le cadre de la résolution des équations mécaniques par l’algorithme de Newton.

La matrice tangente s’obtient par différenciation de l’équation de comportement:

../../../../_images/Shape410.gif

oùon a utilisé

../../../../_images/Shape510.gif

On définit

../../../../_images/Shape610.gif

.

Concentrons-nous sur la quantité

../../../../_images/Shape710.gif

:

../../../../_images/Shape812.gif

En régime dissipatif on a:

../../../../_images/Shape912.gif

Ici, on vient négliger la variation des quantités

../../../../_images/Shape1010.gif

et

../../../../_images/Shape1112.gif

en fonction de

../../../../_images/Shape126.gif

, qui en dépendent via le taux de mixité. Donc on fait l’hypothèse que le taux de mixité de chargement n’évolue pas trop brutalement. Si nécessaire, on veillera à discrétiser suffisamment le chargement pour satisfaire cette hypothèse au cours d’un calcul.

et par ailleurs:

../../../../_images/Shape135.gif

où on a utilisé

../../../../_images/Shape144.gif

Pour finir on obtient la matrice tangente

../../../../_images/Shape154.gif

suivante:

  • en régime dissipatif

    ../../../../_images/Shape164.gif

:

../../../../_images/Shape173.gif

  • en régime linéaire

    ../../../../_images/Shape183.gif

ou

../../../../_images/Shape193.gif

:

../../../../_images/Shape202.gif

avec, on le rappelle:

../../../../_images/Shape219.gif

On constate que la matrice

../../../../_images/Shape223.gif

n’est pas définie en

../../../../_images/Shape232.gif

. On prolongera la fonction à droite en ce point.

Variables internes de la loi CZM_TURON#

La loi CZM_TURON possède variables internes. Du point de vue de l’implémentation du comportement, seule la première est une variable interne au sens thermodynamique. Les autres sont des post-traitements qui fournissent des indications sur l’état du joint à un instant donné.

! V1 : PLUS GRANDE NORME DU SAUT ÉQUIVALENT TOTAL (variable

../../../../_images/Shape242.gif

telle que définie p. 5 )

! V2 : VARIABLE SEUIL (variable

../../../../_images/Shape252.gif

telle que définie p. 5 )

! V3 : VARIABLE D’ENDOMMAGEMENT (variable d)

! V4 : INDICATEUR DE DISSIPATION (0 : SI REGIME LIN, 1 : SI REGIME DISS)

! V5 : INDICATEUR D’ENDOMMAGEMENT (0 : SAIN, 1: ENDOMMAGE, 2: CASSE)

! V6 A V8 : VALEURS DU SAUT DANS LE REPERE LOCAL

! V9 : VALEUR DU SAUT EQUIVALENT TOTAL (variable

../../../../_images/Shape263.gif

)

! V10 : VALEUR DU SAUT EQUIVALENT TANGENTIEL (variable

../../../../_images/Shape272.gif

)

! V11 : POURCENTAGE D’ENERGIE DISSIPEE

! V12 : VALEUR DE L’ENERGIE DISSIPEE

! V13 : TAUX DE MIXITE BETA A T+ (CALCULE PAR LES SAUTS)

! V14 : TAUX DE MIXITE B A T+ (CALCULE PAR LES TAUX DE RESTITUTION D’ENERGIE)

! V15 : SEUIL D’INITIATION DE L’ENDOMMAGEMENT EN MODE MIXTE A T+ (variable

../../../../_images/Shape282.gif

)

! V16 : SEUIL DE PROPAGATION DE LA FISSURE EN MODE MIXTE A T+ (variable

../../../../_images/Shape292.gif

)

Pilotage du chargement par prédiction élastique#

Le pilotage du chargement, est disponible pour les lois de comportement cohésives. De type prédiction élastique: PRED_ELAS (voir [R5.03.80]), il permet de suivre des solutions dissipatives présentant des instabilités. Cette technique est décrite précisément dans la documentation [R7.02.11] concernant d’autres types de lois cohésives disponibles dans code_aster. On pourra également se référer à [ 4 ] pour plus de détails.

Mots-clés#

La loi de comportement CZM_TURON peut être utilisée dans STAT_NON_LINE et DYNA_NON_LINE. Elle est utilisable sur des éléments de joint 2D ou 3D avec les modélisations: PLAN_JOINT, AXIS_JOINT ou 3D_JOINT.

Quatre paramètres sont à saisir dans DEFI_MATERIAUsous les mots clé RUPT_TURON:

GC_N: densité d’énergie de surface critique en mode N pur (Pression x Longueur)

GC_T: densité d’énergie de surface critique en mode T pur (Pression x Longueur)

SIGM_C_N: contrainte critique en mode N pur (Pression)

SIGM_C_T: contrainte critique en mode T pur (Pression)

K:rigidité initiale de l’interface (Pression x Longueur -1)

ETA_BK: paramètre puissance dans la loi de Benzeggagh-Kenane

C_RUPT: paramètre de rigidité résiduelle post-rupture (sans dimension, valant 0.001 par défaut et 0.1 au maximum)

CRIT_INIT: paramètre de sélection du critère d’initiation de l’endommagement (vaut ‘YE’ou ‘TURON’)

Commandes

STAT_NON_LINE COMPORTEMENT/RELATION CZM_TURON

DYNA_NON_LINE COMPORTEMENT/RELATION CZM_TURON

AFFE_MODELE MODELISATION PLAN_JOINT

AXIS_JOINT

3D_JOINT

DEFI_MATERIAU RUPT_TURON GC_N

GC_T

SIGM_C_N

SIGM_C_T

K

ETA_BK

C_RUPT

CRIT_INIT

Cette loi de comportement est testée en statique:voir cas test SSNV110 [V6.04.110].

Références#

  1. TURON A., CAMANHO P. P., COSTA J., DAVILA C. G., A damage model for the simulation of delamination in advanced composites under variable-mode loading. Mechanics of Materials 38 , pp. 1072-1089(2006)

  1. CHARLOTTE M., LAVERNE J. , MARIGO J.-J, Initiation of cracks with cohesive force models: a variational approach. European Journal of Mechanics, A/Solids , v 25, n 4, July/August, 2006, p 649-669

  1. BARENBLATT G. I., The mathematical theory if equilibrium cracks in brittle fracture. Adv. Appl. Mech ., 7, pp. 55-129 (1962)

  1. LAVERNE J.: “Formulation énergétique de la rupture par des modèles de forces cohésives: considérations théoriques et implantations numériques“. Thèse de l’Université Paris XIII, novembre 2004 (sous la direction de J.-J. MARIGO et E. LORENTZ)

  1. LORENTZ E., A mixed interface finite element for cohesive zone models, Comput. Methods Appl. Mech. Engrg. 198 (2008), 302-317

  1. R3.06.13: «Eléments finis d’interface mixte pour des modèles de zone cohésive (xxx_INTERFACEet xxx_INTERFACE_S), documentation de code_aster

  1. R7.02.11 : «Lois de comportement cohésives : CZM_xxx_xxxet pilotage du chargement», documentation de code_aster