r5.03.33 Lois de comportement plastique etviscoplastique à écrouissage isotrope non linéaire VMIS_ISOT_NL et VISC_ISOT_NL

Résumé:

Ce document décrit le modèle de comportement élasto-visco-plastique à écrouissage isotrope non linéaire, disponible en version locale et à gradient de variables internes (GRAD_VARI ou GRAD_INCO). Une modélisation des grandes déformations est accessible via la cinématique GDEF_LOG. La fonction d’écrouissage cumule plusieurs termes (un terme affine, deux termes exponentiels, un terme en puissance, un plateau de Lüders) à même de couvrir la plupart des réponses non linéaires d’un matériau métallique; chacun peut facilement être annulé, ce qui offre une forme de modularité à la relation de comportement (par exemple: écrouissage linéaire ou écrouissage en loi puissance uniquement). Enfin, le comportement peut être indépendant du temps (VMIS_ISOT_NL) ou viscoplastique (VISC_ISOT_NL). Dans ce dernier cas, la plasticité évolue selon un modèle de Norton dans lequel intervient la fonction seuil.

Modèle continu

Équations de comportement : modèle plastique local

Il s’agit d’un modèle de plasticité à seuil de von Mises et écrouissage isotrope classique. La déformation (mécanique) \(\strain\) se décompose de manière additive en une partie plastique \(\strain^{p}\) et une partie élastique \(\strain^{e}\), cette dernière étant liée à la contrainte par la relation élastique suivante:

(2128)\[\sigma = K \trace (\strain^{e})\tensTwoUnit + 2\mu \,\dev(\strain^{e})\quad;\quad 3K=\frac{E}{1-2\nu }\quad;\quad 2\mu =\frac{E}{1+\nu }\]

où \(\trace\) et \(\dev\) désignent respectivement la trace et le déviateur d’un tenseur d’ordre deux.

On introduit maintenant la fonction seuil, où la fonction d’écrouissage \(R(\kappa)\) a déjà été introduite en (5001) :

(2129)\[F(\stress,\kappa ) = \stressCmp_{eq}-R(\kappa)\quad;\quad\stressCmp_{eq}=\sqrt{\frac{3}{2}\dev(\stress):\dev(\stress)}\]

On reste dans le cadre de la plasticité à écrouissage positif si bien qu’on impose que la fonction d’écrouissage \(R(\kappa)\) soit croissante.

L’évolution de la déformation plastique est gouvernée par l’équation d’écoulement. Lorsque la fonction seuil est dérivable par rapport à \(\stress\), la direction d’écoulement est normale à la surface seuil :

(2130)\[\text{Si } \stressCmp_{eq}\ne 0 \text{ alors } {\dot\strain}^{p}=\dot{\kappa} \tensTwo{N} \text{ avec } \tensTwo{N}=\frac{3}{2}\dev\frac{\sigma}{\stressCmp_{eq}}\]

Dans le cas singulier où \(\stressCmp_{eq}=0\), la notion de dérivée est généralisée via des notions d’analyse convexe (sous-gradient). La direction d’écoulement n’est assujettie qu’à la condition suivante:

(2131)\[\text{Si } \stressCmp_{eq}=0 \text{ alors } \sqrt{\frac{2}{3}} \Vert {\dot\strain}^{p} \Vert \le \dot{\kappa} \text{ et } \trace({\dot\strain}^{p})=0\]

Enfin, l’évolution de la variable d’écrouissage est fixée par la condition de cohérence:

(2132)\[\dot{\kappa} \ge 0 \quad ; \quad F(\sigma, \kappa)\le 0 \quad ; \quad \dot{\kappa}F(\sigma, \kappa)=0\]

Prise en compte de la viscosité

En présence de viscosité, l’écoulement plastique est différé. Plus exactement, la vitesse d’évolution de l’écrouissage n’est plus fixée par la condition de cohérence (5006) mais elle est dorénavant fonction de l’intensité de la fonction seuil suivant la loi de Norton:

(2133)\[\dot{\kappa} = {\left(\frac{\langle F(\sigma,\kappa ) \rangle}{K}\right)}^{n}\]

où les crochets de McAuley \(\langle \cdot \rangle\) désignent la partie positive. Les équations (5000) à (5005) restent inchangées. Lorsque \(K \to 0\) , on retrouve le modèle élastoplastique.

Prise en compte du gradient d’écrouissage

Une formulation non locale à gradient de variable interne est adaptée lorsque le gradient d’écrouissage \(\nabla \kappa\) devient important. Elle consiste à introduire un terme supplémentaire \({\Phi}^{\mathit{grad}}\) dans l’énergie libre qui reflète les interactions entre points matériels voisins, voir [R5.04.01]:

(2134)\[{\Phi}^{\mathit{grad}}(\nabla \kappa) = \frac{1}{2} \, c \, {\left(\nabla \kappa \right)}^{2}\]

En utilisant la formulation relaxée augmentée présentée dans [R5.04.01], l’impact sur la relation de comportement se traduit par une modification de la force thermodynamique associée à la variable d’écrouissage donc, en pratique, par une modification de la fonction d’écrouissage :

(2135)\[\tilde{R}(\kappa) = R(\kappa) + r\kappa - \varphi \text{ avec } \varphi = b\lambda + r \, a\]

où \(\lambda\) est le multiplicateur de Lagrange associé à la relaxation, \(r>0\) le coefficient d’augmentation et \(a\) l’interprétation de la variable d’écrouissage à l’échelle de la structure. Ces trois quantités sont des données pour ce qui concerne la relation de comportement. Au final, la fonction d’écrouissage est corrigée d’un terme affine, ce qui ne rajoute donc pas de complexité par rapport au traitement numérique du modèle. En revanche, le signe et l’amplitude de \(\varphi\) n’étant pas fixés, cela explique qu’il faille prendre en compte le cas d’un écoulement singulier pour lequel \(\stressCmp_{eq}=0\) .

Intégration numérique

Pour ce qui concerne l’intégration de la loi de comportement, c’est-à-dire le calcul des variables internes et des contraintes à histoire de déformation (mécanique) donnée, il n’y a plus lieu de distinguer s’il s’agit d’une formulation locale ou à gradient. En effet, l’impact de cette dernière se limite à modifier la fonction d’écrouissage selon (4999). On continuera ici à la noter \(R(\kappa)\) sans distinction.

En revanche, on distinguera les cas plastique et viscoplastique, même si on peut montrer que ce dernier se limite également à une correction de la fonction seuil dans les équations discrétisées.

Enfin, en présence d’un plateau de Lüders, on a choisi de simplifier la résolution en se ramenant à la situation sans plateau de Lüders. En effet, du fait de l’unicité de la solution du problème discrétisé, celle-ci se situe soit sur le plateau, soit au-delà. On commence par la chercher sur le plateau: on est alors ramené au cas de la plasticité parfaite (un cas particulier de fonction \(\widehat{R}\) ). Si la solution obtenue de cette manière s’avère au-delà du plateau, on ne la prend pas en compte (ce n’est pas la solution du problème) et on répète la résolution mais cette fois en cherchant la solution sur la partie \(\widehat{R}(\kappa)\) de la fonction d’écrouissage, c’est-à-dire sans se préoccuper du plateau de Lüders (par unicité, on est alors certain que cette nouvelle solution est effectivement au-delà du plateau). C’est pourquoi on se contente de décrire ci-dessous la méthode de résolution en l’absence de plateau de Lüders.

Équations discrétisées

La discrétisation en temps des équations de comportement s’appuie sur un schéma d’Euler implicite (même en viscoplasticité), c’est-à-dire que les différentes variables du problème sont exprimées à l’instant final du pas de temps considéré. On note \({Q}_{n}\) la valeur d’une quantité \(Q\) au début du pas de temps, \(\Delta Q\) son incrément pendant le pas de temps et (simplement) \(Q\) sa valeur en fin de pas de temps. L’état mécanique au début du pas de temps \((\strain_{n},\strain^{p}_{n},{\kappa}_{n})\) est supposé connu ainsi que l’incrément de déformation \(\Delta \strain\) (et donc également la déformation \(\strain\) ). Il s’agit alors de calculer les incréments de variables internes \(\Delta \strain^{p}_{n}\) et \(\Delta \kappa\), ainsi que la contrainte à la fin du pas de temps \(\stress\) .

La discrétisation de la relation contrainte – déformation (5002) s’écrit comme suit:

(2136)\[\stress = K \trace(\strain - {\strain^{p}}_{n}-\Delta \strain^{p})\tensTwoUnit + 2 \, \mu \,\dev(\strain-{\strain^{p}}_{n}-\Delta \strain^{p})\]

Quel que soit le régime d’écoulement, régulier (5004) ou singulier (5005), on constate que la trace de l’incrément de déformation plastique est nul, si bien qu’il en va de même pour la déformation plastique à tout instant (déformation plastique isochore). La relation (4948) se simplifie et on peut y faire apparaître une contrainte élastique (ou contrainte d’essai élastique), notée \({\stress}^{e}\) et qui est connue:

(2137)\[\stress = {\stress}^{e} - 2 \, \mu \Delta \strain^{p} \text{ avec } {\stress}^{e} = K \trace(\strain)\tensTwoUnit + 2 \, \mu \,\dev( \strain -\strain^{p}_{n} )\]

Dans le cas d’un écoulement régulier, la discrétisation temporelle de (5004) s’écrit:

(2138)\[\Delta \strain^{p} = \frac{3}{2} \, \frac{\Delta \kappa}{\stressCmp_{eq}} \, \dev({\stress}^{e}-2\mu \Delta \strain^{p})\]

Après quelques manipulations algébriques, on en déduit:

(2139)\[\stressCmp_{eq} = \stressCmp_{eq}^{e} - 3\mu \Delta \kappa\]

(2140)\[\Delta \strain^{p} = \frac{3}{2} \, \Delta \kappa \frac{\dev({\stress}^{e})}{\stressCmp_{eq}^{e}}\]

(2141)\[\stress = \frac{1}{3} \, \trace({\stress}^{e})\tensTwoUnit + \left ( \stressCmp_{eq}^{e} - 3\mu \Delta \kappa \right ) \frac{\dev({\stress}^{e})}{\stressCmp_{eq}^{e}}\]

Le régime régulier n’est valide que si \(\stressCmp_{eq}>0\), soit d’après (5010)

(2142)\[3 \mu \Delta \kappa > \stressCmp_{eq}^{e}\]

Dans le cas d’un écoulement singulier, \(\stressCmp_{eq}=0\) . Compte tenu de la définition de la contrainte équivalente de von Mises, c’est équivalent à \(\dev(\stress) = \tensTwoZero\). En particulier, la contrainte s’écrit simplement:

(2143)\[\stress = \frac{1}{3} \trace({\stress}^{e})\tensTwoUnit \text{ et } \stressCmp_{eq}=0\]

Par ailleurs, en injectant \(\dev(\stress)=\tensTwoZero\) dans (5008), on en déduit l’incrément de déformation plastique:

(2144)\[2 \mu \Delta \strain^{p} = \dev({\stress}^{e})\]

Et l’écoulement est assujetti à (5005), c’est-à-dire après discrétisation:

(2145)\[\sqrt{\frac{2}{3}} \Vert {\Delta \strain }^{p}\Vert \le \Delta \kappa\]

En reportant (4846) dans (4847), il vient:

(2146)\[\stressCmp_{eq}^{e} \le 3\mu \Delta \kappa\]

On peut noter que les conditions (5013) et (4848) permettent de discriminer le cas régulier du cas singulier. On introduit en particulier la valeur de séparation \({\kappa}^{s}\) entre les deux régimes:

(2147)\[{\kappa}^{s} = {\kappa}_{n} + \frac{\stressCmp_{eq}^{e}}{3\mu }\]

Résolution du problème discrétisé – Cas plastique sans viscosité

À ce stade, il ne reste plus qu’à déterminer l’incrément \(\Delta \kappa\) via la condition de cohérence (5006) qui s’écrit après discrétisation:

(2148)\[\Delta \kappa \ge 0\quad ; \quad F(\stress,\kappa )\le 0\quad ; \quad\Delta \kappa F(\stress,\kappa )=0\]

Pour cela, on commence par envisager l’hypothèse d’une solution élastique, c’est-à-dire \(\Delta \kappa =0\) et \(\Delta \strain^{p}=\tensTwoZero\) . Dans ce cas, \(\stress={\stress}^{e}\) et cette solution est valide d’après (4850) tant que \({F}^{e} = F({\stress}^{e},{\kappa}_{n}) \le 0\) .

Supposons maintenant que \({F}^{e}>0\) . Il s’agit de trouver \(\Delta \kappa > 0\) tel que \(F(\stress ,\kappa )=0\). Dans l’hypothèse d’un écoulement régulier (\({\kappa}_{n}<\kappa <{\kappa}^{s}\)), on peut injecter (5010) dans l’expression de la fonction seuil qui vaut alors:

(2149)\[F(\stress,\kappa) = \stressCmp_{eq}^{e} - 3\mu \Delta \kappa - R(\kappa) = \widehat{M}(\kappa)\]

Comme la fonction d’écrouissage \(R\) est croissante, la fonction \(\widehat{M}\) est strictement décroissante. Par ailleurs, \(\widehat{M}({\kappa}_{n})={F}^{e}>0\). Il existe une (unique) solution en régime d’écoulement régulier si et seulement si \(\widehat{M}({\kappa}^{s})<0\). Une méthode de Newton à bornes contrôlées permet alors de la déterminer.

Dans le cas d’un écoulement singulier (\({\kappa}^{s} \le \kappa\)), \(\stressCmp_{eq}=0\) et la condition de cohérence s’écrit simplement:

(2150)\[F(\stress,\kappa) = -R(\kappa) = 0\]

On peut remarquer que \(-R\) est une fonction (strictement) décroissante. En outre, Lorsque \(\kappa \to +\infty\) , \(-R\) est négative. En effet, dans le cas local, \(R(+\infty) \ge R(0)>0\). Et dans le cas non-local, \(r>0\) si bien que la fonction \(R\) tend vers l’infini. Il existe donc une (unique) solution en régime d’écoulement singulier si et seulement si \(-R({\kappa}^{s})=\widehat{M}({\kappa}^{s})\ge 0\). Là encore, une méthode de Newton permet de déterminer cette solution.

Finalement, il existe une solution unique au problème d’intégration de la loi de comportement. Le régime dépend des valeurs \(\widehat{M}^{e}=\widehat{M}({\kappa}_{n})\) et \(\widehat{M}^{s}=\widehat{M}({\kappa}^{s})\), où \(\widehat{M}^{e}\ge \widehat{M}^{s}\) :

• si \(0 \ge \widehat{M}^{e} \ge \widehat{M}^{s}\) : solution élastique ;

• si \(\widehat{M}^{e} \ge 0 \ge \widehat{M}^{s}\) : solution plastique à écoulement régulier ;

• si \(\widehat{M}^{e} \ge \widehat{M}^{s} \ge 0\) : solution plastique à écoulement singulier.

Résolution du problème discrétisé : cas viscoplastique

En présence de viscosité, la condition de cohérence (5006) est remplacée par la loi d’évolution de Norton (5007). Une fois discrétisée, cette dernière s’écrit:

(2151)\[\frac{\Delta \kappa}{\Delta t} = \left( \frac{\langle F(\stress,\kappa) \rangle }{K} \right)^{n}\]

Et en inversant cette relation:

(2152)\[{C}_{\Delta t}{\Delta \kappa }^{q} = \langle F(\stress ,\kappa ) \rangle \quad ; \quad{C}_{\Delta t} = \frac{K}{{\Delta t}^{1/n}} \quad ; \quad q=\frac{1}{n}\]

Elle peut encore s’exprimer de manière équivalente sous la forme d’une relation de Kuhn et Tucker:

(2153)\[\Delta \kappa \ge 0\quad ; \quad{F}_{v}(\stress ,\kappa )\le 0\quad ; \quad{F}_{v}(\stress ,\kappa )\,\Delta \kappa =0\]

où

(2154)\[{F}_{v}(\stress,\kappa) = F(\stress,\kappa) - {C}_{\Delta t}{\Delta \kappa }^{q} = \stressCmp_{eq} - {R}_{v}(\kappa)\]

avec

\[{R}_{v}(\kappa) = R(\kappa) + {C}_{\Delta t}{(\kappa -{\kappa}_{n})}^{q}\]

Après discrétisation, la prise en compte de la viscosité se réduit à l’ajout d’un terme (positif et croissant) dans la fonction d’écrouissage.

On pourrait être tenté d’appliquer la même méthode de résolution que précédemment mais, malheureusement, ce nouveau terme n’est pas dérivable en \({\kappa}_{n}\) car \(q<1\) en général (\(n>1\)).C’est pourquoi on propose de commencer à chercher la solution du problème sans viscosité (de manière grossière) et on bascule ensuite sur la solution avec viscosité par une méthode de Newton moyennant le changement de variable \(\Delta \kappa \to {C}_{\Delta t}\,{\Delta \kappa }^{q}\), éliminant ainsi le problème de non dérivabilité (et les fortes pentes qui sont pénalisantes pour une méthode de Newton).

Matrice tangente

Dans le cas local, il s’agit de déterminer le tenseur d’ordre quatre \({H}_{\stressCmp \strainCmp} = \delta \stress / \delta \strain\). Dans le cas non local, il faut également construire les tenseurs d’ordre deux \({H}_{\stressCmp \varphi} =\delta \stress /\delta \varphi\) et \({H}_{\kappa \strainCmp}=\delta \kappa /\delta \strain\) ainsi que le scalaire \({H}_{\kappa \varphi }=\delta \kappa /\delta \varphi\) .

L’opérateur tangent dépend du régime d’écoulement. Dans le cas élastique, il se réduit à:

(2155)\[\delta \stress = K \, \trace(\delta \strain) \tensTwoUnit + 2\mu \dev( \delta \strain ) \text{ et } \delta \kappa =0\]

Dans le cas du régime singulier, l’expression est là aussi très simple:

(2156)\[\delta \stress = K \, \trace(\delta \strain) \tensTwoUnit \text{ et } \delta \kappa = {\left(\frac{\partial {\tilde{R}_{v}}}{\partial \kappa }\right)}^{-1}\delta \varphi\]

où la fonction d’écrouissage contient le cas échéant les termes liés au non local et à la viscosité (d’où le tilde et l’indice pour le rappeler explicitement).

Le régime d’écoulement régulier conduit à un opérateur tangent plus complexe. L’usage d’un critère de von Mises conduit aux dérivées suivantes:

(2157)\[\delta \stressCmp_{eq} = 3 \mu N:\delta \strain\]

et

(2158)\[\delta N = \frac{1}{\stressCmp_{eq}} \left[2\mu \dev(\delta \strain)-\delta \stressCmp_{eq}N\right]\]

En appliquant ensuite des dérivées de fonctions composées, on obtient:

(2159)\[\delta \stress = K \trace(\delta \strain )\tensTwoUnit + 2\mu \left[ 1 - \frac{3\mu \Delta \kappa}{\stressCmp_{eq}} \right] \, \dev(\delta \strain) + \frac{9{\mu}^{2}\Delta \kappa}{\stressCmp_{eq}} \tensTwo{N} \otimes \tensTwo{N} : \delta \strain - 3\mu \tensTwo{N} \delta \kappa\]

(2160)\[\delta \kappa = \left[ 3\mu + \frac{\partial {\tilde{R}_{v}}}{\partial \kappa} \right]^{-1} \, \left[ \delta \stressCmp_{eq} + \delta \varphi \right]\]

En phase de prédiction (\(\strain =\strain_{n}\)), deux spécificités sont à noter :

• En présence de viscosité, \(\delta \kappa=0\) car la dérivée \(\frac{\partial {\tilde{R}_{v}}}{\partial \kappa}\) est infinie.

• Les opérateurs tangents ne sont pas continus aux transitions entre régimes (\(\widehat{M}^{e}\approx 0\) et \(\widehat{M}^{s}\approx 0\) ) et dans ce cas on opte pour le domaine correspondant au régime du pas précédent. En revanche, si l’état mécanique du début du pas de temps n’est pas au voisinage d’une transition pour, alors le régime est fixé comme pour l’intégration de la loi.

Pilotage PRED_ELAS

Dans le cadre de la relation de comportement VMIS_ISOT_NL, le pilotage du chargement PRED_ELAS consiste à déterminer l’intensité du chargement piloté \(\eta\) tel que l’incrément de variable d’écrouissage \(\Delta \kappa\) soit égal à une valeur cible \(\Delta p\) fixée par le pas de temps courant \(\Delta t\) et la valeur du coefficient de pilotage COEF_MULT \(C\) via \(\Delta p = \frac{\Delta t}{C} \times \text{5 %}\). On pourra se reporter à [R5.03.80] pour plus de détails sur le principe de la méthode de pilotage PRED_ELAS. On rappelle en particulier que la notion de pilotage n’a de sens qu’en l’absence d’un comportement visqueux. On se limite donc ici au cas sans viscosité. En outre, seul le pilotage pour le modèle local (i.e. sans prise en compte du gradient de plasticité) est traité à l’heure actuelle.

A l’échelle algorithmique de la loi de comportement, le pilotage introduit une déformation affine par rapport à \(\eta\), à savoir \(\strain=\tensTwo{\strain^0}+\eta\,{\strain^1}\) où \({\strain^0}\) et \({\strain^1}\) sont des quantités (tensorielles) connues à ce stade. Compte tenu de (4851), il vient alors :

(2161)\[\stressCmp^e_{eq}(\eta) = 3\mu\,\Delta p + R(\kappa_n + \Delta p)\]

où le membre de droite est une quantité connue. Comme \(\strain\) est affine en \(\eta\), il en va de même de \({\stress^e}\) d’après (5008), si bien que résoudre (2161) se ramène à la recherche des racines d’un polynôme du second degré. On détermine ainsi les solutions \(\eta\) de l’équation de pilotage (au nombre de 0, 1 ou 2).

Bibliographie

Lorentz E., Andrieux S. (1999) A variational formulation for nonlocal damage models. Int. J. Plas. 15, 119-138.

Zhang Y., Lorentz E., Besson J. (2018) Ductile damage modelling with locking-free regularis2. ed GTN model. Int. J. Numer. Methods Eng. 113, 1871-1903.