# A Study on Compressive Anisotropy and Nonassociated Flow Plasticity of the AZ31 Magnesium Alloy in Hot Rolling

Authors: Guoqiang Wang; Xiaolei Qian; Xuefei Li; Haijing Hou; Yi Liu; Yanshan Lou
Journal: Mathematical Problems in Engineering
 (2014)

---

## Abstract

Effect of anisotropy in compression is studied on hot rolling of AZ31 magnesium alloy with a three-dimensional constitutive model based on the quadratic Hill48 yield criterion and nonassociated flow rule (non-AFR). The constitutive model is characterized by compressive tests of AZ31 billets since plastic deformations of materials are mostly caused by compression during rolling processes. The characterized plasticity model is implemented into ABAQUS/Explicit as a user-defined material subroutine (VUMAT) based on semi-implicit backward Euler's method. The subroutine is employed to simulate square-bar rolling processes. The simulation results are compared with rolled specimens and those predicted by the von Mises and the Hill48 yield function under AFR. Moreover, strip rolling is also simulated for AZ31 with the Hill48 yield function under non-AFR. The strip rolling simulation demonstrates that the lateral spread generated by the non-AFR model is in good agreement with experimental data. These comparisons between simulation and experiments validate that the proposed Hill48 yield function under non-AFR provides satisfactory description of plastic deformation behavior in hot rolling for AZ31 alloys in case that the anisotropic parameters in the Hill48 yield function and the non-associated flow rule are calibrated by the compressive experimental results.


---

## Body

## 1. Introduction

Anisotropic mechanical characteristics of the AZ31 magnesium alloy have a significant influence on plastic deformation and should be considered carefully for accurate simulation of forming processes. For decades many researchers have focused on the anisotropic constitutive models of metallic sheets which are all rolled products, whereas rare anisotropic models are adopted in hot rolling simulations. Montmitonnet et al. [1] simulated the hot rolling processes of steels but did not present experimental validation. Massé et al. [2] studied the lateral spread of pearlitic steel flat wire in cold rolling and concluded that taking into account plastic anisotropy significantly improves the estimation of the final width. Except for reported results above, effect of anisotropy has not been comprehensively investigated during hot rolling processes especially for hot rolling of AZ31 magnesium alloys with strong effect of anisotropy.The foundation of most anisotropic constitutive models has been based on the associated flow rule (AFR) hypothesis which states that the flow rule is associated with the yield criterion. In other words, the yield function is also the potential for plastic strain rate in AFR based models. Accordingly, starting from Hill’s quadratic anisotropy model [3], various yield functions have been proposed to describe the anisotropy of metallic sheets. A detailed review can be found in Mohsen Safaei’s doctoral dissertation [4]. More recently, strength-differential effect was extensively studied and modeled by several yield functions proposed, such as Cazacu and Barlat [5], Cazacu et al. [6], Lou et al. [7], and Yoon et al. [8].Most yield functions proposed are coupled with AFR. However, various studies described the difficulty of the AFR concept in dealing with highly anisotropic materials. For instance, Cvitanić et al. [9] showed that Hill 1948 and Karafillis-Boyce 1993 anisotropic yield functions cannot predict both directionalR-values and yield stresses simultaneously for AA2008-T4 and AA2090-T3. Yoon et al. [10] showed that Hill 1948 is unable to predict the exact numbers of ears in deep drawn cups made of AA2090-T3. Park and Chung [11] reported that Hill 1948 and Yld2000-2d generate poor accuracy for directionalR-values and yield stresses for AA2090-T3 and AA5042. Therefore it can be concluded that it is difficult to describe a highly anisotropic material in terms of both plastic strain rate and yielding behavior with an identical function for yielding criterion and plastic potential.The nonassociated flow rule approach neglects the constraint of equality of plastic potential and yield function enforced by the AFR assumption and can be an alternative but efficient approach to improve the accuracy of anisotropic yield functions even though there are some disadvantages related to non-AFR, such as issues about stability and thermodynamic consistency. Accordingly, the yield function and plastic potential function are assumed to be independent and used to describe the elastic limit and plastic strain rate direction, respectively. This assumption leads to various advantages and flexibility such that a larger number of experimental data are used for calibration of the parameters of yield and plastic potential functions, which results in a better agreement between simulation and experimental data, for example, better prediction of yield stresses and Lankford coefficients (width to thickness plastic strain ratios) at multiple in-plane orientations. WhenR-values and yield stresses from a highly anisotropic material are fit to one function, large gradients in the curvature of that function may occur and result in numerical convergence problems. These numerical convergence problems are solved by the non-AFR approach because these nonphysical curvature gradients are reduced by non-AFR [12].Stoughton [13] first accurately predicted theR-value distribution and uniaxial and biaxial yield behaviors based on a non-AFR Hill48 model. After that, the non-AFR based approach has attracted increasing attentions for metal plasticity. Cvitanić et al. [9] developed a non-AFR model based on Hill 1948 quadratic and Karafillis and Boyce nonquadratic yield functions which showed an improved prediction of the heights of deep drawn cups made of an aluminum alloy. Stoughton and Yoon [14] proposed a material model based on non-AFR to describe anisotropic hardening of materials under proportional loading. The anisotropic hardening model under non-AFR was applied to five experimental results of aluminum alloys and stainless steels for the verification purpose. A nonassociated, mixed hardening model proposed by Taherizadeh et al. [15] significantly improves the prediction of earing in the cup drawing process and the prediction of springback in the sidewall of drawn channel sections. Recently, Park and Chung [11] proposed a new formulation leading to the symmetric stiffness modulus for the nonassociated flow rule, which was validated for earing in circular cup drawing of AA2090-T3 and AA5042 sheets. Safaei et al. [16] proposed an evolutionary anisotropic model based on non-AFR that excellently predicted distortional hardening and evolution of instantaneousR-values in seven uniaxial directions as well as balanced biaxial loading condition.In this paper, effect of anisotropy in compression on hot rolling of AZ31 magnesium alloy is studied using a non-AFR constitutive model based on the Hill 1948 yield function. The Hill48 yield function under non-AFR is calibrated byR-values and yield stresses from compression tests since materials mainly undergo compressive loading during hot rolling. The calibrated constitutive model is implemented into the commercial FE code ABAQUS/Explicit [17] using semi-implicit backward Euler's algorithm since the semi-implicit backward Euler method proposed by Moran et al. [18] reduces the equation solving effort during the state update to the solution of a single scalar nonlinear equation for the incremental effective plastic strain, regardless of the material model. The Hill48 yield function under non-AFR is applied to simulate square-bar rolling process. The hot rolling processes are also simulated by the von Mises yield function and the Hill48 yield function under AFR for the comparison purpose. Moreover, the proposed Hill48 yield function under non-AFR is applied to strip rolling for the further evaluation of the proposed yield function on hot rolling processes.


## 2. Constitutive Model

### 2.1. Anisotropic Yield Function

The three-dimensional yield function proposed by Hill [3] is defined as
(1)
f
y
=
F
y
(
σ
22
-
σ
33
)
2
+
G
y
(
σ
33
-
σ
11
)
2
+
H
y
(
σ
11
-
σ
22
)
2
+
2
L
y
σ
23
2
+
2
M
y
σ
31
2
+
2
N
y
σ
12
2
-
σ
-
y
2
,
where F
y
,
G
y
,
H
y
,
L
y
,
M
y
,
N
y are used to describe the direction-dependent yield stresses. To calibrate these coefficients, consider a uniaxial compression at an angle θ from the rolling direction and denote the uniaxial compressive yield stress σ
θ. The stress components in the Cartesian axis system (orthotropic axes 1, 2, and 3 with 1-axis in the rolling direction, 2-axis in the transverse direction, and 3-axis in the thickness direction) are
(2)
σ
11
=
σ
θ
cos
⁡
2
θ
,
σ
22
=
σ
θ
si
n
2
θ
,
σ
12
=
σ
θ
sin
θ
cos
⁡
θ
.The substitution of (2) into (1) provides
(3)
σ
θ
=
(
-
2
H
y
si
n
2
θ
co
⁡
s
2
θ
+
2
N
y
si
n
2
θ
co
⁡
s
2
θ
)
-
1

(
1
)
(
(
F
y
+
H
y
)
si
n
4
θ
+
(
G
y
+
H
y
)
cos
⁡
4
θ
-
2
H
y
si
n
2
θ
cos
⁡
2
θ
+
2
N
y
si
n
2
θ
cos
⁡
2
θ
)
-
1
)
1
/
2
σ
-
y
.
Specifically, the uniaxial compressive yield stresses are formulated for the rolling direction (0°), diagonal direction (45°), and transverse direction (90°) as follows:
(4)
σ
0
=
(
1
G
y
+
H
y
)
1
/
2
σ
-
y
,
σ
45
=
(
4
F
y
+
G
y
+
2
N
y
)
1
/
2
σ
-
y
,
σ
90
=
(
1
F
y
+
H
y
)
1
/
2
σ
-
y
.The uniaxial compression yield stress in the vertical direction denoted asσ
z can be easily expressed as below:
(5)
σ
z
=
(
1
F
y
+
G
y
)
1
/
2
σ
-
y
.Four anisotropic parameters are formulated by solving equations of (4) and (5) as follows:
(6)
F
y
=
1
2
(
1
σ
90
2
-
1
σ
0
2
+
1
σ
z
2
)
σ
-
y
2
,
(7)
G
y
=
1
2
(
1
σ
0
2
-
1
σ
90
2
+
1
σ
z
2
)
σ
-
y
2
,
(8)
H
y
=
1
2
(
1
σ
0
2
+
1
σ
90
2
-
1
σ
z
2
)
σ
-
y
2
,
(9)
N
y
=
(
2
σ
45
2
-
1
2
σ
z
2
)
σ
-
y
2
.The through-thickness shear-related anisotropic parametersL
y and M
y are assumed to be identical with the shear-related anisotropic parameter of N
y and calculated as below:
(10)
L
y
=
M
y
=
N
y
.


### 2.2. Nonassociated Flow Rule

Assuming isotropic linear elasticity and additive decomposition of the strain increment tensord
ε into elastic part d
ε
e and plastic part d
ε
p, the stress increment tensor reads
(11)
d
σ
=
C
e
:
d
ε
e
=
C
e
:
(
d
ε
-
d
ε
p
)
,
where d
σ is the stress increment tensor, d
ε is the total strain increment tensor, d
ε
e is the elastic strain increment tensor, d
ε
p is the plastic strain increment tensor, and C
e is the tensor of elastic moduli. In non-AFR based models, the two independent formulations of yield stress and plastic flow direction determine the yielding and directional plastic strain rate, respectively. The plastic flow direction is determined by the plastic potential f
p as below:
(12)
d
ε
p
=
d
λ
∂
f
p
(
σ
)
∂
σ
.The plastic potential functionf
p is assumed to have a Hill48 form of
(13)
f
p
=
F
p
(
σ
22
-
σ
33
)
2
+
G
p
(
σ
33
-
σ
11
)
2
+
H
p
(
σ
11
-
σ
22
)
2
+
2
L
p
σ
23
2
+
2
M
p
σ
31
2
+
2
N
p
σ
12
2
-
σ
-
p
2
.
Here f
p is a continuously differentiable function and d
λ is a nonnegative scalar called plastic multiplier or consistency parameter. By partial derivatives, the specific plastic strain increments are derived:
(14)
d
ε
11
=
2
d
λ
[
H
p
(
σ
11
-
σ
22
)
-
G
p
(
σ
33
-
σ
11
)
]
,
d
ε
22
=
2
d
λ
[
F
p
(
σ
22
-
σ
33
)
-
H
p
(
σ
11
-
σ
22
)
]
,
d
ε
33
=
2
d
λ
[
G
p
(
σ
33
-
σ
11
)
-
F
p
(
σ
22
-
σ
33
)
]
,
d
ε
12
=
d
ε
21
=
4
d
λ
N
p
σ
12
,
d
ε
23
=
d
ε
32
=
4
d
λ
L
p
σ
23
,
d
ε
31
=
d
ε
13
=
4
d
λ
M
p
σ
31
.For the sample compressed at the direction ofθ
° from the rolling direction, the Lankford coefficient is
(15)
R
θ
=
d
ε
θ
+
9
0
°
d
ε
33
=
-
(
(
∂
f
p
∂
σ
11
+
∂
f
p
∂
σ
22
)
-
1
(
∂
f
p
∂
σ
11
si
n
2
θ
+
∂
f
p
∂
σ
22
cos
⁡
2
θ
-
∂
f
p
∂
σ
12
sin
θ
cos
⁡
θ
)
×
(
∂
f
p
∂
σ
11
+
∂
f
p
∂
σ
22
)
-
1
)
=
(
(
G
p
+
F
p
ta
n
2
θ
)
-
1
(
H
p
+
H
p
ta
n
2
θ
si
n
2
θ
+
2
N
p
si
n
2
θ
-
(
F
p
+
G
p
+
3
H
p
)
si
n
2
θ
)
×
(
G
p
+
F
p
ta
n
2
θ
)
-
1
)
.
Based on the experimental data of R
0
,
R
45 and R
90 as well as the uniaxial compressive yield stress along the through-thickness direction of σ
z, the anisotropic parameters of the Hill48 plastic potential function are obtained:
(16)
F
p
=
R
0
(
R
0
+
R
90
)
σ
z
2
,
G
p
=
R
90
(
R
0
+
R
90
)
σ
z
2
,
H
p
=
R
0
R
90
(
R
0
+
R
90
)
σ
z
2
,
N
p
=
2
R
45
+
1
2
σ
z
2
.
The anisotropic parameters of L
p and M
p are assumed to be identical with the value of N
p. It should be noticed that the difference between the stress-shaped yield function andR-value shaped plastic potential function of Hill 1948 is the type of experimental inputs for parameter identification.


## 2.1. Anisotropic Yield Function

The three-dimensional yield function proposed by Hill [3] is defined as
(1)
f
y
=
F
y
(
σ
22
-
σ
33
)
2
+
G
y
(
σ
33
-
σ
11
)
2
+
H
y
(
σ
11
-
σ
22
)
2
+
2
L
y
σ
23
2
+
2
M
y
σ
31
2
+
2
N
y
σ
12
2
-
σ
-
y
2
,
where F
y
,
G
y
,
H
y
,
L
y
,
M
y
,
N
y are used to describe the direction-dependent yield stresses. To calibrate these coefficients, consider a uniaxial compression at an angle θ from the rolling direction and denote the uniaxial compressive yield stress σ
θ. The stress components in the Cartesian axis system (orthotropic axes 1, 2, and 3 with 1-axis in the rolling direction, 2-axis in the transverse direction, and 3-axis in the thickness direction) are
(2)
σ
11
=
σ
θ
cos
⁡
2
θ
,
σ
22
=
σ
θ
si
n
2
θ
,
σ
12
=
σ
θ
sin
θ
cos
⁡
θ
.The substitution of (2) into (1) provides
(3)
σ
θ
=
(
-
2
H
y
si
n
2
θ
co
⁡
s
2
θ
+
2
N
y
si
n
2
θ
co
⁡
s
2
θ
)
-
1

(
1
)
(
(
F
y
+
H
y
)
si
n
4
θ
+
(
G
y
+
H
y
)
cos
⁡
4
θ
-
2
H
y
si
n
2
θ
cos
⁡
2
θ
+
2
N
y
si
n
2
θ
cos
⁡
2
θ
)
-
1
)
1
/
2
σ
-
y
.
Specifically, the uniaxial compressive yield stresses are formulated for the rolling direction (0°), diagonal direction (45°), and transverse direction (90°) as follows:
(4)
σ
0
=
(
1
G
y
+
H
y
)
1
/
2
σ
-
y
,
σ
45
=
(
4
F
y
+
G
y
+
2
N
y
)
1
/
2
σ
-
y
,
σ
90
=
(
1
F
y
+
H
y
)
1
/
2
σ
-
y
.The uniaxial compression yield stress in the vertical direction denoted asσ
z can be easily expressed as below:
(5)
σ
z
=
(
1
F
y
+
G
y
)
1
/
2
σ
-
y
.Four anisotropic parameters are formulated by solving equations of (4) and (5) as follows:
(6)
F
y
=
1
2
(
1
σ
90
2
-
1
σ
0
2
+
1
σ
z
2
)
σ
-
y
2
,
(7)
G
y
=
1
2
(
1
σ
0
2
-
1
σ
90
2
+
1
σ
z
2
)
σ
-
y
2
,
(8)
H
y
=
1
2
(
1
σ
0
2
+
1
σ
90
2
-
1
σ
z
2
)
σ
-
y
2
,
(9)
N
y
=
(
2
σ
45
2
-
1
2
σ
z
2
)
σ
-
y
2
.The through-thickness shear-related anisotropic parametersL
y and M
y are assumed to be identical with the shear-related anisotropic parameter of N
y and calculated as below:
(10)
L
y
=
M
y
=
N
y
.


## 2.2. Nonassociated Flow Rule

Assuming isotropic linear elasticity and additive decomposition of the strain increment tensord
ε into elastic part d
ε
e and plastic part d
ε
p, the stress increment tensor reads
(11)
d
σ
=
C
e
:
d
ε
e
=
C
e
:
(
d
ε
-
d
ε
p
)
,
where d
σ is the stress increment tensor, d
ε is the total strain increment tensor, d
ε
e is the elastic strain increment tensor, d
ε
p is the plastic strain increment tensor, and C
e is the tensor of elastic moduli. In non-AFR based models, the two independent formulations of yield stress and plastic flow direction determine the yielding and directional plastic strain rate, respectively. The plastic flow direction is determined by the plastic potential f
p as below:
(12)
d
ε
p
=
d
λ
∂
f
p
(
σ
)
∂
σ
.The plastic potential functionf
p is assumed to have a Hill48 form of
(13)
f
p
=
F
p
(
σ
22
-
σ
33
)
2
+
G
p
(
σ
33
-
σ
11
)
2
+
H
p
(
σ
11
-
σ
22
)
2
+
2
L
p
σ
23
2
+
2
M
p
σ
31
2
+
2
N
p
σ
12
2
-
σ
-
p
2
.
Here f
p is a continuously differentiable function and d
λ is a nonnegative scalar called plastic multiplier or consistency parameter. By partial derivatives, the specific plastic strain increments are derived:
(14)
d
ε
11
=
2
d
λ
[
H
p
(
σ
11
-
σ
22
)
-
G
p
(
σ
33
-
σ
11
)
]
,
d
ε
22
=
2
d
λ
[
F
p
(
σ
22
-
σ
33
)
-
H
p
(
σ
11
-
σ
22
)
]
,
d
ε
33
=
2
d
λ
[
G
p
(
σ
33
-
σ
11
)
-
F
p
(
σ
22
-
σ
33
)
]
,
d
ε
12
=
d
ε
21
=
4
d
λ
N
p
σ
12
,
d
ε
23
=
d
ε
32
=
4
d
λ
L
p
σ
23
,
d
ε
31
=
d
ε
13
=
4
d
λ
M
p
σ
31
.For the sample compressed at the direction ofθ
° from the rolling direction, the Lankford coefficient is
(15)
R
θ
=
d
ε
θ
+
9
0
°
d
ε
33
=
-
(
(
∂
f
p
∂
σ
11
+
∂
f
p
∂
σ
22
)
-
1
(
∂
f
p
∂
σ
11
si
n
2
θ
+
∂
f
p
∂
σ
22
cos
⁡
2
θ
-
∂
f
p
∂
σ
12
sin
θ
cos
⁡
θ
)
×
(
∂
f
p
∂
σ
11
+
∂
f
p
∂
σ
22
)
-
1
)
=
(
(
G
p
+
F
p
ta
n
2
θ
)
-
1
(
H
p
+
H
p
ta
n
2
θ
si
n
2
θ
+
2
N
p
si
n
2
θ
-
(
F
p
+
G
p
+
3
H
p
)
si
n
2
θ
)
×
(
G
p
+
F
p
ta
n
2
θ
)
-
1
)
.
Based on the experimental data of R
0
,
R
45 and R
90 as well as the uniaxial compressive yield stress along the through-thickness direction of σ
z, the anisotropic parameters of the Hill48 plastic potential function are obtained:
(16)
F
p
=
R
0
(
R
0
+
R
90
)
σ
z
2
,
G
p
=
R
90
(
R
0
+
R
90
)
σ
z
2
,
H
p
=
R
0
R
90
(
R
0
+
R
90
)
σ
z
2
,
N
p
=
2
R
45
+
1
2
σ
z
2
.
The anisotropic parameters of L
p and M
p are assumed to be identical with the value of N
p. It should be noticed that the difference between the stress-shaped yield function andR-value shaped plastic potential function of Hill 1948 is the type of experimental inputs for parameter identification.


## 3. Experiments

As the billet during hot rolling mainly experiences compressive loading, the compression tests are carried out to calibrate anisotropic parameters in the yield function and non-AFR potential function. The thickness of the original billet is 300 mm and reduces to 12 mm after fifteen rolling passes. The rolling temperature declines from the initial 370°C to 265°C. Specimens of 12 mm× 12 mm × 12 mm are machined and compressed along the planar directions of 0°, 30°, 45°, 60°, 90°, and the vertical direction, as illustrated in Figure 1. Specimens are heated up to 265°C and kept at this temperature for 5 minutes before specimens are tested. The specimens are compressed by the reduction of 40% at the speed of 0.5 mm/minute. Specimens after compression are demonstrated in Figure 2.Figure 1
Sampling schematic.Figure 2
The compressed samples along different directions.The yield stresses and Lankford coefficients at different directions are measured and presented in Table1. It is obvious that the compressive yield stress along the transverse direction is much larger than those along other directions, while theR-value of vertically compressed specimen is a little higher than the others. Then these experimental data in Table 1 are utilized to calibrate the anisotropic parameters in the yield function and plastic potential function of Hill48, respectively with (6)–(10) and (16). T herein the yield stress at the vertical direction σ
z is used to represent the effective yield stress σ
-
y in (1) and σ
-
p in (13). The anisotropic parameters calibrated are summarized in Tables 2 and 3.Table 1
Directional compressive yield stress and Lankford coefficients of hot rolled AZ31.
Direction
0
30
45
60
90
Vertical
Yield stress
61.29
61.96
63.89
66.56
87.46
63.43
R-value
0.4910
0.5087
0.5256
0.5385
0.5613
1.4401Table 2
Material parameters in the Hill48 yield function for hot rolled AZ31.
F
y
G
y
H
y
L
y
M
y
N
y
0.2275
0.7725
0.2985
1.4713
1. 4713
1. 4713Table 3
Material parameters in the Hill48 plastic potential function for hot rolled AZ31.
F
p
G
p
H
p
L
p
M
p
N
p
0.4666
0.5334
0.2619
1.0256
1.0256
1.0256The anisotropies in compressive yield stresses andR-values by the AFR and non-AFR Hill48 models are predicted and compared with experimental results in Figures 3 and 4, respectively. The results demonstrate that the proposed non-AFR Hill48 model provides sufficient predictability of the directionality both in compressive yield stresses andR-values. The advantage of the non-AFR Hill48 model lies in that a more comprehensive consideration of the effect of anisotropy on deformations is correctly introduced.Figure 3
Comparison of the compressive yield stress directionality for hot rolled AZ31.Figure 4
Comparison of theR-value directionality for compressed hot rolled AZ31.Figure5 compares the shape of plastic potential and yield function of Hill48 at zero vertical stress. It is observed that the potential function is much different from the yield surface and the non-AFR Hill48 yield function is accurate enough to describe the anisotropy in both the yield stress and plastic flow.Figure 5
Comparison of yield function and plastic potential of Hill48 inx
-
y stress plane.The stress-plastic strain curves are compared for the uniaxial compressions along different directions in Figure6. The comparison obviously indicates that the compressive stress along the rolling direction is much larger than those along other directions. To be consistent with the yield stress at the vertical direction σ
z which is used for parameters identification above, the strain-hardening behavior at the vertical direction is incorporated into the user-defined material subroutine VUMAT. Unlike common exponential or power forms, there is a concave transition between yielding stress and the maximum compressive stress on each stress-plastic strain curve; thus in this paper a quintic polynomial is utilized to describe the strain-hardening behavior:
(17)
σ
-
=
65360
*
ε
-
p
5
+
(
-
72170
)
*
ε
-
p
4
+
29880
*
ε
-
p
3
+
(
-
5325
)
*
ε
-
p
2
+
398.5
*
ε
-
p
+
52.27
,
where σ
- and ε
-
p are the equivalent stress and the equivalent plastic strain, respectively. The coefficients in (17) are calibrated by experimental stress-plastic strain curve from uniaxial compression test along the vertical direction.Figure 6
Stress-plastic strain curves of compressed AZ31 at different directions.


## 4. Numerical Integration Algorithm

Finite element implementation of the constitutive model requires numerical integration of the constitutive equations over the time increment. Here the semi-implicit backward Euler method proposed by Moran et al. [18] is adopted. The dominant feature of the semi-implicit backward Euler method is that it is implicit in plasticity parameter and explicit in the plastic flow direction and plastic moduli. That means the increments in the plasticity parameter are calculated at the end of the step while the plastic flow direction and plastic moduli are calculated at the beginning of the step. To avoid drift from the yield surface, the yield condition is enforced at the end of the step. The integration scheme is summarized as follows [19].(1)
Give the set(
σ
t
,
ε
t
,
ε
-
t
p
)  at time t and the strain increment Δ
ε
t
+
Δ
t.
(2)
Suppose the strain increment at the current time step is purely elastic, a trial stress or an elastic predictor is calculated:(18)
σ
t
+
Δ
t
trial
=
σ
t
+
C
:
Δ
ε
t
+
Δ
t
.
(3)
Substituteσ
t
+
Δ
t
trial into (1) to check the yield condition:
(19)
f
t
+
Δ
t
e
=
f
y
(
σ
t
+
Δ
t
trial
,
ε
-
t
p
)
.
If f
t
+
Δ
t
e
≤
0, it means an elastic process; then σ
t
+
Δ
t
=
σ
t
+
Δ
t
trial, ε
-
t
+
Δ
t
p
=
ε
-
t
p; else, the yield surface is expanded and the following plastic correction is needed.
(4)
Plastic correction:
(a)
Initialization: set the initial equivalent plastic strain to the converged values at the end of previous time-step, zero the increment in plasticity parameterλ, and evaluate the elastic trial stress:
(20)
k
=
0
:
ε
-
(
0
)
=
ε
-
t
p
,
Δ
λ
(
0
)
=
0
,
σ
(
0
)
=
σ
t
+
Δ
t
trial
-
Δ
λ
(
0
)
C
:
r
t
,
where r
=
∂
f
p
/
∂
σ represents the plastic flow direction and Δ
λ
(
0
)
C
:
r
t is the plastic corrector.
(b)
Check yield condition and convergence at thekth iteration:
(21)
f
y
(
k
)
=
f
y
(
σ
(
k
)
,
ε
-
p
(
k
)
)
.
If f
y
(
k
)
<
TOL, converged, return (1); else, go to (c).
(c)
Compute the increment in plasticity parameter:(22)
δ
λ
(
k
)
=
f
(
k
)
∂
f
(
k
)
:
A
(
k
)
:
r
~
t
,
where
(23)
[
A
(
k
)
]
=
[
C
0
0
-
1
]
(
k
)
,
{
r
~
t
}
=
{
r
t
1
}
,
[
∂
f
(
k
)
]
=
[
∂
f
(
k
)
∂
σ
i
j
∂
f
(
k
)
∂
ε
-
p
]
.
(d)
Obtain the increments in stress and equivalent plastic strain:(24)
{
Δ
σ
(
k
)
Δ
ε
-
p
(
k
)
}
=
-
δ
λ
(
k
)
[
A
(
k
)
]
{
r
~
t
}
.
(e)
Update the stress and equivalent plastic strain:(25)
ε
-
p
(
k
+
1
)
=
ε
-
p
(
k
)
+
Δ
ε
-
p
(
k
)
,
Δ
λ
(
k
+
1
)
=
Δ
λ
(
k
)
+
δ
λ
(
k
)
,
σ
(
k
+
1
)
=
σ
(
k
)
+
Δ
σ
(
k
)
.
If k
←
k
+
1, go to (b).The experimental and simulated deformations of the sample after uniaxial compression are compared in Figure7 to validate the anisotropic constitutive model built according to the above theories and identified parameters. The Hill48 yield function under AFR is also used to simulate these deformation processes for the comparison purpose. It is obvious that the dimensions of all the three directional samples simulated by the non-AFR Hill48 model are in much better agreement with the experimental data than those predicted by the Hill48 yield function under AFR.Figure 7
Comparison between experimental and numerical flattening of compressed samples.


## 5. Application to the Hot Rolling Simulation

### 5.1. Square-Bar Rolling Pass

The calibrated Hill48 yield function under non-AFR is applied to simulate the elastoplastic deformation of AZ31 magnesium alloy during hot rolling process. For the comparison purpose, the von Mises yield function and the Hill48 yield function under AFR are also used in the simulation of square-bar hot rolling process. As shown in Figure8, the square-bar billet with a cross-section of 20 mm × 20 mm is rolled by a thickness reduction of 40%. The elastic modulus of AZ31 is 44.9 GPa while the roller is assumed to be rigid. Rolling velocity is set to 1 m/s.Figure 8
FE model adopted in the square-bar rolling simulation.For simulation results predicted by both the isotropic and anisotropic material models, the feeding end of the billet is warped to some degrees as a result of initial interaction with the roller. However, it has little effect on the flatness of the remainder part. Besides, the lateral displacement near both ends is larger than that at the middle of the billet as depicted in Figure9. Thus, the contour at the middle of the billet after rolling in Figure 9 is extracted from simulations with three plasticity models as compared in Figure 10. The comparison apparently reveals that the lateral spread at the middle of the billet is strongly affected by anisotropy of the material. It is observed that the lateral spread at the middle of the billet is 21% for the isotropic von Mises model, 39.9% for the AFR Hill48 function, and for the non-AFR Hill48 model the predicted lateral spread is 28.7%, which is between those predicted by the von Mises and Hill48 under AFR models.Figure 9
Comparison of the square-bar billet’s lateral displacements for the three material models.Figure 10
Comparison of the cross-section profiles at the middle of the rolled square-bar billet.


### 5.2. Strip Rolling Pass

Simulations and experiments of a strip rolling pass are carried out to evaluate the predictability of the non-AFR Hill48 model in predicting lateral spread. The rolling is operated on a reversible four-roller hot rolling mill for magnesium alloy plates. The hot rolling mill consists of two symmetrically arranged backup rollers and two working rollers as illustrated in Figure11 (left). The billet thickness in experiments is identical with that of compression specimens so that the calibrated yield function in Section 3 can be utilized directly in simulation of strip rolling pass. The finite element model is described in Figure 11 (right) for an AZ31 strip billet with the 12 mm × 100 mm cross-section.Figure 11
The four-roller hot rolling mill for magnesium alloys and the FE model for strip rolling.The experimental section at the middle of the deformed billet is compared with those predicted by the isotropic von Mises yield function, the Hill48 yield function under AFR, and the non-AFR Hill48 yield in Figure12. The comparison demonstrates that the lateral spread of strips is severely underestimated by the isotropic von Mises yield function, while the Hill48 yield function under AFR overestimates the lateral spread. The Hill48 yield function under non-AFR predicts the lateral spread of the strip after rolling with very high accuracy.Figure 12
Comparison of the cross-section profiles at the middle of the rolled strip billet.


## 5.1. Square-Bar Rolling Pass

The calibrated Hill48 yield function under non-AFR is applied to simulate the elastoplastic deformation of AZ31 magnesium alloy during hot rolling process. For the comparison purpose, the von Mises yield function and the Hill48 yield function under AFR are also used in the simulation of square-bar hot rolling process. As shown in Figure8, the square-bar billet with a cross-section of 20 mm × 20 mm is rolled by a thickness reduction of 40%. The elastic modulus of AZ31 is 44.9 GPa while the roller is assumed to be rigid. Rolling velocity is set to 1 m/s.Figure 8
FE model adopted in the square-bar rolling simulation.For simulation results predicted by both the isotropic and anisotropic material models, the feeding end of the billet is warped to some degrees as a result of initial interaction with the roller. However, it has little effect on the flatness of the remainder part. Besides, the lateral displacement near both ends is larger than that at the middle of the billet as depicted in Figure9. Thus, the contour at the middle of the billet after rolling in Figure 9 is extracted from simulations with three plasticity models as compared in Figure 10. The comparison apparently reveals that the lateral spread at the middle of the billet is strongly affected by anisotropy of the material. It is observed that the lateral spread at the middle of the billet is 21% for the isotropic von Mises model, 39.9% for the AFR Hill48 function, and for the non-AFR Hill48 model the predicted lateral spread is 28.7%, which is between those predicted by the von Mises and Hill48 under AFR models.Figure 9
Comparison of the square-bar billet’s lateral displacements for the three material models.Figure 10
Comparison of the cross-section profiles at the middle of the rolled square-bar billet.


## 5.2. Strip Rolling Pass

Simulations and experiments of a strip rolling pass are carried out to evaluate the predictability of the non-AFR Hill48 model in predicting lateral spread. The rolling is operated on a reversible four-roller hot rolling mill for magnesium alloy plates. The hot rolling mill consists of two symmetrically arranged backup rollers and two working rollers as illustrated in Figure11 (left). The billet thickness in experiments is identical with that of compression specimens so that the calibrated yield function in Section 3 can be utilized directly in simulation of strip rolling pass. The finite element model is described in Figure 11 (right) for an AZ31 strip billet with the 12 mm × 100 mm cross-section.Figure 11
The four-roller hot rolling mill for magnesium alloys and the FE model for strip rolling.The experimental section at the middle of the deformed billet is compared with those predicted by the isotropic von Mises yield function, the Hill48 yield function under AFR, and the non-AFR Hill48 yield in Figure12. The comparison demonstrates that the lateral spread of strips is severely underestimated by the isotropic von Mises yield function, while the Hill48 yield function under AFR overestimates the lateral spread. The Hill48 yield function under non-AFR predicts the lateral spread of the strip after rolling with very high accuracy.Figure 12
Comparison of the cross-section profiles at the middle of the rolled strip billet.


## 6. Conclusions

The Hill48 yield function is combined with a non-AFR in a Hill48 form to describe the plastic behavior of AZ31 magnesium alloy. The non-AFR Hill48 yield function is calibrated by uniaxial compression tests since the proposed plasticity model is applied to hot rolling process during which materials mainly experience compression. The calibrated non-AFR Hill48 yield function was implemented in a user-defined material subroutine of ABAQUS/Explicit based on a semi-implicit backward Euler method. The user subroutine is utilized to simulate plastic deformation of uniaxial compression tests to verify the correct implementation of the non-AFR Hill48 yield function. The implemented non-AFR Hill48 yield function is utilized to simulate both square-bar rolling pass and strip rolling pass. The simulated lateral spreads are compared with experimental results and those predicted by the von Mises yield function and the Hill48 yield function under AFR. The comparison validates that the non-AFR Hill48 yield function calibrated by uniaxial compression tests can accurately model the plastic behavior of AZ31 magnesium alloy in hot rolling processes. Hence, the proposed non-AFR Hill48 yield function is suggested to be calibrated by the uniaxial compression tests for the accurate simulation of various rolling processes in design of rolling tools and factor optimization in rolling.



---

*Source: 256194-2014-03-13.xml*