r7.01.47 Loi d’endommagement ENDO_LOCA_TC

Résumé:

Ce document décrit le modèle de comportement isotrope quasi-fragile ENDO_LOCA_TC. Il reprend les fondamentaux des modèles ENDO_ISOT_BETON [R7.01.04] et ENDO_LOCA_EXP [R7.01.42] et les enrichit en introduisant :

• un seuil d’endommagement rationnalisé dans l’optique d’optimiser la robustesse numérique ;

• un effet de l’endommagement en compression, de sorte à y borner les contraintes, sans pour autant chercher à modéliser la phase post-pic. La contrainte seuil en compression tend de manière monotone vers une limite, la résistance en compression.

Ce document synthétise les informations principales du modèle pour son utilisation dans code_aster en s’appuyant sur la note interne [Lorentz-2025]. L’utilisateur désirant avoir plus d’informations sur le modèle (descriptif complet des équations, propriétés mécaniques ou méthématiques…) est invité à se reporter à cette référence.

Phénoménologie & Utilisation

Phénomènes couverts et choix de modélisation

Le modèle d’endommagement ENDO_LOCA_TC est isotrope. Il ne rend donc compte ni d’une anisotropie initiale – l’élasticité initiale y est isotrope de même que le seuil d’endommagement – ni d’une anisotropie induite : l’endommagement s’accumule sans garder la mémoire des directions de sollicitations qui conduisent à son évolution, si bien que l’élasticité et le seuil d’endommagement courants restent isotropes.

Ce cadre isotrope permet néanmoins de distinguer des effets de traction vs. compression. En effet, parmi les directions de sollicitations (les directions propres des déformations), il est possible d’identifier celles correspondant à des déformations principales positives (respectivement négatives) : les directions de traction (resp. compression). Cette distinction permet de séparer la contrainte élastique (celle qu’on observerait sous la même déformation en l’absence d’endommagement) en deux termes, l’un relatif à la traction et l’autre la compression, qui seront affectés différemment par l’endommagement pour fournir la contrainte effectivement subie par le point considéré. A endommagement fixé, le modèle est donc élastique non linéaire.

Le modèle rend compte de deux contributions distinctes à l’endommagement. Une première, dite endommagement de traction, est pilotée par les déformations de traction et ne pondère que la contrainte élastique de traction. Pour des déformations de traction croissantes, ce coefficient de pondération (noté A), initialement égal à 1, le reste jusqu’à un certain seuil de résistance en traction, puis décroît jusqu’à devenir nul (endommagement ultime total) pour une déformation suffisamment élevée mais finie si bien que seule la contrainte élastique de compression contribue encore à la contrainte. Cela signifie qu’à endommagement de traction ultime, le point ne résiste plus à des sollicitations en traction mais uniquement en compression.

Une seconde composante d’endommagement est pilotée par les déformations de compression (plus précisément la contrainte élastique de compression). C’est pourquoi elle est appelée endommagement de compression, même si elle affecte non seulement la contrainte élastique de compression mais aussi celle de traction (elle pondère donc toute la contrainte élastique). Pour des déformations de compression croissantes (au signe près), ce coefficient de pondération (noté C), initialement égal à 1, le reste jusqu’à un certain seuil de perte de linéarité en compression, puis décroît et tend asymptotiquement vers 0 mais sans jamais atteindre la valeur nulle. Son évolution en fonction de la déformation est telle que les contraintes restent bornées. Mais à la différence de la traction, elles continuent de croître (pas d’adoucissement en compression), jusqu’à atteindre une valeur ultime, la limite de résistance en compression.

Le modèle ENDO_LOCA_TC ne présente pas d’autres mécanismes dissipatifs que les endommagements de traction et de compression. En particulier, il ne rend compte d’aucune déformation irréversible. Enfin, il est possible de tenir compte de l’influence sur la réponse mécanique d’autres phénomènes physiques comme la thermique ou les transferts hydriques au sein du béton. D’une part, la température peut influer sur les paramètres caractéristiques du modèle. Et d’autre part, ce dernier est exprimé en déformations mécaniques, c’est-à-dire que les déformations d’origine thermique ou de séchage viennent se soustraire à la déformation totale.

Paramètres du modèles

En sus des paramètres élastiques (module Young et coefficient de Poisson), les paramètres du modèle sont au nombre de cinq :

• \(f_{t}\) la résistance en traction avec \(f_{t} > 0\)

• \(p\) un coefficient d’écrouissage en traction, avec en pratique \(p\in\left[-0.2,9\right]\)

• \(\bar{\omega}\) l’énergie à rupture en traction confinée normalisée par l’énergie élastique à l’initiation de l’endommagement avec \(\bar{\omega}> \frac{4}{3\pi}\left(p+2\right)^{3/2}\)

• \(f_{c}\) la résistance en compression (moyenne) avec \(f_{c} > 2\,f_{t}\)

• \(\sigma^{0}\) la contrainte seuil d’amorçage en compression avec \(2\,f_{t} \leq \sigma^{0} < f_{c}\)

En l’absence d’essais qui permettraient un recalage plus précis, le fib Model Code ([fib-2010]) fournit des ordres de grandeur sur la base de corrélations pour le module de Young E (en MPa), le coefficient de Poisson \(\nu\), l’énergie de fissuration \(G_{f}\) (en N/m) et la résistance en traction \(f_{t}\) (en MPa) en fonction de la résistance en compression \(f_{c}\) (en MPa) :

(4189)\[ E = 21500\,\left(\frac{f_{c}}{10}\right)^\left(1/3\right) \quad ; \quad \nu = 0.2 \quad ; \quad G_{F}=73\,f_{c}^{0.18} \quad ; \quad f_{t}=0.3\,\left(f_{c}-8\right)^\left(2/3\right)\]

En cohérence avec l’estimation du module de Young sécant proposée dans Eurocode 2 ([Eurocode-2-2005]), la valeur du seuil d’amorçage en compression peut être choisie par :

(4190)\[ \sigma^{0}=0.4\,f_{c}\]

La distance inter-fissures moyenne \(L_{F}\) peut être estimée à priori en fonction des armatures et du schéma d’homogénéisation simpliste illustré en Fig. 280 :

Fig. 280 Illustration du principe d’équivalence entre modèle homogène et modèle cohésif

Ce schéma d’homogénéisation 1D simpliste, décrit dans [Lorentz-2023], permet de passer de la réponse cohésive d’une fissure (CZM, déduite du modèle non local) à la réponse moyenne du VER constitué d’une fissure et d’une zone élastique environnante. Les réponses des deux modèles sont construites de sorte que sous une contrainte \(\sigma\) donnée, l’allongement du VER, noté \(\triangle U\), est identique. On voit apparaître de manière eplicite le rôle de la distance moyenne inter-fissures \(L_{F}\) (i.e. la taille de la cellule VER) qui établit le lien entre l’allongement \(\triangle U\) et la déformation moyenne \(\frac{\triangle U}{L_{F}}\) qui est la variable d’entrée du modèle local.

Avec cette valeur de distance inter-fissures moyenne \(L_{F}\), l’énergie consommée normalisée \(\bar{\omega}\) s’exprime en fonction des autres paramètres comme suit :

(4191)\[ \bar{\omega}=\frac{2E_{C}\,G_{F}}{L_{F}\,f_{t}^2}\]

Avec le module Young en traction confiné défini par :

(4192)\[ E_{c}=\lambda+2\mu=\frac{\left(1-\nu\right)\,E}{\left(1-2\nu\right)\left(1+\nu\right)}\]

Enfin, pour ce qui concerne le coefficient d’écrouissage \(p\), une valeur de - 0.2 conduit à un écrouissage (négatif) quasi-linéaire, tandis que des valeurs comprises entre 1 et 4 sont plus représentatives du creusement des courbes de traction présentées dans [Peterson-1981], voir [Lorentz-2017]. Notons d’ailleurs qu’un coefficient d’écrouissage \(p=1.5\) correspond à la réponse la plus proche du modèle cohésif bilinéaire « moyen » introduit dans [Peterson-1981].

Sur la base de ces suggestions, un jeu de valeurs typiques pour un béton de résistance \(f_{c}=40\,MPa\) et une distance inter-fissure \(L_{F} =30\,cm\) est fourni dans Tableau 67.

Tableau 67 Jeu de paramètres caractéristiques

Nom	Symbole	Valeur	Unité
Module de Young	\(E\)	34000	MPa
Coefficient de Poisson	\(\nu\)	0.2	–
Résistance en traction	\(f_{t}\)	3	MPa
Energie à rupture normalisée	\(\bar{\omega}\)	3.5	–
Coefficient d’écrouissage en traction	\(p\)	1	–
Résistance en compression	\(f_{c}\)	40	MPa
Seuil d’amorsage en compression	\(\sigma^{0}\)	16	MPa
Distance inter-fissure	\(L_{F}\)	30	cm

Identification des paramètres

Sans compter les paramètres élastiques et l’activation de viscosité interne de la loi, le modèle ENDO_LOCA_TC dépend des 5 paramètres détaillés dans le paragraphe précédent. Ces paramètres peuvent être identifiés via les réponses du modèle en traction simple (\(f_{t}\) résistance en traction, \(p\) coefficient d’écrouissage en traction et \(\bar{\omega}\) énergie à rupture normalisée) et en compression simple (\(f_{c}\) résistance en compression, \(\sigma^{0}\) seuil d’amorsage en compression). Les valeurs dans Tableau 67 sont en partie accessibles directement depuis les courbes ci-dessous illustrant le comportement du béton en traction et compression simple.

Fig. 281 Illustration de la réponse du modèle en traction simple

Ainsi, \(f_{t}\) est la contrainte maximale atteinte tandis que :math`p` et \(\bar{\omega}\) doivent être calibrés pour retrouver la courbe post-pic en traction simple.

Fig. 282 Illustration de la réponse du modèle en compression simple

Dans l’essai de compression simple, \(\sigma^{0}\) correspond à la contrainte de perte de linéarité de la réponse du béton et \(f_{c}\) est la contrainte maximale atteinte (contrainte au pic).

On rappelle également la relation (4193) entre l’énergie de fissuration \(G_{F}\) et l’énergie volumique à rupture en traction confinée \(\omega\) (relié à \(\bar{\omega}\) par (4191)), dans laquelle \(L_{F}\) désigne la distance moyenne inter-fissure :

(4193)\[ G_{F}=L_{F}\,\omega\]

Il est possible d’utiliser DEFI_MATER_GC [U4.42.07] pour trouver les valeurs des paramètres de ENDO_LOCA_TC à l’aide de données caractéristiques du béton ou via les recommandations du FIB model code.

Illustraction du comportement mécanique

Quelques caractéristiques du modèle sont illustrées dans cette section en utilisant les valeurs de Tableau 67. Dans un premier temps les réponses en situation de traction et compression simples en enchaînant traction puis compression (Fig. 283) ou bien compression puis traction (Fig. 285) sont présentées :

Fig. 283 Réponse du modèle sous sollicitations uniaxiale. Traction puis compression simple.

Fig. 284 Réponse du modèle sous sollicitations uniaxiale avec zoom sur la contrainte en traction. Traction puis compression simple.

Fig. 285 Réponse du modèle sous sollicitations uniaxiale. Compression simple puis traction.

Fig. 286 Réponse du modèle sous sollicitations uniaxiale avec zoom sur la contrainte en traction. Compression simple puis traction.

En traction (Fig. 283), le modèle conduit à un adoucissement de la réponse après franchissement d’un pic de contrainte (à la valeur de la résistance en traction \(f_{t}\)). Il n’y a pas de phase d’endommagement pré-pic. L’endommagement ultime est total, c’est-à-dire que la contrainte devient nulle. Il est atteint pour une valeur finie de la déformation.

Les décharges conduisent à une réponse élastique sécante, c’est-à-dire sans déformations irréversibles. En enchaînant en compression, on observe une restauration d’une grande partie de la rigidité et une réponse relativement proche de celle qu’aurait le matériau sans phase préalable d’endommagement en traction.

Vis-à-vis de la compression (Fig. 285), une perte de linéarité au-delà de la contrainte seuil \(\sigma^{0}\) (au signe près) est apparente. La contrainte continue de croître et tend vers une valeur asymptotique qui est la résistance en compression \(f_{c}\). Il n’y a pas de phase d’adoucissement (l’écrouissage reste positif). Par ailleurs, les décharges conduisent là-aussi à une réponse élastique sécante : il n’y a pas non plus de déformations irréversibles en compression. Enfin, en enchaînant par un chargement de traction, on observe que l’endommagement en compression a un effet significatif sur les propriétés en traction : la rigidité, le seuil d’endommagement et l’énergie dissipée (aire sous la courbe déformation-contrainte) sont nettement abaissés par rapport à ce qu’ils auraient été sans compression préalable.

Pour illustrer le rôle de chacun des paramètres, leur impact sur les réponses en traction simple et en compression simple est examiné. La courbe noir correspond à la réponse nominale avec le jeu de paramètres de Tableau 67. La courbe en bleu correspond à une valeur plus basse de l’un des paramètres (les autres conservant leur valeur nominale) et la courbe rouge correspond à une valeur plus haute de ce paramètre.

Fig. 287 Influence de la limite en traction sur la courbe en traction simple

Fig. 288 Influence de l’énergie à rupture normalisée en traction confinée sur la courbe en traction simple

Fig. 289 Influence du coefficient d’écrouissage en traction sur la courbe en traction simple

Fig. 290 Influence de la limite en compression sur la courbe en compression simple

Fig. 291 Influence de la contrainte seuil d’amorçage en compression sur la courbe en compression simple

On illustre maintenant l’ensemble des contraintes atteignables, représenté dans l’espace des contraintes principales (Fig. 292) puis en coupe dans le plan des contraintes planes (Fig. 293).

Fig. 292 Domaine des contraintes atteignables dans l’espace des contraintes principales.

Fig. 293 Domaine des contraintes atteignables dans l’espace des contraintes planes.

L’examen de l’ensemble des contraintes atteignables en Fig. 292 montre que le modèle ne conduit à aucune contrainte supérieure à \(f_{t}\) ni inférieure à \(-f_{c}\) : les contraintes sont bornées quelle que soit le trajet de chargement.

Dans le domaine des tractions dans toutes les directions (toutes les contraintes principales sont positives), la surface ultime suit la forme d’un critère de Rankine, ce qui est conforme aux observations expérimentales pour le béton [Lee-et-al-2004].

Dans le domaine des compressions dans toutes les directions (toutes les contraintes principales sont négatives), la surface ultime procède également d’un critère de Rankine. En contraintes planes (Fig. 293), c’est une approximation acceptable de la réalité expérimentale du béton. En revanche, l’effet de confinement en tri-compression n’est pas modélisé, c’est-à-dire qu’il n’y a pas de ductilité plus élevée dans cette direction de chargement.

Aide à la convergence

Les calculs d’endommagement peuvent être difficiles à faire converger car il existe un risque de propagation brutale de l’endommagement liée au caractère quasi-fragile du béton en endommagement de traction. Deux stratégies numériques permettent de retrouver une forme de continuité de l’évolution de l’endommagement par rapport au chargement appliqué.

La première consiste à adapter l’intensité du chargement indépendamment de sa chronologie réelle ; on parle alors de pilotage du chargement (par une méthode de continuation). Une fonction de pilotage bien adaptée, baptisée pilotage par prédiction élastique dans [Lorentz-et-Badel-2004], consiste à contrôler l’incrément maximal d’endommagement dans la structure pendant le pas de temps courant. Cette méthode, désignée par PRED_ELAS dans code_aster, peut être utilisée en l’absence de variables de commande (telles qu’un transitoire thermique) et de chargements dépendant du temps physique, voir document [R5.03.80].

Une stratégie numérique alternative consiste à introduire une part de viscosité dans l’évolution de l’endommagement. Elle est caractérisée par un paramètre noté \(\tau\) qui est homogène à un temps : plus il est élevé, plus l’effet de viscosité est prononcé. La viscosité induira une évolution continue de l’endommagement ce qui aidera à la convergence mais elle entraînera un retard de l’endommagement. Le temps caractéristique de viscosité :math:`tau`correspond à TAU_REGU_VISC dans DEFI_MATERIAU. Il est facultatif et requiert une connaissance avancée du modèle pour être employé.

Variables internes

Les deux composantes de l’endommagement, A et C, sont stockées dans les variables internes du modèle (plus précisément, ENDOTRAC correspond à \(1-A\) et ENDOCOMP à \(1-C\), pour garder l’usage courant d’un endommagement initialement nul qui évolue jusqu’à 1). Ce sont les deux seules « vraies » variables internes.

D’autres variables sont également rangées dans la liste des variables internes mais elles correspondent à des post-traitements des deux précédentes. Ainsi, les variables A et C peuvent être combinées pour refléter la perte de rigidité relative dans la direction de sollicitation courante (pour laquelle se superposent généralement des sollicitations à la fois en traction et en compression) : c’est la variable « interne » ENDO. Dans le même registre, ENDOTOT traduit la perte de rigidité la plus pénalisante en considérant toutes les directions de sollicitation ; elle vaut \(1\,-\,A\,C\).

Par ailleurs, deux autres mesures (équivalentes) de l’endommagement (HISTTRAC et HISTCOMP) sont également rangées dans la liste des variables internes. Il s’agit des endommagements « bruts » directement pilotés par les seuils en traction et en compression (ce sont les variables qu’on notera ultérieurement \(a\) et \(c\) où HISTTRAC correspond à \(a\) et HISTCOMP à \(c-1\)), conservés essentiellement pour des raisons numériques pratiques.

Enfin, les énergies de déformation correspondant aux contraintes élastiques de traction (ENERTRAC) et de compression (ENERCOMP) sont également stockées parmi les variables « internes » du modèle. En présence de viscosité, une dernière variable interne (SIGMVISC) fournit une estimation du retard pris par l’endommagement pour le niveau de déformation courant. Il s’agit de l’écart entre la contrainte courante et la contrainte qu’on observerait en l’absence de viscosité, c’est-à-dire pour un endommagement sans retard. L’écart est évalué en termes de contrainte principale maximale.

Le tableau suivant synthétise les noms des variables internes :

Tableau 68 Variables internes de ENDO_LOCA_TC

Nom	Numéro	Signification
ENDO	V1	Endommagement correspondant à la perte de rigidité relative dans la direction de sollicitation actuelle
ENDOTRAC	V2	Rigidité en traction égal à \(1-A\) avec A la composante d’endommagement de traction
HISTTRAC	V3	Endommagement en traction correspondant à a
ENERTRAC	V4	Energie de déformation correspondant aux contraintes élastiques de traction
ENDOCOMP	V5	Rigidité en compression égal \(1-C\) avec C la composante d’endommagement de compression
HISTCOMP	V6	Endommagement en compression correspondant à \(c-1\)
ENERCOMP	V7	Energie de déformation correspondant aux contraintes élastiques de compression
SIGMVISC	V8	Ecart entre la contrainte courante et la contrainte observable en l’absence de viscosité
ENDOTOT	V9	Endommagement le plus pénalisant en examinant toutes les directions de sollicitation

Equations du modèles

Relation déformation-contrainte

L’état mécanique est décrit par la déformation \(\tensTwo{ \varepsilon}\), un endommagement de traction \(a\in\left[0,1\right]\) qui évolue de 0 (sain) à 1 (rompu) et un endommagement de compression \(c\) qui évolue de 1 (sain) à \(+\infty\) (saturé). L’énergie de déformation est définie par :

(4194)\[ w \left(\tensTwo{\varepsilon},a,c\right)=C\left(c\right)\left[A\left(a\right) w^{+}\left(\tensTwo{ \varepsilon}\right)+ w^{-}\left(\tensTwo{ \varepsilon}\right) \right] \quad ; \quad w^{+/-}\left(\tensTwo{\varepsilon}\right)=\frac{1}{2}\left[\lambda\langle trac\left(\tensTwo{\varepsilon}\right)\rangle^{2}_{+/-}+2\mu\sum_1^3\langle\varepsilon_{i}\rangle^{2}_{+/-}\right]\]

\(w^{+}\) et \(w^{-}\) correspondent à une partition de l’énergie élastique \(w_{e}\) entre les effets de traction et de compression tel que :

(4195)\[ w_{e}=w^{+}+w^{-}\]

La relation contrainte – déformation dérive du potentiel précédent :

(4196)\[ \sigma=\frac{\partial w}{\partial \tensTwo{ \varepsilon}}=C\left(c\right)\left[A\left(a\right)\tensTwo{ s}^{+}+\tensTwo{ s}^{-} \right]\]

Où

(4197)\[ \tensTwo{ s}^{+}=\frac{\partial w^{+}}{\partial \tensTwo{ \varepsilon}}\left(\tensTwo{ \varepsilon}\right) \quad ; \quad \tensTwo{ s}^{-}=\frac{\partial w^{-}}{\partial \tensTwo{ \varepsilon}}\left(\tensTwo{ \varepsilon}\right)\]

\(A\left(a\right)\) est la fonction de rigidité normalisée en traction confinée et elle est choisie identique à celle du modèle ENDO_LOCA_EXP ([R7.01.42]).

(4198)\[ A\left(a\right)=\frac{1-a}{e\left(a\right)} \quad ; \quad e\left(a\right) = 1-a+\frac{\bar{\omega}}{2}\left[m_{0}a+\left(D_{1}-m_{0}\right)a^{r}\right]\]

\(m_{0}\), \(D_{1}\) et r s’exprime uniquement à l’aide du paramètre \(p\) de la loi ENDO_LOCA_TC :

(4199)\[ m_{0}=\frac{3\pi}{2}\left(p+2\right)^{-3/2}, D_{1}=\frac{3\pi}{4}\sqrt{1+p}\quad ; \quad r= \frac{2\left(D_{1}-1\right)-m_{0}}{2-m_{0}}\]

Ainsi, les paramètres introduits se limitent au coefficient \(p\) qui influe sur l’allure de la réponse post-pic, d’autant plus creusée que \(p\) est élevé et à l’énergie à rupture normalisée en traction confinée définie par

(4200)\[ \bar{\omega}=\frac{\omega}{w_{c}}\]

Où \(\omega\) est l’énergie volumique à rupture du béton dans un essai de traction confinée et \(w_{c}\) est l’énergie élastique à l’initiation de l’endommagement en traction confinée.

Le graphe de la fonction de rigidité A pour le jeu de valeurs de Tableau 67 est présenté en Fig. 294.

Fig. 294 Fonction de rigidité de traction A en fonction de l’endommagement traction a.

\(C\left(c\right)\) est la rigidité résiduelle normalisée de compression. Elle affecte toute la rigidité sous l’effet de la compression. Elle vaut :

(4201)\[ C\left(c\right)=\frac{1}{c}\left[1+\left(\frac{f_{c}}{\sigma^{0}}-1\right)\tanh\left(\frac{c-1}{\frac{f_{c}}{\sigma^{0}}-1}\right)\right]\quad ; \quad f_{c}>\sigma^{0}\]

Deux nouveaux paramètres sont introduits : la résistance en compression du béton \(f_{c}\) et le seuil d’amorçage en compression \(\sigma^{0}\). Le graphe de la fonction C pour le jeu de valeurs de Tableau 67 est présenté en Fig. 295.

Fig. 295 Fonction de rigidité de traction C en fonction de l’endommagement de compression c.

Evolution de l’endommagement de traction

La déformation relative à l’endommagement de traction est définie par :

(4202)\[ \chi_{t}=\langle\frac{s_{3}}{f_{t}}\rangle_{+} \quad ; \quad s_{3}=\lambda_{3}\left(\tensTwo{ s}\right)\]

\(f_{t}\) est un paramètre de la loi correspondant à la limite en traction du béton. \(s_{3}=\lambda_{3}\left(\tensTwo{ s}\right)\) désigne la valeur propre maximale du tenseur des contraintes élastiques \(\tensTwo{ s}\) défini comme la dérivée de la densité d’énergie élastique ((4203)).

(4203)\[ \tensTwo{ s} =\frac{\partial w_{e}}{\partial \tensTwo{\varepsilon}}\left(\tensTwo{\varepsilon}\right)=s^{+}+s^{-}\]

La mesure de déformation \(\chi_{t}\) permet de construire la fonction seuil pour l’endommagement de traction, ainsi que l’équation d’évolution correspondante. L’équation régissant la fonction seuil en traction est (4204) :

(4204)\[ g_{t}\left(\tensTwo{{\varepsilon}},a\right) = \chi_{t}\left(\tensTwo{{\varepsilon}}\right)^{2}-e\left(a\right)^{2}\]

avec \(e\left(a\right)\) le seuil d’endommagement défini par (4198). La loi d’évolution de l’endommagement en traction est (4205) :

(4205)\[ g_{t}\left(\tensTwo{{\varepsilon}},a\right) \leq 0\quad ; \quad \dot{a}\geq0\quad ; \quad \dot{a}\,g_{t}\left(\tensTwo{{\varepsilon}},a\right) = 0\]

Pour un état initial non endommagé, on a \(a\left(0\right)=0\)

Evolution de l’endommagement de compression

L’évolution de la rigidité résiduelle relative en compression C est gouvernée par un mécanisme similaire à celui de la traction. La mesure de déformation relative à l’endommagement en compression est définie comme suit :

(4206)\[ \chi_{c}=\langle-\,\frac{s_{1}^{-}}{\sigma^{0}}\rangle_{+} \quad ; \quad s_{1}^{-}=\lambda_{1}\left(\tensTwo{ s}^{-}\right)\]

Où \(\sigma^{0}\geq 2\,f_{t}\) afin de garantir certaines propriétés du modèle comme l’absence d’initiation d’endommagement en tration dans une sollicitation en compression. \(\lambda_{1}\left(\tensTwo{ s}^{-}\right)\) est la la plus petite valeur propre du tenseur \(\tensTwo{ s}^{-}\) défini en (4197). L’évolution de l’endommagement est alors gouvernée par l’intermédiaire de la fonction seuil suivante :

(4207)\[ g_{c}\left(\tensTwo{\varepsilon},c\right) = \chi_{c}\left(\tensTwo{\varepsilon}\right)-c\]

La loi d’évolution de l’endommagement en compression s’écrit :

(4208)\[ g_{c}\left(\tensTwo{\varepsilon},c\right) \leq 0\quad ; \quad \dot{c}\geq0\quad ; \quad \dot{c}\,g_{c}\left(\tensTwo{\varepsilon},c\right) = 0\]

Pour un état initial non endommagé, on a \(c\left(0\right)=1\). Contrairement à la mesure d’endommagement en traction a qui évolue de 0 à 1, l’endommagement en compression c évolue de 1 à \(+\infty\).

Implémentation numérique

De manière générale, l’intégration d’une relation de comportement consiste à déterminer l’histoire fonction du temps des variables internes, ici \(t\rightarrow a\left(t\right)\) et \(t\rightarrow c\left(t\right)\), et des contraintes \(t\rightarrow \tensTwo{\sigma}\left(t\right)\) pour une sollicitation en déformation \(t\rightarrow \tensTwo{\varepsilon}\left(t\right)\). Traditionnellement, cette opération est réalisée en procédant à une discrétisation en temps des équations différentielles qui caractérisent le comportement puis au traitement des systèmes algébriques (non linéaires) qui en résultent.

Equations discrétisées en temps

Pour la discrétisation en temps de la relation de comportement ENDO_LOCA_TC, on adopte un schéma d’Euler implicite, d’ordre 1. Pour un pas de temps \(\left[t_{n},t_{n+1}\right]\) et une quantité quelconque \(X\), on note \(X_{n}=X\left(t_{n}\right)\), \(X_{n+1}=X\left(t_{n+1}\right)\) et \(\triangle X =X_{n+1}-X_{n}\). La discrétisation de la relation contrainte – déformation (4196) via le schéma retenu conduit à :

(4209)\[ \tensTwo{\sigma}_{n+1}=C\left(c_{n+1}\right)\left[A\left(a_{n+1}\right)\tensTwo {s}_{n+1}^{+}+\tensTwo {s}_{n+1}^{-} \right]\]

Où

(4210)\[ \tensTwo {s}_{n+1}^{+}=\frac{\partial w^{+}}{\partial \tensTwo{ \varepsilon}}\left(\tensTwo{ \varepsilon}_{n+1}\right) \quad ; \quad \tensTwo {s}_{n+1}^{-}=\frac{\partial w^{-}}{\partial \tensTwo{ \varepsilon}}\left(\tensTwo {\varepsilon}_{n+1}\right)\]

La discrétisation de l’équation d’évolution en traction (4205) s’écrit :

(4211)\[ g_{t}\left(\tensTwo{\varepsilon}_{n+1},a_{n+1}\right) \leq 0\quad ; \quad \triangle a\geq0\quad ; \quad \triangle a \,g_{t}\left(\tensTwo{\varepsilon}_{n+1},a_{n+1}\right) = 0\]

Enfin, l’équation d’évolution en compression (4208) s’écrit après discrétisation :

(4212)\[ g_{c}\left(\tensTwo{\varepsilon}_{n+1},c_{n+1}\right) \leq 0\quad ; \quad \triangle c\geq0\quad ; \quad \triangle c_{n+1}\, g_{c}\left(\tensTwo{\varepsilon}_{n+1},c_{n+1}\right) = 0\]

Résolution du système algébrique

L’intégration du comportement consiste à déterminer \(a_{n+1}\) et \(c_{n+1}\) connaissant \(a_{n}\), \(c_{n}\), \(\varepsilon_{n}\) et \(\triangle \varepsilon\) puis à calculer \(\sigma_{n+1}\) via (4209). Il est intéressant de noter que le découplage entre les lois d’évolution permet de déterminer \(a_{n+1}\) et \(c_{n+1}\) de manière indépendante. Commençons par le calcul de l’endommagement de traction en fin de pas de temps \(a_{n+1}\) à la l’aide de la définition du seuil (4204):

(4213)\[ g_{t}\left(\tensTwo{\varepsilon}_{n+1},a_{n+1}\right) = \chi_{t}\left(\tensTwo{\varepsilon}_{n+1}\right)^{2}-e\left(a_{n+1}\right)^{2}\]

Comme la fonction \(g_{t}\) est strictement décroissante sur [0,1], (4211) admet une solution unique :

• si \(g_{t}\left(\tensTwo{\varepsilon}_{n+1},a_{n}\right)\leq 0\) alors \(a_{n+1}=a_{n}\).

• si \(g_{t}\left(\tensTwo{\varepsilon}_{n+1},1\right)\geq 0\) alors \(a_{n+1}=1\) car l’endommagement en traction est borné.

• sinon \(a_{n+1}\) est l’unique solution de \(g_{t}\left(\tensTwo{\varepsilon}_{n+1},a_{n+1}\right)= 0\) et compris strictement entre \(a_{n}\) et 1.

Pour résoudre le dernier cas, il est nécessaire de passer par la résolution d’une équation scalaire. Comme \(g_{t}\) est continûment dérivable, on peut mettre en œuvre une méthode de Newton avec contrôle des bornes.

Pour ce qui concerne l’évolution de l’endommagement en compression, le problème (4212) admet une solution unique. En effet, comme la fonction \(c\rightarrow g_{c}\left(\tensTwo{\varepsilon},c\right)=\chi_{c}\left(\tensTwo{\varepsilon}\right)-c\) est affine, la solution s’écrit simplement :

(4214)\[ c_{n+1}=max\left(c_{n},\chi_{c}\left(\tensTwo{\varepsilon}_{n+1}\right)\right)\]

Matrice tangente

ENDO_LOCA_TC est une loi de comportement formulée en déformation mécanique. Pour bien illustrer le rôle de ce choix, on considère ici à titre d’exemple l’effet de la température \(T\) à travers une déformation thermique \(\tensTwo{\varepsilon^{th}}\) et une dépendance des paramètres du modèle à la température. La contrainte issue du processus d’intégration du comportement du paragraphe précédent peut s’écrire formellement:

(4215)\[ \tensTwo{\sigma^{num}}\left(\tensTwo{\varepsilon},T\right)\]

Où la déformation mécanique est définie par :

(4216)\[ \tensTwo{\varepsilon}=\triangledown^{s}\vector{u}-\tensTwo{\varepsilon^{th}}\]

En suivant [Filiot-et-al-2017], l’algorithme de résolution des équations d’équilibre de Code_Aster s’appuie sur un développement limité à l’ordre 1 des contraintes en chaque point d’intégration. Dans le cas présent, il s’écrit :

(4217)\[ \tensTwo{\sigma^{sol}}_{n+1}=\tensTwo{\sigma^{num}}\left(\tensTwo{\varepsilon^{sol}}_{n+1},T_{n+1}\right)\approx\tensTwo{\sigma^{num}}\left(\tensTwo{\varepsilon},T_{n+1}\right)+\frac{\partial \tensTwo{\sigma^{num}}}{\partial \tensTwo{\varepsilon}}\left(\tensTwo{\varepsilon},T_{n+1}\right)\bf{:}\left(\tensTwo{\varepsilon^{sol}}_{n+1}-\tensTwo{\varepsilon}\right)\]

où \(\tensTwo{\varepsilon^{sol}}_{n+1}\) et \(\tensTwo{\sigma^{sol}}_{n+1}\) désignent respectivement la déformation et la contrainte solutions du problème à la fin du pas de temps. La déformation mécanique \(\tensTwo{\varepsilon}\) au voisinage de laquelle est effectuée le développement limité est la déformation mécanique au début du pas de temps \(\tensTwo{\varepsilon}_{n}\) en phase de prédiction de l’algorithme et le dernier itéré de la déformation mécanique en phase de correction.

Lors de l’intégration du comportement cette déformation mécanique \(\tensTwo{\varepsilon}\) est une donnée d’entrée, à l’instar de la température en fin de pas de temps \(T_{n+1}\) (chaînage et pas couplage de la thermique et de la mécanique). Quant à l’opérateur tangent (cohérent), il s’agit du facteur du terme d’ordre 1 dans le développement limité, à savoir :

(4218)\[ \tensFour{L}=\frac{\partial \tensTwo{\sigma^{num}}}{\partial \tensTwo{\varepsilon}}\left(\tensTwo{\varepsilon},T_{n+1}\right)\]

Pour calculer cet opérateur tangent, on procède par différentiation, dans l’optique d’identifier \(\tensFour{L}\) via \(d\tensTwo{\sigma}=\tensFour{L}\bf{:}d\tensTwo{\varepsilon}\). En partant de l’expression de la contrainte (4209) où l’indice n+1 est omis pour alléger les notations, on obtient l’équation suivante :

(4219)\[ d\tensTwo{\sigma}=dC\left[A\,\tensTwo{s}^{+}+\tensTwo{s}^{-}\right]+C\,dA\,\tensTwo{s}^{+}+C\left[A\frac{\partial^{2} w^{+}}{\partial \tensTwo{\varepsilon}^{2}}+\frac{\partial^{2} w^{-}}{\partial \tensTwo{\varepsilon}^{2}}\right]\bf{:}\,d\tensTwo{\varepsilon}\]

Dans cette expression, on exploite le caractère \(C^{2}\) de \(w^{+}\) et \(w^{-}\) en mettant en oeuvre une régularisation de ces dernières introduite dans [Lorentz-2025]. Les variations de A et C se calculent ensuite en différentiant les équations de cohérence correspondantes. En régime d’endommagement de traction actif (sinon \(dA=0\)), l’équation \(g_{t}\) issue de (4211) et (4213) conduit à :

(4220)\[ dg_{t}=2\chi_{t}d\chi_{t}-2e\left(a\right)e^{'}\left(a\right)da=0\]

et de là :

(4221)\[ dA=A^{'}\left(a\right)da\quad ; \quad da=\frac{\chi_{t}d\chi_{t}}{e\left(a\right)e^{'}\left(a\right)}\]

On procède de même pour évaluer \(dC\) en régime d’endommagement de compression actif (sinon \(dC=0\)), via la différentiation de l’équation de cohérence \(g_{c}=0\) issue de (4214).

(4222)\[ dC=C^{'}\left(c\right)dc\quad ; \quad dc=d\chi_{c}\]

L’opérateur tangent \(\tensFour{L}\) se déduit alors directement des expressions (4219) à (4222).

Pilotage PRED_ELAS

De manière générale, le pilotage consiste à décomposer le chargement extérieur (forces et déplacements imposés) en deux termes, l’un fixe (FIXE_CSTE, dans le vocabulaire de Code_Aster) et l’autre (FIXE_PILO) proportionnel au paramètre de pilotage \(\eta\) (ETA_PILOTAGE). On cherche alors à déterminer \(\eta\) de sorte à respecter une équation de pilotage qui dépend des degrés de liberté du problème et d’un paramètre qui contrôle la « vitesse » d’avancée le long de la branche de solutions.

Dans le cas des modèles d’endommagement, le choix du pilotage PRED_ELAS introduit dans [Lorentz-et-Badel-2004] est bien adapté : le paramètre d’avancée est un incrément d’endommagement cible noté \(\triangle a_{cible}\) et l’équation de pilotage assure qu’un point (d’intégration) de la structure atteint effectivement cet endommagement cible tandis que tous les autres points s’endommagent moins (ou pas du tout). De cette manière, on ajuste le chargement pour que la structure s’endommage progressivement.

L’algorithme de résolution proposé dans [Lorentz-et-Badel-2004] s’appuie sur une méthode de Newton pour la résolution des équations d’équilibre à l’échelle de la structure et une technique d’intervalles emboités pour l’équation de pilotage (qui requiert que le domaine de réversibilité soit convexe). La seule spécificité relative au modèle de comportement retenu se traduit par la résolution de l’équation scalaire suivante en \(\eta\), au niveau de chaque point d’intégration :

(4223)\[ g_{t}\left(\tensTwo{\varepsilon_{fixe}}+\eta\,\tensTwo{\varepsilon_{pilo}},a_{n}+\triangle a_{cible}\right)=0\]

Où \(\tensTwo{\varepsilon_{fixe}}\) et \(\tensTwo{\varepsilon_{pilo}}\) sont des déformations connues à ce stade et associées respectivement aux chargements fixes et pilotables, \(g_{t}\) la fonction seuil en traction et enfin \(a_{n}\) dénote à nouveau l’endommagement au début du pas de temps courant. Compte tenu de la forme de la fonction seuil, l’équation locale de pilotage s’écrit encore :

(4224)\[ \chi_{t}\left(\tensTwo{\varepsilon_{fixe}}+\eta\,\tensTwo{\varepsilon_{pilo}}\right)=e\left(a_{n}+\triangle a_{cible}\right)\]

[Lorentz-2025] détaille le calcul de la fonction de pilotage et des bornes associées.

Cas-tests

La loi ENDO_LOCA_TC est vérifiée et validée par les tests suivants :

Tableau 69 Cas tests de vérfication et validation

Code	Nom	lien doc
SSNV272	essai de traction - compression homogène et confinée sans viscosité	[V6.04.272]
SSNV273	essai de traction homogène et confinée avec viscosité	[V6.04.273]
SSNP180	Modélisation des essais SCIENCE avec ENDO_LOCA_TC	[V6.03.180]

Références

[Lorentz-2025] (1,2,3)

Lorentz, E., 2025. Construction, justification théorique et implantation numérique du modèle d’endommagement local de béton ENDO_LOCA_TC, Note EDF R&D 6125-1723-2024-03235-FR.

[Lorentz-2023]

Lorentz, E., 2023. Modélisation de la Rupture Quasi-fragile. Dans: J. Besson, F. Lebon & E. Lorentz, éds. Modélisation numérique en mécanique fortement non linéaire : contact et rupture. London: ISTE editions, pp. 189-273.

[Lee-et-al-2004]

Lee, S.-K., Song, Y.-C. & Han, S.-H., 2004. Biaxial behavior of plain concrete of nuclear containment building. Nuclear Engineering and Design, Volume 227, pp. 143-153.

[Lorentz-et-Badel-2004] (1,2,3)

Lorentz, E. & Badel, P., 2004. A new path-following constraint for strain-softening finite element simulations. International Journal for Numerical Methods in Engineering, Volume 60, pp. 499-526.

[fib-2010]

fib, 2010. Model Code for Concrete Stuctures. Fédération Internationale du Béton. éd. Berlin: Ernst & Sohn.

[Eurocode-2-2005]

Eurocode 2, 2005. Calcul des structures en béton. NF EN 1992-1-1, p. 32.

[Peterson-1981] (1,2)

Peterson, P.-E., 1981. Crack growth and development of fracture zones in plain concrete and similar materials, s.l.: Report TVBM—1006, Lund Institute of Technology.

[Lorentz-2017]

Lorentz, E., 2017. A nonlocal damage model for plain concrete consistent with cohesive fracture. International Journal of Fracture, Volume 207, pp. 123-159.

[Filiot-et-al-2017]

Filiot, A., Abbas, M. & Lorentz, E., 2017. Projet de développements pour l’amélioration de la prédiction dans code_aster. CR-T66-2017-131.