Code_Aster ®
Version
5.0
Titrate:
Integration of the viscoelastic relations in STAT_NON_LINE
Date:
02/02/01
Author (S):
P. of BONNIERES
Key:
R5.03.08-A
Page:
1/14
Organization (S): EDF/MTI/MN
Handbook of Référence
R5.03 booklet: Nonlinear mechanics
Document: R5.03.08
Integration of the relations of behavior
viscoelastic in operator STAT_NON_LINE
Summary
This document describes the viscoelastic behaviors in the case of the ingredients necessary to the setting in
work of non-linear algorithm STAT_NON_LINE described in [R5.03.01]. Data input of all them
viscoelastic relations of behavior integrated in Aster have in a general way the same form. Only
the way of introducing the principal data (the function speed of viscous deformation) varies: it is presented
according to the various key words which make it possible the user to choose the relation of behavior wished.
These quantities are calculated by a method of integration semi-implicit. From the initial state, or from
the moment of preceding calculation, one calculates the stress field resulting from an increment of deformation.
Handbook of Référence
R5.03 booklet: Nonlinear mechanics
HI-75/01/001/A
Code_Aster ®
Version
5.0
Titrate:
Integration of the viscoelastic relations in STAT_NON_LINE
Date:
02/02/01
Author (S):
P. of BONNIERES
Key:
R5.03.08-A
Page:
2/14
Contents
1 Introduction ............................................................................................................................................ 3
2 Relation continues ................................................................................................................................... 3
3 Nature of the function G for each relation of behavior ................................................. 4
3.1 Relation LEMAITRE ......................................................................................................................... 4
3.2 Relations ZIRC_CYRA2 and ZIRC_EPRI .......................................................................................... 4
4 Integration of the relation of behavior ........................................................................................... 7
4.1 Establishment of the scalar equation for the implicit scheme and with elastic coefficients
constant ......................................................................................................................................... 7
4.2 Resolution of the scalar equation: principle of routine ZEROF2 .................................................. 8
4.3 Calculation of the constraint at the end of the step of current time ............................................................... 10
4.4 Semi-implicit diagram .................................................................................................................. 11
4.5 Taking into account of the variation of the elastic coefficients according to the temperature .......... 12
5 Calculation of the tangent operator ............................................................................................................... 13
6 Bibliography ........................................................................................................................................ 14
Handbook of Référence
R5.03 booklet: Nonlinear mechanics
HI-75/01/001/A
Code_Aster ®
Version
5.0
Titrate:
Integration of the viscoelastic relations in STAT_NON_LINE
Date:
02/02/01
Author (S):
P. of BONNIERES
Key:
R5.03.08-A
Page:
3/14
1 Introduction
The tubes of sheath in Zircaloy of the fuel pin of power stations REP present a behavior
strongly viscous mechanics.
Within the framework of the chaining enters Code_Aster and the code of the fuel pin CYRANO3, two
non-linear viscoelastic models specific to Zircaloy were introduced into Code_Aster.
One of them is the model used in Cyrano2 and Cyrano3. Other was developed by the EPRI.
In addition, a model much more general and correspondent with other materials that Zircaloy has
also introduced summer. It is about the non-linear viscoelasticity of Lemaître, which can be brought back for
certain particular values of the parameters to a relation of viscoelastic behavior of
Norton.
For these three models, one supposes that the material is isotropic. They can be used in 3D, in
plane deformations (D_PLAN) and into axisymmetric (AXIS).
One presents in this note the equations constitutive of the models and their establishment in
Code_Aster.
2 Relation
continuous
One places oneself on the assumption of the small disturbances and one divides the tensor of the deformations into
an elastic part, a thermal part, a anelastic part (known) and a viscous part.
equations are then:
early = E + HT + has + v
= A (T) E
3 ~
!v = (
G eq, T) 2 eq
with:
2
: cumulated viscous deformation! =
! :
v!
3
v
~
1
: diverter of the constraints ~ = - Tr () I
3
3
~ ~
eq: equivalent constraint eq =
:
2
WITH (T): tensor of elasticity
Handbook of Référence
R5.03 booklet: Nonlinear mechanics
HI-75/01/001/A
Code_Aster ®
Version
5.0
Titrate:
Integration of the viscoelastic relations in STAT_NON_LINE
Date:
02/02/01
Author (S):
P. of BONNIERES
Key:
R5.03.08-A
Page:
4/14
3
Nature of the function G for each relation of
behavior
3.1 Relation
LEMAITRE
In this case, G is expressed explicitly (is a scalar here):
N
(
1
1
1
G, T) =
with
,
0
,
0 N 0
/
K 1 m
>
K
m
The data of the material characteristics are those provided under the key words factors
LEMAITRE or LEMAITRE_FO of operator DEFI_MATERIAU.
1
1
/LEMAITRE:
NR:
N UN_ SUR_ K: UN_ SUR_ M:
K
m
The Young modulus E and the Poisson's ratio are those provided under the key words factors
ELAS or ELAS_FO.
3.2 Relations
ZIRC_CYRA2 and ZIRC_EPRI
For these relations, G is not expressed explicitly. The behavior is represented by one
unidimensional creep test, with constant constraint, which utilizes time passed since
the moment when the constraint is applied. The relation of behavior is defined here by the data of
four functions F, G, F, G
1
1 2
2 describing the evolution of the viscous deformation in the course of time:
v = = f1 (T) g1 (, T) + f2 (T) g2 (, T)
éq 3.2-1
The function G is calculated then numerically by eliminating time T in the following way:
1) for a given triplet (, T), one solves in T the equation [éq 3.2-1] by the method of
Newton (see [bib2]). An approximation of the solution T is found (, T),
2) one obtains the value of the function G in (, T) by deriving compared to time the equation
[éq 3.2-1] (see [bib1]):
! =!
v
= (
G, T) = F “1 (T) g1 (, T) + F”
2 (T) g2 (, T)
and in substituent in this new equation the value of T (, T) found previously.
One finds the formulation uniaxial following:
! =!
v
= (
G, T) = F “1 (T (, T) g1 (, T) + F”
2 (T (, T) g2 (, T)
For each of two relations ZIRC_CYRA2 and ZIRC_EPRI, the form of the four functions
F, G, F, G
1
1 2
2 is preset and the user introduces only some parameters into the file of
order.
Handbook of Référence
R5.03 booklet: Nonlinear mechanics
HI-75/01/001/A
Code_Aster ®
Version
5.0
Titrate:
Integration of the viscoelastic relations in STAT_NON_LINE
Date:
02/02/01
Author (S):
P. of BONNIERES
Key:
R5.03.08-A
Page:
5/14
Thus, for ZIRC_CYRA2, one a:
F
- T
1 (T) = [Cth fab (1 - E
) +t] Frec
(0,00266 - 0,413)
2
2nd
g1 (, T) =
With ex
HT
(p K/(T +,
273)
15
3
3
F
- K T
irr
2 (T) = [Cirr fab (1 - E
) +t] Frec
4
G2 (, T) = A
3 irr
with: C = 4450
HT
= 4,
- 3
1
5 10 H
To = 9,
17
529 10
HT
K = 39000° K
F
=,
- 4
1.816 10 ex
rec
(p6400/(T +,
273
rec
)
15
C
= 4000
irr
- 3
1
K
= 3 10 H
irr
Airr =
-
1 2 10 22
,
Positive parameters fab, R
T EC. and are those provided under the key word factor ZIRC_CYRA2 of
operator DEFI_MATERIAU:
/ZIRC_ CYRA2: (EPSI_ FAB: TEMP_ REHEATED:T
FLUX_ PHI:
fab
rec
)
Handbook of Référence
R5.03 booklet: Nonlinear mechanics
HI-75/01/001/A
Code_Aster ®
Version
5.0
Titrate:
Integration of the viscoelastic relations in STAT_NON_LINE
Date:
02/02/01
Author (S):
P. of BONNIERES
Key:
R5.03.08-A
Page:
6/14
In the same way, for ZIRC_EPRI, one a:
F
5
With
1 (T) = T
With
G
A3
4
A6 - A7 T+273 15
1 (, T) = A (
1 S (
H2
With)
/(
,)
R E
p
F
B2
2 (T) = T
B
G
B3 B4 - B5 T+273 15
B6
7
2 (, T)
/(
,)
= 1
B E
RP (cos max
)
with:
1
With =,
8
1.603 10
2
To = 4,
5
567 10
3
To = 2,28
4
With =,
0 997
5
With =,
0 77
6
With =,
0 956
7
To = 2,
4
3 10
1
B =,
- 21
3.296 10
B2 =,
0 811
B3 =,
0 595
B4 = 1 352
,
5
B = 22 9
, 1
B6 = 1 5
, 8
B7 = 2 228
,
Positive parameters, RP and max 0 max
2 are those provided under the key word factor
ZIRC_EPRI of operator DEFI_MATERIAU:
/ZIRC_ EPRI: (FLUX_ PHI:
R_ P:R THETA_ MAX:
p
max
)
It will be noted that, for the two relations of behavior and all the functions:
T is expressed in hours
T is expressed in °C
express yourself in MPa
The document of Young E and the Poisson's ratio are those provided under the key words factors
ELAS or ELAS_FO.
Handbook of Référence
R5.03 booklet: Nonlinear mechanics
HI-75/01/001/A
Code_Aster ®
Version
5.0
Titrate:
Integration of the viscoelastic relations in STAT_NON_LINE
Date:
02/02/01
Author (S):
P. of BONNIERES
Key:
R5.03.08-A
Page:
7/14
4
Integration of the relation of behavior
4.1 Establishment of the scalar equation for the implicit scheme and
with constant elastic coefficients
One indicates by early the total deflection at the moment T + T
and by early the variation of deformation
total during the step of current time. One calls O the deformation imposed at the moment T + T
and
O variation of deformation imposed during the step of current time.
This imposed deformation results from thermal dilation and the anelastic deformations. One has
thus:
O = (T + T)
[
((Tt+t) - rTef) - (T) ((Tt) - rTef)]I3+a (t+t) - has (T)
where I3 is the tensor identity of command 2 in dimension 3.
One poses = early - O
As it is supposed here that µ is constant, one with the following relation between the diverters of
and:
~ = 2µ (~
- v)
éq 4.1-1
However, the law of flow is written, an implicit way:
~
v
3
-
= G eq, +
,
éq 4.1-2
2
(v) T
eq
T
eq
One thus has, by eliminating v between [éq 4.1-1] and [éq 4.1-2]:
~
-
2µ~ = ~ + 3µt G,
eq +
(v), T
eq
eq
-
G
,
éq 4.1-3
eq +
(v)
(
, T
~
eq
- + 2µ~) = 1+ 3µ
~
T
eq
By posing ~e
~-
~
= + 2µ, one thus has:
E
-
eq = 3µ T
G eq, +
(v), T
éq 4.1-4
eq
eq
+
Handbook of Référence
R5.03 booklet: Nonlinear mechanics
HI-75/01/001/A
Code_Aster ®
Version
5.0
Titrate:
Integration of the viscoelastic relations in STAT_NON_LINE
Date:
02/02/01
Author (S):
P. of BONNIERES
Key:
R5.03.08-A
Page:
8/14
However, one has according to [éq 4.1-2]:
(
-
v)
= T G, +
, T
eq
eq
(v) eq
From where:
eeq = 3µ (v) +
eq
eq
(
1
E
v)
=
-
eq
3 (eq eq)
µ
In substituent this last expression in [éq 4.1-4], one a:
1
E
-
E
eq = 3µ T
G eq, +
,
3 (eq - eq) T + eq
µ
If one poses, eeq, -, T and T
being known:
1
F (X)
µ T
G X, -
3
3 (E - X), T + X
E
=
+
eq
- eq
µ
one can then calculate the quantity
-
eq = (+) as being the solution of the scalar equation:
eq
F (X) = 0 where X = eq, convention adopted for the following paragraphs.
4.2
Resolution of the scalar equation: principle of routine ZEROF2
One easily shows that, if the requirements in the paragraph [§3] on the characteristics of
materials are checked, the function F is strictly increasing and the equation F (X) = 0 admits one
single solution.
If E
E
eq = 0, then the solution is X = 0. If not, one a: F ()
0 = - eq < 0
The problem thus consists in finding for a function F unspecified the solution of the equation
F (X) = 0 knowing that this solution exists, that F ()
0 < 0 and that F is strictly increasing.
The algorithm adopted in ZEROF2 is as follows:
· one leaves A = 0 and B = X
where X
0
0
ap
ap is an approximation of the solution. If it is
necessary (i.e. if F (b0) < 0), one brings back oneself by the method of the secants
year F (bn) - bn F (year)
(Zn =
then has 1 = B and B 1 = Z) in one or more iterations with the case
F (B
n+
N
n+
N
N) - F (year)
where F (A) < 0 and F (b) > 0:
Handbook of Référence
R5.03 booklet: Nonlinear mechanics
HI-75/01/001/A
Code_Aster ®
Version
5.0
Titrate:
Integration of the viscoelastic relations in STAT_NON_LINE
Date:
02/02/01
Author (S):
P. of BONNIERES
Key:
R5.03.08-A
Page:
9/14
F (b1)
= B has
1
0
0
b1
F (b0)
F ()
0
(In the case of the figure above, this first sentence was done in an iteration:
has = B and F (B
1
0
1) > 0).
· one
calculate
Nd = left whole (Nmax) or Nmax is the maximum number of iterations that
one was given. One then solves the equation by the method of the secants while using however
method of dichotomy with each time N is multiple of Nd:
1)
If Nd divides N
+ B has
Z
N
N
N =
2
if not
year F (bn) - bn F (year)
Zn =
F (bn) - F (year)
finsi
N = N + 1
if F (Z) >
if F (Z) < 0
1 = Z has
B 1 = B
n+
N
n+
N
if not
has 1 = has
B 1 = Z
n+
N
n+
N
finsi
to go into 1)
if not
The solution is: X = Z FIN
N
finsi
This second part of the algorithm makes it possible to treat in a reasonable iteration count the cases
where F is very strongly non-linear, whereas the method of the secants would have converged too much
slowly. These cases of strong non-linearity meet in particular with the law of LEMAITRE, for
N
values of
large.
m
Handbook of Référence
R5.03 booklet: Nonlinear mechanics
HI-75/01/001/A
Code_Aster ®
Version
5.0
Titrate:
Integration of the viscoelastic relations in STAT_NON_LINE
Date:
02/02/01
Author (S):
P. of BONNIERES
Key:
R5.03.08-A
Page:
10/14
4.3
Calculation of the constraint at the end of the step of current time
According to [éq 4.1-3], if X is the solution of the scalar equation, while posing:
(
1
X
B X, E
eq) =
= E
1
-
G X, + 3 (E
- X
eq
eq
), T
µ
1 + 3µ T
X
one a:
~ = B (X, E E
eq) ~
éq 4.3-1
If eeq = 0, which is equivalent according to the scalar equation to X = 0, one prolongs B by
1
continuity. For that, one poses (
y X)
-
(E
=
+
- X
G X = G (X, there X, T)
eq
)
µ
and ()
(). The derivative of G
3
express yourself according to the derivative partial of G at the point (X (
, there X), T):
G
1 G
G' (X) =
(X (
, there X), T) -
(X (
, there X), T)
X
µ
3 y
The prolongation of B by continuity gives then:
(
1
B,
0)
0 = 1+ 3µtG' () 0
and one has, always if E
~
eq = 0, = 0.
Once one calculated ~
, one obtains by the relation (K is supposed to be constant here):
= - + = ~
1
+ Tr (-
) + KTr ()
I
éq 4.3-2
3
3
Handbook of Référence
R5.03 booklet: Nonlinear mechanics
HI-75/01/001/A
Code_Aster ®
Version
5.0
Titrate:
Integration of the viscoelastic relations in STAT_NON_LINE
Date:
02/02/01
Author (S):
P. of BONNIERES
Key:
R5.03.08-A
Page:
11/14
4.4 Diagram
semi-implicit
With an implicit numerical diagram [éq 4.1-2], in the case, for example, where G does not depend on,
only intervenes by the calculation of the v value of the constraint at the end of the step of time. It can in
to result from the important numerical errors if the constraint strongly varies in the course of time (see
[bib2]).
To cure that and to improve the resolution, one discretizes the law of flow in way
semi-implicit:
~
~-
+
v
3
-
(v) eq
T
2
= G +
-
,
, T -
éq 4.4-1
T
2
2
+
+
2
2 -
eq
+
2 eq
To transform in the most economic way what was programmed previously (while following
implicit formulation [éq 4.1-2]), it is enough to divide each member of the equation [éq 4.4-1] by 2:
-
(v) eq
-
-
T
G +
,
, T
- ~
~
(
2
+
+
2
2
+
v/2)
eq
3
2
=
T
2
2
-
+
2 eq
and to make the same thing with the relation [éq 4.1-1]:
~
~
= 2µ
-
v
2
2
2
It is noted that this system is same form as that consisted the equations [éq 4.1-1] and
[éq 4.1-2], the data being
instead of, unknown factors being respectively
and
v with
2
2
2
G
place of and v and the function replacing the function G.
2
One can thus use the resolution of the paragraphs [§4.1] with [§4.3] as well as the corresponding algorithm
while introducing
and by dividing the function G by 2. It then remains to multiply the results
2
and
v
by 2 to obtain the increments of calculated constraint and viscous deformation
2
2
by the semi-implicit diagram (it and the v of the equation [éq 4.4-1]).
Handbook of Référence
R5.03 booklet: Nonlinear mechanics
HI-75/01/001/A
Code_Aster ®
Version
5.0
Titrate:
Integration of the viscoelastic relations in STAT_NON_LINE
Date:
02/02/01
Author (S):
P. of BONNIERES
Key:
R5.03.08-A
Page:
12/14
It will be noticed that the calculation of the tangent operator is not affected by this modification of the diagram
numerical. Indeed, one has obviously:
()
2
=
2
4.5 Taking into account of the variation of the elastic coefficients in
function of the temperature
One has, if A is the tensor of elasticity:
= + (- 1
v
With)
with:
(A-1) A-1 (-
T
T) (-) A-1
(- T) -
=
+
+
-
This is translated in the equations of [§4.4] by:
-
(v)
-
eq
-
T
G +
,
,
T
-
~
~
~
2
+
+
2
2
+
eq
-
-
2
2µ + 2µ
2µ
~
3µ
~-
T
-
2 -
+
2 =
2
-
-
4µ
+
2 eq
While posing:
-
~
~
2µ
2µ
E
+
~
-
=
+
2µ
-
4µ
2
and
3K- + 3K
-
Tr (E
) =
Tr
-
() +
3K Tr
6K
2
one is reduced exactly to the preceding case [§4.4].
Handbook of Référence
R5.03 booklet: Nonlinear mechanics
HI-75/01/001/A
Code_Aster ®
Version
5.0
Titrate:
Integration of the viscoelastic relations in STAT_NON_LINE
Date:
02/02/01
Author (S):
P. of BONNIERES
Key:
R5.03.08-A
Page:
13/14
5
Calculation of the tangent operator
If eeq = 0 and X = 0, one take the tensor of elasticity as tangent operator.
If not, one obtains this operator by deriving the equation [éq 4.3-1] compared to:
~
~
(
B X, E
E
eq)
~
~e
E
=
=
+ (
B X, eq)
then while also deriving [éq 4.3-2] compared to:
~
Tr () ~
=
+ K 3
=
+ K T
I
I I
3 3
It will be noted that, in these equations, the tensors of command 2 and command 4 are respectively compared to
vectors and with matrices. I3 is here a tensor of a nature 2, compared to a vector:
T I3 = (11,1,0, 0,)
0
One has moreover:
(bx, E
E
eq) B
X
=
X E
+
X E
eq
X (
B
, eq)
E (, eq)
eq
X
It is thus necessary to calculate. For that, one derives the scalar equation implicitly compared to.
To simplify, one will omit thereafter in the writing of G and his derivative the parameter T.
One has then:
E
E
[
X
G
eq
3µt G (X) +]
1
+ T
(X, y)
eq
=
y
From where:
G
1 - T
(X, y)
X
y
E
eq
=
1+ µ
3
T
G (X)
G
1 - T
(X, y)
X
y
µ
3 T ~e
=
1+ µ
3
T
G (X) E
eq
with the expression of G' (X) obtained with [§4.3].
Handbook of Référence
R5.03 booklet: Nonlinear mechanics
HI-75/01/001/A
Code_Aster ®
Version
5.0
Titrate:
Integration of the viscoelastic relations in STAT_NON_LINE
Date:
02/02/01
Author (S):
P. of BONNIERES
Key:
R5.03.08-A
Page:
14/14
One obtains finally the following expression of the tangent operator:
= K T
I I
~
~
3
3 + 2µ E T E + B X, E
With
[
(eq)]
with:
~
1
With =
= J6 -
T
I I where J is the matrix identity of row 6.
3
3
6
3
G
1 -
,
3
T
(X y)
E
=
y
(
3
eq
- X
2nd
1 3
G
eq)
+ µt (X)
Note:
In the case of laws ZIRC_CYRA2 and ZIRC_EPRI, it is checked easily that:
1
2
2
G (X) =
'
G G F '
1 1 1 -
'
''
'
'
''
'
'
'
F F
1 1 + G G F
2 2
2 -
F F
2 2 + G G
1 2 (F F
1 2 - F F
1 2)
F G
1 1+ F G
2 2
'
''
'
1
'
'
+ G G
2 1 (F F
1 2 - F F
1 2) -
(F g11+f g22)
3µ
G (
F '' G
F 'G
X, y, T)
1 1+
=
2 2
F 'G
1 1+ F '2 g2
where F, F “F '' F, F” F '
1
,
1
,
1 2 2, 2 indicate the values of F and F
1
2 and of their derivative at the point T (X, y, T) and
where G, g', G, g'
1
1
2
2 indicate the values of G and G
1
2 and of their derivative compared to at the point
(X, T) (see [bib1]).
6 Bibliography
[1]
BONNIERES P.: Writing in generalized standard form of the laws of behavior
viscoplastic of Zircaloy, notes EDF-DER HI-71/7940-Index A, 1992
[2]
BONNIERES P., ZIDI Mr.: Introduction of viscoplasticity into the module of
thermomechanics of Cyrano3: principle, description and validation, note EDF-DER HI-71/8334,
1993
Handbook of Référence
R5.03 booklet: Nonlinear mechanics
HI-75/01/001/A
Outline document