Plasticit 3D
Formalisme gnral,
dtail des modles ECL, EI, ECNL
et bases des schmas dintgration
Marc Franois
1
Les spcificits de
la plasticit 3D
Linda Besemer
1.1
La dformation plastique
Les expriences montrent que la transformation plastique
est un phnomne isochore ( volume constant) :
V
p
trace(" ) =
=0
V
Le tenseur des dformations plastiques
est un dviateur.
1
pD
p
p
p
"
="
trace" I = "
3
ou, de manire quivalente
H
p
P :" = 0
PD : "p = "p
De plus, on a, comme en 1D, la
[Link]
partition des dformations:
(si on soccupe de ces eets).
1.2
La plasticit vit dans lespace dviatorique
Les expriences de P.W. Bridgman (Nobel, 1946) montrent que
seule la partie dviatorique du tenseur des contraintes a un rle
en plasticit. La pression hydrostatique :
na aucune influence (pour les pressions modres). (revoir le
cours 1D).
Donc toute la plasticit ne dpendra que de , le dviateur de
. On rappelle :
Donc, les causes et les eetssont dviatoriques: les
autres variables de la plasticitle seront aussi (ou scalaires).
1.3
La loi dlasticit en lastoplasticit
Quel que soit le modle, la relation dlasticit sera du type:
ou, depuis la partition des dformations:
p
= C : (" " )
Ce qui revient dire que touts les potentiels dtat
commenceront comme:
car (TPI) :
@
@
@
=
= e =
@"
@"
@(" "p )
La loi dlasticit peut scrire sous les trois formes:
Young:
ELASY
Lam:
ELASL
Kelvin:
ELASK
Bien que moins classique, cette dernire forme est la plus
adapte cause de la sparation H/D et sera utilise par dfaut
dans ce cours. Et donc, lnergie libre de Gibbs commence
toujours comme:
Elle sinverse facilement :
KELAS
1.4
Le critre dlasticit
En 1D on avait des critres dlasticit scalaires :
Saint Venant
ECL, ECNL
EI
En 3D, on gnralise ces critres en bornant la norme
euclidienne naturelle du dviateur des contraintes. Pour le
modle trivial de Saint Venant, cela correspond au critre de
Von Mises(voir cours dalgbre tensoriel).
Von Mises :
8
Pour lcrouissage cinmatique (linaire ou non) on introduit
une variable tensorielle X qui dcrit le dplacement du centre
du domaine lastique dans lespace des dviateurs:
Esp. dv. (5D)
ECL, ECNL :
Pour lcrouissage isotrope, on conserve la variable de surcrouisage R qui dcrit laugmentation du domaine dlasticit.
Cela correspond une dilatation dans lespace des dviateurs :
EI :
Esp. dv. (5D)
Quelques mots sur Von Mises et Henki et Huber !
On trouve le critre en gnral crit comme suit :
[Link] Mises
Mais le sens physique de J2 est moins clair que celui de la
norme
En base canonique, on a:
En base de tenseurs:
Le critre de Von Mises en fonction des contraintes
principales( viter car cest trs couteux en CPU):
10
Le coecient sert caler sur lessai 1D:
Si lon regarde les essais trs finement le critre de Tresca colle
mieux mais il nest pas drivable
S!
[0 0 0 0 0 1]
[0 0 0 0 0 1]
Tresca
8!
8"
2/ 6
1/ 6
1/ 6 0 0 0
Essais sur aluminium Au4G-T4
Thse M. Rousset LMT 1985.
Modle M. Franois IJP 1995.
S"
2/ 6
1/ 6
1/ 6 0 0 0
11
2
Ltat
[Link]
12
2.1
Les variables dtat
Rappels de TPI :
En plasticit on a forcment la premire variable dtat qui est
la dformation plastique:
Le autres modles introduisent dautres variables dtat et leurs
forces thermodynamiques associes. Le choix nest pas unique:
dans ce cours on a retenu:
Vk
Ak
ECL
"p
D
Vk
Ak
EI
"p
D
r
R
ECNL
Vk "p
D
Ak
X
13
2.1.1
crouissage Cinmatique Linaire
Pour lECL on introduit une force thermodynamique
associe avec la dformation plastique
Vk
Ak
"
D
On se donne lnergie libre de Gibbs suivante:
GIBBS
On obtient la contrainte par drivation:
laide de :
on retrouve bien ELASK
14
La contrainte dcrouissage est obtenue par la drivation
suivante :
ETAT
Le 2/3 sert retrouver le X scalaire de ce modle un 1D
15
T.D. de cours :
On donne :
= 0
0
0 0
0 0
0 0
0
"p /2
0
"p
p
" = 0
0
0
0
p
" /2
Calculer et
en dduire et et
retrouver les quations 1D
+
Calculer lnergie libre de Gibbs
et identifier les aires correspondantes
sur le graphe 1D ci-contre. Que
reprsente lnergie libre ?
O
+
16
2.1.2
crouissage isotrope
Dans ce cas on introduit la force thermodynamique r associe
avec la force R, le (sur-crouissage) (voir T8).
Vk
Ak
"p
D
r
R
On se donne lnergie libre de Gibbs suivante :
GIBBS
Pour la contrainte on retrouve ELASK
Et pour le sur-crouissage:
On verra ensuite que
ou
et que, comme en 1D:
ETAT
17
2.1.3
crouissage cinmatique non linaire
Dans ce modle, il faut introduire une nouvelle variable
tensorielle associe X, en plus de la dformation plastique:
Vk
Ak
"
GIBBS
Lnergie libre de Gibbs est
Et la loi dlasticit est toujours :
ELASK
Mais bien sr on dfinit ici
ETAT
mais le rapport avec le 1D est encore lointain
18
3
La surface de charge
Mercedes-Benz
19
Cest une fonction convexe des forces thermodynamiques.
Elle est ngative dans le domaine dlasticit et nulle lors de la
transformation plastique. Elle ne peut pas tre positive
(domaine interdit). En relation avec ce qui prcde on a :
ECL, ECNL:
SDC
SDC
EI
SDC
20
4
La normalit
21
En 1D on na quune direction dcoulement possible ; seul le
signe de la dformation plastique est dfinir. Par exemple en
ECL:
En 3D, il faut connaitre la direction du tenseur :
sera-t-il colinaire ? ?
Cest lobjet de ce chapitre. Pour lintensit, ce sera au chapitre
suivant.
Deux approches direntes conduisent la rgle de normalit.
22
4.1
Les essais et la normalit la surface de
charge
Des campagnes dessais trs fines ont permis didentifier que
la direction dcoulement plastique, tait normale la surface
de charge (mme si celle-ci tait loin de celle de Von Mises):
Rousset, LMT-Cachan, 1985
23
Cete normalit a lieu dans lespace des contraintes. On en
dduit lquation de normalit de lcoulement :
@f
d" =
d
@
p
NORM
dans laquelle f est la surface seuil, cest dire que llasticit est
dfinie par f0 et un multiplicateur de Lagrange positif qui
sera dfini plus [Link] nos modles nous avons :
EI :
Esp. dv. (5D)
24
ECL, ECNL :
Esp. dv. (5D)
O n correspond aussi la normale (tensorielle et norme) la
surface de charge.
25
Une consquence macroscopique de la rgle de normalit :
aprs une traction, lapplication dune torsion (contrainte de
cisaillement) commence par allonger la pice !
d"
2
3
L coulement plastique a lieu dans la direction de la normale
la surface de charge ; elle ne change qu long terme, pas
instantanment.
26
Cette rgle ne sapplique pas sur tous les matriaux : par
exemple, pour un matriau granulaire (sable), on a :
Mais inapproprie pour les matriaux granulaires
(sables):
f( ) = 0
d"
12
p
2 11 / 6
Critre de Coulomb
coulement volume (approximativement) constant:
nest pas suivant la normale au critre
27
Cours 1D
Remarque : on peut vrifier que cette condition est cohrente
avec les critures du cas 1D.
Ex. ECL
1D
Depuis la condition prcdente :
+
@f
d" =
d
@
p
28
4.2
Condition dvolution plastique
Comme on la vu en 1D, tre f=0 ne sut pas avoir une
transformation plastique. Il faut que le chargement tende
sortir du domaine dlasticit.
1D
3D
Esp. dv. (5D)
+
En 3D cette condition scrit :
et
ou encore :
29
4.3
La TPI et le principe de dissipation
maximale
La TPI impose la positivit de la dissipation intrinsque:
@
qi1 =
dVk > 0
@Vk
Mais on na pas la direction des tenseurs
Les rgles dvolution montres dans le cours de TPI doivent
tre adaptes au cas de la plasticit.
On postule lexistence dun pseudo-potentiel de dissipation
donnant
les relations forces-flux:
'(V k , ~q /T )
!
@'
@'
Ak =
gradT =
@~q /T
@ V k
30
'
0, . . .) = 0
est
suppos convexe et nulle zro '(0,
positive et continue
Pour une telle fonction, on a la proprit:
~x
O
!
@f
= grad(~x)
@~x
@f
.~x > 0
@~x
f (~x) = k > 0
soit, ici: ~x $ V k
@'
Vk = Ak V k = qi1 > 0
@ V k
ce qui revient garantir la positivit de la dissipation.
31
Mais en plasticit la dissipation nintervient que si le critre est
atteint. Or, il est exprim en fonction des forces
thermodynamiques, pas des flux.
'(V k , ~q /T )
On doit alors passer du pseudo-potentiel en flux
un autre en forces
!
F (Ak , gradT ) 6 0
Ceci est faite grace la transforme de Legendre-Fenchel
(opration topologique).
!
!
!
~q .gradT
~q
F (Ak , gradT ) = sup Ak Vk
' Vk ,
T
T
V k ,~
q /T
32
Transforme de Legendre-Fenchel en 1D (on ne considre que
"p
A1 et V1 soit et ).
F ( ) = sup
: "p
' "p
"p
33
On peut dmontrer que si F est drivable, la proprit de
normalit est conserve et:
@F
V k =
@Ak
~q
@F
=
!
T
@ gradT
Il faut, pour que le second principe soit vrifi, que F soit
convexe, positive et nulle lorigine.
La seconde quation (thermique) correspond au terme de
dissipation thermique(cours TPI) :
!
~q .grad(T )
qi2 =
dt
T
34
En plasticit le temps physique nintervient pas (voir cours). On
remplace
@F
dVk =
dt
@Ak
par
dVk
@F
d
@Ak
NORG
d
o est
un multiplicateur de Lagrange (multiplicateur
plastique) qui sera dtermin ultrieurement. Il est commun
toutes les Vk.
Cette quation est nomme normalit gnralise.
35
4.4
Lien ente les deux approches
La normalit vue au 3.2 impose que le pseudo-potentiel F est
li la surface de charge f par certaines relations :
Pour respecter la normalit il faut que, comme V1 = "p
@f
@F
p
d" =
d =
d
@
@A1
En rgle gnrale, f et F sont identiques par rapport :
A1 =
Les conditions ci-dessus sont respectes si (ex.
EI) ou
A1
si =+
(ex/ ECL).
36
3.5
Classes de modles standard
Un matriau obissant aux rgles prcdentes est dit matriau
standard gnralis (Halphen et Coq Son). De plus:
Matriau standard gnralis associ
F = f + cste
Matriau standard gnralis non associ: f et F ne sont
identique quenpour respecter la normalit.
F, = f,
37
3.6
cr. Cinmatique Linaire
Le modle ECL est associ donc :
r
3 D
D
F(
X) =
k
Xk
2
Vk
Ak
"
D
Lcoulement plastique est donn en direction par la
normalit :
@F
p
d" =
d
@( D X)
Cest un calcul que lon a dj fait:
EVO1
On dfinit la plasticit
cumule p :
EVO2
38
3.7
crouissage Isotrope
Vk
Ak
"p
D
r
R
Le modle EI est associ donc :
Pour la premire V.E. on retrouve automatiquement la normalit
et la plasticit cumule:
EVO1
EVO2
Lautre volution est donne par la drivation / -R
39
p
V
"
3.8
cr. Cinmatique Non Linaire k D
Ak
X
Le modle ECNL est non associ : la surface seuil est ajoute
un terme complmentaire qui ne dpend pas de afin de
respecter la normalit la surface de charge.
r
3 D
3
D
F(
X) =
k
Xk +
X:X
2
4C
Lcoulement plastique est similaire lECL :
r
D
X
3
p
n=
d" =
nd
k D Xk
2
EVO1
EVO2
la dformation plastique
cumule :
40
On drive F par rapport la seconde force thermo.:
et depuis ltat (2.1.3)
r
3
En tenant compte de :on
trouve la relation
d" =
nd
2
suivante, similaire celle du 1D.
p
2
p
dX = Cd"
3
Xdp EVO3
On remarque que lon pratiquement se passer de la VEqui
est ncessaire pour la TPI mais directement associe X.
41
3.9
La dformation plastique cumule
Tous les modles ont en commun la valeur de la dformation
plastique cumule :
Esp. dv. (5D)
On a donc :
ce qui, au coecient prs, correspond
la longueur de la courbe dcrite par
dans lespace des dviateurs.
En plasticit pure, p ne sert au lEI. Mais il sert aussi en
lastoplasticit endommageable, etc
42
5
La consistance
43
5.1
Lquation de consistance
On connait maintenant la direction de lcoulement plastique
(normalement la surface de charge) et des autres variables
internes, grce au pseudo-potentiel de dissipation (qui est trs
proche de la surface seuil).
dVk
@F
d
@Ak
Mais il reste dfinir le multiplicateur plastique
Pour cela on utilise lquation de consistance, qui impose
lors de la transformation plastique (ce qui revient dire que lon
emmne la surface de charge):
@f
df = 0
dAk = 0
@AK
44
5.2
Consistance ECL
Vk
Ak
"
D
Depuis SDC on a:
Depuis EVO2
CONS
Cette quation reprsente lentrainement de la surface de
charge, suivant n.
Esp. dv. (5D)
45
5.3
Consistance EI
Vk
Ak
"p
D
r
R
Depuis SDC on a:
Depuis EVO2
CONS
Cette quation reprsente lentrainement de la surface de
charge, suivant n.
Esp. dv. (5D)
46
5.4
Consistance ECNL
Vk
Ak
"p
D
Depuis SDC on a:
Depuis EVO2
CONS
Cette quation reprsente lentrainement de la surface de
charge, suivant n.
Esp. dv. (5D)
47
6
Lvolution
48
6.1
Le pilotage
Dun point de vue pratique il faut connaitre quelque chose pour
faire le calcul ! Selon les cas on peut avoir:
Le pilotage en contraintes :
on connait
on cherche
Cest utile pour les calculs analytiques et souvent plus simple en
plasticit.
Le pilotage en dformation :
on connait
on cherche
Cest plus dlicat en plasticit mais cest le pilotage utilis dans
les codes de calcul : le code envoie une dformation et rcuper
la contrainte. Ex. UMAT dAbqus.
49
6.2
Le calcul explicite
Cest la programmation la plus simple : les VE sont prises
linstant n et on calcule les incrments de n vers n+1.
On ne ractualise pas les donnes: les pentes (ou gradients)
sont pris ltat prcdent.
On peut commettre beaucoup derreur : il faut prendre des pas
de calculs assez petits.
Pente:
cart
50
6.3
LECL pilot en contraintes
Vk
"
D
Ak
X
En plus de ltat actuel, lincrment en contraintesest
connu. Il nous faut calculer le multiplicateur plastique en
fonction de celui-ci.
ETAT et EVO1:
CONS :
CONSS
51
Lalgorithme de lECL pilot en contraintes
Donnes :
lasticit?
SDC
NON
OUI
Plasticit
lasticit
n=
X
Xk
EVO2
CONSS
EVO1
ETAT
KELAS
52
function [S,X,Ep] = ECNL3D(E);
% Modele crouissage cinmatique non linaire pilote en deformations
% d'Armstrong-Fredericks
% Mthode explicite pure
% M. Franois - Univ de Nantes - octobre 2011
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% Inits
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Constantes matriaux
global C gamma Y nu sy
% En Kelvin
K = Y/(3*(1-2*nu));
mu = Y/(2*(1+nu));
% Projecteurs
PH = zeros(6);
PH(1:3,1:3) = 1/3;
PD = eye(6) - PH;
% Init des rsultats
S = zeros(size(E)); % Contrainte
X = zeros(size(E)); % Backstress
Ep= zeros(size(E)); % Deformation plastique
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% Programme principal
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Methode explicite pure : l'etat est pris au pas n-1
% Corps du programme
for k=2:size(E,2)
% Prdicteur lastique
X(:,k) = X(:,k-1);
Ep(:,k) = Ep(:,k-1);
S(:,k) = (2*mu*PD+3*K*PH)*(E(:,k)-Ep(:,k));
% si plasticit : correcteur plastique
if sqrt(3/2)*norm( PD*(S(:,k)-X(:,k)) )-sy>0
% Normale a la surface de charge
n = PD*(S(:,k-1)-X(:,k-1));
n = n/norm(n);
% Multiplicateur plastique
dp = sqrt(3/2)*2*mu*(n'*(E(:,k)-E(:,k-1)))/...
(3*mu + C - sqrt(3/2)*gamma*n'*X(:,k-1));
% Actualisation des variables internes
dEp = sqrt(3/2)*n*dp;
Ep(:,k) = Ep(:,k-1) + dEp;
X(:,k) = X(:,k-1) + 2*C*dEp/3 - gamma*X(:,k-1)*dp;
% Actualisation de la contrainte
S(:,k) = (2*mu*PD+3*K*PH)*(E(:,k)-Ep(:,k));
end
53
6.4
LECL pilot en dformations
Il nous faut calculer le multiplicateur plastique en fonction
de lincrment de dformation connu.
Depuis ELASK :
n dviatorique
CONSS et EVO1 :
CONSE
54
Lalgorithme de lECL pilot en dformations
Donnes :
lasticit?
ELASK
ETAT
SDC
OUI
NON
Plasticit
lasticit
n=
X
Xk
EVO2
CONSE
EVO1
ELASK
6.5
LEI pilot en contraintes
Vk
Ak
"p
D
r
R
Il nous faut en fonction de.
CONS et ETAT:
CONSS
56
Lalgorithme de lEI pilot en contraintes
Donnes :
lasticit?
NON
OUI
SDC
Plasticit
lasticit
EVO2
CONSS
EVO1
KELAS
57
6.6
LEI pilot en dformations
Il nous faut calculer le multiplicateur plastique en fonction
de lincrment de dformation connu.
Depuis ELASK :
n dviatorique
CONSS
et EVO1 :
CONSE
58
Lalgorithme de lEI pilot en dformations
Donnes :
lasticit?
ELASK
SDC
NON
OUI
Plasticit
lasticit
EVO2
CONSE
EVO1
ELASK
56
59
6.7
LECNL pilot en contraintes Vk
Ak
Le multiplicateur plastique en fonction de:
"p
D
EVO1 et EVO3:
CONS :
CONSS
60
Lalgorithme de lECNL pilot en contraintes
Donnes :
lasticit?
SDC
OUI
NON
Plasticit
lasticit
n=
X
Xk
EVO2
CONSE
EVO1
2
p
dX = Cd"
3
Xdp
EVO3
KELAS
61
6.8
LECNL pilot en dformation Vk
Ak
Le multiplicateur plastique en fonction de:
"p
D
Depuis ELASK :
n dviatorique et
EVO1
et
CONSS
do :
CONSE
62
Lalgorithme de lECNL pilot en dformations
Donnes :
lasticit?
ELASK
SDC
OUI
NON
Plasticit
lasticit
n=
X
Xk
EVO2
CONSE
EVO1
2
p
dX = Cd"
3
Xdp EVO3
ELASK
63
7
Complments
[Link]
64
7.1
Module tangent
Le module tangent est dfini comme la drive droite suivante:
et son inverse:
qui intgre donc lvolution de toutes les variables et qui dpend
de lvolution.
Pour mmoire, la drive partielle est calcule en bloquant les
autres variables et :
Cest un tenseur du quatrime ordre, de mme structure que le
tenseur dlasticit. Il est trs souvent demand dans les codes
car il permet un calcul rapide.
65
ELASK
EVO1
(commune tous les modles ECL, EI, ECNL).
On rsume CONSE par :
O h est un modle dcrouissage qui dpend du modle :
ECLEIECNL
il vient :
66
7.2
Module tangent inverse
Il sagit dinverser un tenseur du quatrime ordre.
Heureusement on peut utiliser la proprit des projecteurs
. Il faut pour cela transformer lexpression prcdente
en projecteurs orthogonaux pondrs:
67
Autre preuve pour lexpression de H depuis H-1 (ALV)
depuis H-1 :
Ouf !
68
7.3
Cas gnrique
Si lon regroupe toutes les quations donnes par la TPI, la
normalit, la consistance on peut crire un framework complet
pour la plasticit. Universel pour tous les modles standards
gnraliss (associs ou non) peut tre utile pour des modles
plus labors.
69
Vk avec
Variables : dtat,
V1 = "
Ak
associes aux forces.
(", Vk )
nergie libre (potentiel
dtat) et quations dtat:
Ak
p
@
=
@"
@
=
@Vk
Pseudo-potentiel de dissipation
F (Ak ) convexe et nul lorigine.
do relations dcoulement dVk = @F/@Ak d
@f
@F
Seuilconvexe,
avec condition
f (Ak ) 6 0
=
@
@A1
@f
Condition dcoulement:&
f =0
:d >0
@
Consistance df
donne
(avec)
(Ak ) = 0
H = d /d"
d
@f @ 2
@Ak @Vk @" : d"
@f
@2
@F
@Ak @Vk @Vl @Al
H=C
@2
@f @ 2
@F
@"@Vk @Ak @Al @"@Vl
@f
@2
@F
@Ak @Vk @Vl @Al
p
qi1 = Ak dVk
Dissipationplast.
cumule dp = 2/3kd"p k
70
7.4
Bilan nergtique
On considre lquation suivante (cours de TPI, T15)
Pour lentropie, on a :
dans les modles prsents, on na pas introduit de dpendance
en temprature
On a aussi :lnergie thermique dissipe.
Et
lnergie fournie de lextrieur
se rpartit en une partie dissipe et une stocke
71
Bilan nergtique ECL :
EVO1 et EVO2
SDC
En direntiant GIBBS
on retrouve bien, en 1D, le jeu daires
tudi au cours de Plasticit 1D.
Les deux nergies sont rcuprables ; lnergie lastique lest
plus directement.
72
Bilan nergtique EI :
EVO1 et EVO2
SDC
En direntiant GIBBS
on retrouve bien, en 1D, le jeu daires
tudi au cours de Plasticit 1D.
Lnergie stocke Rdp, strictement croissante, est irrcuprable
tandis que lnergie lastique est rcuprable.
73
p
V
"
k
Bilan nergtique pour lECNL :
D
Ak
En utilisant ETAT, EVO1, EVO2, EVO3 et SDC
Pour GIBBS
Pas dinterprtation simple sous forme daire Les nergies
sont rcuprables, lnergie lastique lest directement.
74
8
Intgration numrique (bases)
75
chaque tape gobale du calcul E.F. le code propose un champ
de dplacement
. Ltape locale consiste dterminer,
partir du champ de dformation qui en rsulte, de calculer
les contraintes associes
, les variables internes
et les forces associes .
Soit en rsum :
Donnnes :
Inconnues :
76
8.1
Le calcul explicite
La mthode la plus simple programmer.
Dans un calcul explicite, les gradients sont calculs au point n
connu.
cart
La contrainte est donne directement par le module tangent
(c.f. cours de plasticit 3D):
Pour tous les modles classiques (ECL, EI, ECNL) on a la
forme suivante pour H (T. 66) :
HNN
77
On commence par tester une volution lastique, ce qui
correspond :
On regarde si cette contrainte est eectivement dans le
domaine dlasticit (en conservant les VI prcdentes) :
Esp. dv. (5D)
si
alors on a eectivement une
transformation plastique et donc :
78
On considre maintenant le cas oppos o lon a plasticit :
Il sagit donc dune transformation plastique. On doit valuer le
module tangent complet au pas n.
ex. ECL
T. 66
Esp. dv. (5D)
Ensuite on calcule la contrainte dans le cas
de lvolution plastique, qui apparait comme
une correction de la prcdente :
(depuis HNN)
Correcteur plastique
79
Les variables internes sont ractualises. On calcule le
multiplicateur plastique grce lquation de consistance en
dformation (CONSE), toujours avec les valeurs au point n de la
forme :
ex. ECL
T. 54
Ensuite on peut ractualiser les variables internes grce la
normalit gnralise NORG (T.35), en pratique les quations
EVOi que lon en a dduit :
ex. ECL
@F
T. 38
d
NORG dVk =
EVO1
@Ak
80
Si besoin on calcule les forces thermodynamiques depuis la
valeur des variables dtat avec les relations dtat:
ex. ECL
T. 15
ETAT
On obtient donc les
et lon peut passer ltape suivante
81
Remarque 1 :
En gnral la contrainte calcule en explicite ne vrifie en
gnral pas exactement f=0.
Esp. dv. (5D)
Remarque 2 :
La mthode explicite ncessite donc de raliser de petits
incrments. En gnral on doit raner le pas de calcul du code
E.F. de manire conserver une prcision acceptable et donc
raliser de nombreuses fois le problme global.
Cest donc simple mais coteux en temps CPU.
82
8.2
Mthode du retour radial (Ortiz et Simo)
La mthode explicite prcdente ne garantit pas exactement :
Ce qui peut entraner des instabilits dans le calcul. La
mthode de retour radial (Ortiz et Simo, 63) le garantit.
Le prdicteur lastique et le cas de la transformation lastique
sont calculs comme prcdemment.
Pour le cas de la transformation plastique
on considre, depuis cet tat *,
une transformation dformation
globale nulle (puisque est (dj) impos)
qui ramne exactement f=0 selon un
processus de Newton.
83
Les d indiqus ici sont donc relatifs la sous itrations de
Newton * ** qui est donc intermdiaire entre le pas n et le
pas n+1.
Pour cette transformation on a donc :
et, depuis EVO1 (commune ECL, EI, ECNL) et EVO2
r
3
p
d" =
nd
2
retour radial : correction suivant n
84
On applique un algorithme de Newton : on suppose la solution
f=0 atteinte au sous-pas suivant (ces sous-pas sont inclus dans la
transformation).
connu (page prcdente).
connu.
Il nous manque juste :
85
et est donn par la normalit :
au bilan on a :
soit le multiplicateur plastique :
laide de ce multiplicateur plastique on ractualise les
variables internes (comme au T5 du prsent cours).
et
on actualise les et les
@F
dVk =
d
@Ak
Et on itre jusqu ce qu convergence avec un critre en
86
Remarque 1 :
On a dj disposition la pluspart des expressions requises par
Exemple : pour lECL :
Pour llasticit :
D
X
n
=
EVO2 :
k D Xk
Vk
Ak
"p
D
depuis SDC et T38
Depuis
Au final (bon exercice !)
87
Remarque 2:
Cette mthode est intensivement employe. Elle permet un
calcul de prcision acceptable mme lorsque le pas
est assez grand, ce qui permet davoir un moins grand nombre
dtapes globales (inversion de la matrice de rigidit du
systme).
Remarque 3 :
Il existe beaucoup dvolutions de cette mthode. Par exemple,
la rigidit tangente cohrente, etc Se rfrer la bibliographie
si besoin.
88
8.3
La -mthode
Algorithme propos par Benallal et coll. dans :
An integration algorithm and the corresponding consistent
tangent operator for fully coupled elastoplastic and damage
equations, Comm. in Appl. and Num. Meth. 4, pp. 731-740.
Particulirement adapt au formalisme des matriaux standart
gnraliss.
Existe aussi pour les comportements visqueux (dpendance au
temps physique).
Pas des plus rapides mais trs prcis.
89
Formalisme
Le pas de temps prcdent connu est not n, le suivant, que lon
calcule est n+1 (tapes globales du code).
Le vecteur regroupe des V.I. ou des V.E. (selon la
commodit, puis quelles sont relies par le potentiel dtat). On
peut avoir des scalaires et/ou des tenseurs.
ECL :
EI :
ECNL:rem. : on aurait pu rendre
En utilisant une base de tenseurs utilise 6 composantes, et
une seule donc q comte 7 composantes pour lEI.
Les quations de plasticit se mettent sous la forme de 3
quations :
Surface seuil :
coulement :
lasticit :
ECL :
SDC, ETAT
EVO1, 2, ETAT
ELASK
91
Le code E.F. demande :
la contrainte
les variables internes
en fonction des donnes
la dformation teste
et des variables internes
92
Prdicteur lastique
On commence par tester une volution en lasticit car cest un
calcul direct donc trs rapide:
Esp. dv. (5D)
si
alors
93
Cas de la transformation plastique
Si le prdicteur lastique donne une contrainte en dehors du
domaine dlasticit:
Si
Alors on a faire avec les quations de plasticit. Dun point de
vue analytique, pour une transformation infinitsimale, on a:
Seuil
coulement
lasticit
On doit discrtiser ces quations. Le gradient g (=dq/d) sera
pris la valeur intermdiaire entre les points et
.
Cest la -mthode.
94
Calcul explicite:
Calcul implicite :
La -mthode:
Le calcul de la pente
pour chaque mthode
Le thorme de la valeur moyenne montre quil existe un
pour lequel la pente entre les deux points est
trouve exactement.
On montre que pourun schma dintgration est
inconditionnellement stable.
95
Donc la discrtisation du systme dquations devient :
que lon rsume sous la forme
dans laquelle
On remarque que G comme Q possdent 7 dimensions de plus
que q (donc 13 pour lECL). Le 0 de G, au point n+1, est trouv
par un algorithme de Newton. On sous-discrtise donc le
problme en des points s, compris entre n et n+1 :
96
6,7 ou 12 Dim.
6 Dim.
im
Tous calculs faits (), ce systme revient la rsolution de ce
problme matriciel:
et de le rsoudre jusqu convergence (dtecte par la norme du
correcteur complet ou simplement de C).
97
Les calculs ncessaires sont relativement simples.
Pour lECL, depuis T16 de ce cours :
Rem. : cest aussi
Remarque : est le projecteur sur
98
8.4
Des tests
Prsentation de MatSGen, logiciel de pilotage en -mthode,
pour un lment de volume unique (ou pour une structure
champ homogne). Cest une aide pour dvelopper des modles
ou des UMAT.
Problme de lessai de traction qui prsente un chargement la
fois en contraintes et en dformations imposes.
Benchmark 1
Trajet 1D TCS
Chargement en traction +/-1500 MPa
modle ECNL :
99
Choix de
1 0 0 0 0 0
1500
S
S
Contrainte
1000
500
0
500
1000
1500
0.01 0.005
=0.0 : Efinal = 0.00103
=0.5 : Efinal = 0.00104
=1.0 : Efinal = 0.00105
0
0.005
Deformation
0.01
0.015
100
Nombre de pas (=0.5)
1 0 0 0 0 0
1500
Contrainte
1000
500
0
500
1000
1500
0.01 0.005
3 pas : Efinal = 0.00104
10 pas : Efinal = 0.00104
100 pas : Efinal = 0.00104
0
0.005
Deformation
0.01
0.015
101
0.7
Benchmark 2
Rochet bi-dimensionnel
Chargement en traction 1000 MPa puis torsion +/-1%
modle ECNL, mme3constantes dlasticit
0
x 10
5
0
5
10
1 0 0 0 0 0
15
20
3
x 10
=0.5 : E11 final = 0.0163 ; 200 pas : 0.0163 ; 3 pas : 0.0167
=1.0 : E11 final = 0.0164
102
0.7
Dpendance importante au pas de temps :
3
x 10
4
2
0
2
4
6
EE
10
1 0 0 0 0 0
15
20
3
x 10
Mme calcul de rochet bi-dimensionnel pour 3, 20 et 200 pas
de temps par cycle.
103
9
Bibliographie
104
Bibliographie
D. Dureisseix, Mthodes numriques appliques la conception par
lments finis. [Link]
PDF/[Link].
M. Ortiz and J. C. Simo. An analysis of a new class of
integration algorithms for elasto-plastic constitutive relations.
International Journal for Numerical Methods in Engineering,
23(3) :353366, 1986.
M. Ortiz and E. P. Popov. Accuracy and stability of integration
algorithms for elastoplastic constitutive relations. International
Journal for Numerical Methods in Engineering, 21(9) :1561
1576, 1985.
K.J. Bathe. Finite Element Procedures.
105
Abaqus eory Manual
J.L. Batoz
106
