ZK60 alloy deformation and damage at different temperatures using unified constitutive model

Abstract

The deformation behaviour and damage mechanism of ZK60 magnesium alloy under temperatures of 25–140 °C and strain rates of 0.01–0.0001 s⁻¹ were investigated by uniaxial hot tensile and characterization experiments. The results suggested that the main factor for damage in ZK60 magnesium alloy was the nucleation, growth and coalescence of microdefects such as microvoids and cracks around the second-phase particles and/or grain boundaries. Then, a unified constitutive model including damage evolution, work hardening and recovery mechanisms was established. The model, verified by experimental data, shows good predictability. It was implemented into ABAQUS software through the VUMAT subroutine, which then predicted and analyzed the microstructure evolution.

Keywords

ZK60 magnesium alloy deformation behaviour damage evolution unified constitutive model

Introduction

ZK60 magnesium alloy has excellent performance of low density and high strength and is widely used in automotive, aviation and aerospace. In the process of heat treatment to improve the formability and mechanical performance, the ZK60 magnesium alloy complex component undergoes dynamic recrystallization (DRX), dislocation slip, twinning and other complex deformation mechanisms,¹ which also cause microcracks that lead to material damage and ductile fracture^2–4 severely limiting its application. However, the attention and research on the damage mechanism of low-symmetry hexagonal close-packed magnesium alloys are still insufficient, one reason is that the fracture stress-strain relationship of face-centred cubic metals (such as aluminum alloys and dual phase steels) is not applicable to magnesium alloys.⁵ Currently, the demand for lightweight components in military and civil fields is huge, which makes it a key issue for manufacturing development to elucidate the deformation behaviour and damage evolution mechanism of ZK60 magnesium alloy at different temperatures.

At present, there have been some research on the evolution of the flow behaviour and microstructures in ZK60 magnesium alloy deformation. Hadadzadeh and Wells⁶ conducted Gleeble thermal simulation tests on cast and extruded (along with the squeezing direction and cross-sectional direction) ZK60 Magnesium alloy and established a constitutive model using hyperbolic sine constitutive equation and Ludwig equation. They analyzed the constitutive parameters, activation energy, strain hardening and strain rate sensitivity of the material. Xu et al.⁷ conducted in-situ Digital Image Correlation (DIC) tensile tests on ZK60 magnesium alloy rolling sheet and observed the micro damage volume fraction and morphology using X-ray computed tomography and scanning electron microscope (SEM).

In order to better describe the damage behaviour, Cockcroft⁸ proposed an empirical damage criterion (C-L model) based on energy accumulation theory, which assumes that the maximum tensile stress causes fracture. Similarly, Johnson-Cook (J-C model),⁹ McClintock,¹⁰ Oyane¹¹ and Rice and Tracey¹² successively proposed the criteria of damage. Subsequently, Gurson¹³ proposed the Gurson model based on the micromechanics of voids and established the plastic potential equation to describe the damage void effect, but did not consider the inter-action between the voids. Later, Tvergaard and Needleman¹⁴ replaced the void volume fraction with a parameter to enhance the prediction accuracy.

The most critical problem in the damage model is the coupling of damage and plastic deformation, so the model based on continuum damage mechanics (CDM) has been widely developed. Hayhurst et al.¹⁵ studied the failure of welding parts at high temperature and established a relation between the internal variables and the effective strain rate of damage. Lin et al.¹⁶ proposed a unified constitutive model to describe the formation and development of damage. Huo et al.¹⁷ established a multi-axial model of microstructure and ductile damage of a high-speed railway axle steel during cross wedge rolling. This model can effectively predict the evolution of grain size and ductile damage. Xiao et al.¹⁸ analyzed the flow behaviour and fracture damage of AA7075 aluminum alloy by establishing a DRX, grain size and damage evolution coupling. Using high-temperature tensile tests, Huo et al.¹⁹ examined the flow stress, microstructural evolution and fracture morphology of forged TC4 titanium alloy. They introduced internal state variables such as DRX volume fraction, grain size and stress triaxiality to improve the GTN damage model. Liu et al.²⁰ analyzed the evolution of damage during the decompression preprocessing process based on tensile tests. For magnesium alloys, Feng et al.²¹ studied the tensile and impact test of AZ31B magnesium alloy at different strain rates and temperatures and established a fracture model based on the J-C model. Kim et al.²² proposes a fracture prediction method for AZ31B magnesium alloy sheet based on a micro-mechanical void model and an asymmetric plasticity constitutive law, which considers the material anisotropy and damage evolution. Zhang et al.²³ developed a fully coupled damage model which is based on the classical CDM approach, for Mg alloy sheet forming simulations at elevated temperatures. Xu et al.²⁴ used a modified GTN model with the representative volume element method to investigates the anisotropic damage mechanism of the LA103Z Mg-Li alloy rolling sheet.

Currently, the constitutive models of the hot deformation behaviour of ZK60 magnesium alloy, such as the hyperbolic sine constitutive equation and the Ludwig equation, usually only consider the effects of temperature, strain and strain rate and other factors. Internal variables such as dislocation density are often neglected in material deformation and damage. It has been found that the nucleation, growth and coalescence of microdefects such as voids and cracks are the important causes of damage and fracture in magnesium alloys,²⁵ but few constitutive models have been reported that consider the damage evolution of microdefects in ZK60 magnesium alloy. This limits the ability of existing constitutive models to describe the flow behaviour of ZK60 magnesium alloy at different temperatures and makes it difficult to fully capture the variation of flow stress with strain. Therefore, it is of great significance to establish a constitutive model that can accurately describe the deformation behaviour and damage evolution of ZK60 magnesium alloy at different temperatures, for optimizing the hot processing parameters and controlling the microstructure of the material.

In this study, it is assumed that the main factor for thermal deformation damage is the microvoids and cracks generated around the second-phase particles and/or grain boundaries, which coalesce to form macrocracks. The relative importance of these two factors for damage varies with strain rate, temperature and dislocation density. The aim of this paper is to propose a unified damage constitutive model based on coupling dislocation density and hardening rate, which describes the relationship between flow stress, dislocation density and damage rate in ZK60 magnesium alloy during hot processing. The accuracy and applicability of the model are verified by experimental data, and then the model is implemented into ABAQUS through the VUMAT subroutine, to predict the evolution of internal variables during uniaxial tensile process.

Experimental materials and procedures

The experiment used hot-rolled ZK60 magnesium alloy with a thickness of 2 mm, which mainly consists of 4.8–6.2% zinc and zirconium with a content greater than 0.45%.

Uniaxial tensile tests of ZK60 magnesium alloy at different temperatures were carried out on a DDL-100 high-low temperature tensile testing machine. The uniaxial tensile specimens were cut from the small-sized rolled plate along the rolling direction, where the width direction of the specimens was the transverse direction of the plate, and the thickness direction was the normal direction of the plate. The gauge section size of the specimens was 50 mm × 8 mm × 2 mm, as shown in Figure 1. The surface of the specimens was polished with 400–1000 grit sandpaper. The specimens were ultrasonically cleaned with acetone solvent. The deformation temperatures were room temperature, 60 °C, 100 °C and 140 °C, respectively, and the initial strain rates were 0.01 s⁻¹, 0.001 s⁻¹ and 0.0001 s⁻¹, respectively. The experimental process was carried out according to the national standard GB/T 228.2–2015. The experimental procedure is shown in Figure 1.

Polish the specimens with sandpaper and use a polishing machine to achieve a mirror finish. The specimens were etched using the mixed solution (1.5 g picric acid +1.25 mL acetic acid + 25 mL absolute ethanol +5 mL deionized water) at 25 °C for 35 s. Clean the etched surface with ethanol and dry it with an air dryer. Microstructures of specimens were observed by using a DM4000 M optical microscope (OM) and a TESCAN S8000 SEM.

Analysis of the deformation behaviour of ZK60 magnesium alloy

ZK60 Magnesium alloy's uniaxial tensile process true stress–strain response is shown in Figure 2. Plastic deformation depends on temperature and strain rate, which is a typical viscoelastic deformation.

The effect of different temperature on deformation behaviour

The plasticity of the alloy varies significantly with different temperatures. From Figure 3(a), it can be seen that when the strain rate is 0.01 s⁻¹, as the temperature increases, the fracture strain have increased, indicating that its plastic deformation ability increases. The hot formability of ZK60 magnesium alloy increases with the rise of temperature. Therefore, choosing hot stamping process for parts with complex shapes and large stamping depths is beneficial for material formability.

It can be seen in Figure 3(c) and (d) that under the same temperature and strain rate conditions, the ZK60 samples occur in work hardening and dynamic softening during thermal deformation. As shown in Figure 3(c), under the conditions of 140 °C and 0.01 s⁻¹, when the strain is less than 0.05, the work hardening rate gradually decreases, but the stress–strain curve is still on the rise, indicating that the hardening effect of the material at this time is greater than the softening effect. The material enters a significant plastic stage when the strain exceeds 0.05. When the strain is between 0.05 and 0.27, the work hardening and stress–strain curve are basically stability, indicating that the softening effect of the material at this time is enhanced due to dynamic recovery, and the hardening effect reaches a balanced state. Figure 3(d) compares the effect of different temperatures on the work hardening rate at a strain rate of 0.01 s⁻¹. It shows that the softening effect of the material at this time is enhanced due to dynamic recovery, and the effect of temperature on dynamic recovery is more obvious; when the stress reaches the peak, as the strain continues to increase, the work hardening rate has begun to decline, and negative values will appear until the fractures occurs. As illustrated in Figure 3(b), the max stress exhibits a negative correlation with the temperature, indicating a thermally induced softening behaviour.

The effect of different strain rates on deformation behaviour

The fracture strain in relation to the temperature and strain rate is depicted in Figure 4(a). It is evident that the fracture strain decreases as the strain rate increases.

The plastic deformation of magnesium alloy under warm forming mainly relies on twinning mechanism and dislocations.²⁶ Twinning leads to stress concentration and crack initiation, reducing the fracture strain of magnesium alloy.²⁷ The flow stress of magnesium alloy under warm forming increases with the increase of strain rate, and high flow stress promotes the occurrence and expansion of twinning, further reducing the fracture strain of magnesium alloy. High strain rate reduces the time for dynamic recovery, which weakens the softening effect, enhances the deformation resistance of the material and leads to earlier fracture. In summary, within the temperature range of warm forming, the plasticity is not sensitive to strain rate, which is conducive to forming complex parts at higher forming speed, to improve productivity and reduce heat loss.

Figure 4(b) indicates that when ZK60 magnesium alloy is deformed at 140 °C, work hardening becomes obvious with the increase of strain rate, which is mainly due to less static recovery occurring during deformation at higher strain rate under the same temperature. During hot working, excessive strain hardening may reduce the ductility of the material and increase the risk of brittle fracture during plastic deformation.

The effect of strain rate on flow stress at a strain of 0.05 is shown in Figure 4(c). The true stress–strain rate curves obtained from experiments were fitted, indicating that the viscoplastic response of the material follows power law. With the increase of strain rate, flow stress increases, showing strain rate hardening characteristics. The strain rate hardening exponent is represented by the slope ‘m’ of the trend line, corresponding to the power exponent in the equation. This is a key factor affecting the deformation uniformity in hot stamping process. At 140 °C, the m value of magnesium alloy is slightly larger than that at 100 °C. Within a experimental temperature range, temperature has little effect on strain rate hardening.

Microstructure characterization of damage failure process

Fracture analysis

The fracture morphologies of the deformed specimens under different conditions are illustrated in Figure 5. In the damage process, the fracture surfaces of the specimens are usually covered with dimples, voids and microcracks. When the material undergoes large plastic deformation, the dimples and voids grow, and those with large enough diameters tend to coalesce and form microcracks, which eventually lead to material failure.²⁸ The deformation parameters, such as strain rate and deformation temperature, have significant effects on the nucleation and growth of dimples and voids. At lower deformation temperature and lower strain rate, as shown in Figure 5(a), there are fewer dimples and voids, and some of the microcracks have larger diameters. Combined with the study of Xiao et al.²⁹ on the fracture of ZK60 magnesium alloy, it can be found that with the increase of strain rate, the number of dimples gradually increases, and the dimples become shallow and stop growing. Combined with the study of Tang et al.,³⁰ it can be found that when the deformation temperature increases, the number of dimples decreases, but the dimples become deeper, and the crack propagation accelerates. Therefore, with the increase of strain rate and the decrease of deformation temperature, the nucleation rate of dimples, the coalescence rate of dimples and the expansion rate of cracks all increase.

Damage mechanism analysis

On the cross-section close to the fracture surface of the specimens after deformation, the microvoids and microcracks are revealed in Figure 6(b). The cracks in ZK60 magnesium alloy first initiate at the voids, which are distributed at the polycrystalline interaction and perpendicular to the grain boundaries. During tensile process, the voids are elongated and expanded, and connected with each other, eventually leading to crack formation. In addition, the cracks propagate along with the shortest path through adjacent grain boundaries, or intersect with other cracks to form cracks. The growth, coalescence and crossing directions of voids and cracks are consistent with the tensile direction.

The initial microstructure of the material was observed by OM, as shown in Figure 6(a). The grains were stretched along RD direction, and the grain distribution showed a strong rolling texture. To further observe the morphology of the second-phase particles, SEM electron micrographs were taken, as shown in Figure 6(b). The results show that the second-phase particles are widely distributed at the grain boundaries, but the second-phase particles at the grain boundaries are significantly larger than those in the grains. The length and width of the second-phase particles at the grain boundaries are about 200 nm, much smaller than the grain size of 1 μm, but much larger than the size of the second-phase particles in the grains. Based on these observations, the possible fatigue damage mechanism is controlled by the microcrack initiation at the second-phase particles at the grain boundaries, which are likely to coalesce and lead to final failure. This is similar to the results of Jordon et al.³¹

Unified constitutive model

Formulation of unified constitutive model

Damage evolution

The damage value is defined as the ratio of the sum of the areas of microvoids and cracks in the typical cross-section of the material to the total area of the cross-section, expressed as follows³²:

D = \frac{A - \bar{A}}{A}

(1)

where A is the total area of the typical cross-section;

\bar{A}

is the effective area of the typical cross-section (

0 \leq \bar{A} \leq A

). Studies have shown that when D reaches the threshold of 0.7, the material can be considered as failed.³³

This paper studies the strain rate effect on damage evolution using the improved GTN model and accounts for damage by microcrack and microvoid formation and coalesced crack growth. The plastic damage model is:

{\dot{D}}_{T} = {\dot{D}}_{N} + {\dot{D}}_{G}

(2)

where

{\dot{D}}_{T}

is the total plastic damage growth rate;

{\dot{D}}_{N}

is the plastic damage nucleation rate;

{\dot{D}}_{G}

is the plastic damage growth rate.

Yang et al.³⁴ suggested that the void nucleation rate of TC6 alloy is a function of the strain rate difference between $β$ phase and $α$ phase, which is suppressed by the total damage value, and is related to temperature. Considering the effect of recovery on reducing the dislocation density during metal thermal deformation, the plastic damage nucleation rate is modified as follows:

{\dot{D}}_{N} = π (1 - D_{N}) \dot{\bar{ρ}}

(3)

where

π

is a temperature-dependent constant.

Considering the damage mechanism of alloy uniaxial tensile deformation, it is assumed that the damage evolution rate can be expressed by the plastic strain rate and the total plastic damage value as follows:

{\dot{D}}_{G} = f (D_{T}) β_{d} η {\dot{ε}}_{p} = \frac{β_{d} {| {\dot{ε}}^{p} |}^{δ} \cosh (γ_{4} ε^{p})}{{(1 - D_{T})}^{γ_{3}}}

(4)

where

β_{d}

δ

γ_{3}

γ_{4}

are material constants.

Equation (4) represents the effect of plastic strain on the damage change rate, where the four material constants control the rate of microcrack growth and coalescence under different deformation conditions, reflecting the strength of the influence of deformation conditions on damage evolution. The introduction of $(1 - D_{T})^{- γ_{3}}$ describes the rapid increase of damage value at the eve of material fracture (after necking).

Plastic flow law

According to the classical Norton equation, the flow stress equation is given by:

{\dot{ε}}_{p} = {⟨ \frac{σ - H - k}{K} ⟩}^{n}

(5)

where

{\dot{ε}}_{p}

is the plastic strain rate;

σ

is the total flow stress; k is the initial yield stress. K is the resistance coefficient; n is the viscous exponent;

⟨ \cdot ⟩

is the McCauley bracket.

The effective stress after removing the damage in the plastic deformation process of the material should be:

\bar{σ} = \frac{σ}{1 - D_{T}}

(6)

Substituting the effective stress into equation (5), the modified flow equation is obtained as follows:

{\dot{ε}}^{p} = {⟨ \frac{\frac{σ}{1 - D_{T}} - H - k}{K} ⟩}^{n} \frac{1}{{(1 - D_{T})}^{γ_{1}}}

(7)

where

γ_{1}

are material constants.

Equation (7) introduces $(1 - D_{T})^{- γ_{1}}$ to regulate the plastic deformation process of different metal materials, making the equation more applicable.

The deformation process of metal materials consists of two stages: elastic deformation and plastic deformation. The relationship between elastic strain and flow stress is obtained by Hooke's law as follows:

σ = E (1 - D_{T}) (ε_{T} - ε_{p})

(8)

where E is the elastic modulus;

ε_{T}

is the total strain;

ε_{p}

is the plastic strain.

Dislocation density evolution

Considering the effect of recovery on reducing dislocations, Lin et al.³⁵ introduced a rate equation for normalized dislocations as follows:

\dot{\bar{ρ}} = A (1 - \bar{ρ}) | {\dot{ε}}_{p} | - C {\bar{ρ}}^{n_{2}}

(9)

where A, C, n₂ are material constants. The above equation can be rewritten as:

\dot{\bar{ρ}} = A | {\dot{ε}}_{p} | - A \bar{ρ} | {\dot{ε}}_{p} | - C {\bar{ρ}}^{n_{2}}

(10)

The equation has three terms: the first item represents the accumulated dislocations due to plastic deformation, and parameter A is temperature dependent; the second term represents the decrease in dislocation caused by dynamic recovery, which bounds the normalized dislocation density to 1; the third term represents the decrease in dislocation caused by static recovery, which can be ignored for cold forming, with C as another temperature-dependent parameter.

Hardening is directly related to dislocation density, which can be expressed by a dislocation function:

H = B {\bar{ρ}}^{\frac{1}{2}}

(11)

Rewrite the above equation as an evolution equation, that is

\dot{H} = 0.5 B {\bar{ρ}}^{- 0.5} \dot{\bar{ρ}}

(12)

where B is the material constant.

Determination of material constants

Based on the above analysis, the unified damage constitutive equation for different temperatures can be written as follows:

{\begin{matrix} {\dot{ε}}^{p} = {⟨ \frac{\frac{σ}{1 - D_{T}} - H - k}{K} ⟩}^{n} \frac{1}{{(1 - D_{T})}^{γ_{1}}} \\ \dot{H} = 0.5 B {\bar{ρ}}^{- 0.5} \dot{\bar{ρ}} \\ \dot{\bar{ρ}} = A (1 - \bar{ρ}) | {\dot{ε}}^{p} | - C {\bar{ρ}}^{γ_{2}} \\ \begin{matrix} {\dot{D}}_{N} = π (1 - D_{N}) \dot{\bar{ρ}} \\ {\dot{D}}_{G} = \frac{β_{d} {| {\dot{ε}}^{p} |}^{δ} \cosh (γ_{4} ε^{p})}{{(1 - D_{T})}^{γ_{3}}} \\ {\dot{D}}_{T} = {\dot{D}}_{N} + {\dot{D}}_{G} \\ σ = E (1 - D_{T}) (ε^{t} - ε^{p}) \end{matrix} \end{matrix}

(13)

There are 24 constants to be determined in the constitutive model, and the temperature-dependent material constants and their expressions are listed in Table 1. Where T is the absolute temperature, R is the gas constant, Q, Qn, Qc, Q_E and Q_D are activation energy-related constants.

Table 1.

List of temperature-dependent parameters.

$n = n_{0} \exp (Q_{n} / R T)$	$π = π_{0} \exp (Q_{J} / R T)$
$k = k_{0} \exp (Q_{k} / R T)$	$β_{d} = β_{0} \exp (Q_{D} / R T)$
$K = K_{0} \exp (Q_{K} / R T)$	$γ_{3} = γ_{30} \exp (Q_{γ 3} / R T)$
$B = B_{0} \exp (Q_{B} / R T)$	$γ_{4} = γ_{40} \exp (Q_{γ 4} / R T)$
$C = C_{0} \exp (- Q_{C} / R T)$	$E = E_{0} \exp (Q_{E} / R T)$

The remaining parameters were determined based on the experimental stress–strain flow curve, and the forward Euler iteration was carried out in three steps using Matlab software.

Step 1:

For each equation that has been constructed, the initial value of all variables given in time is 0. According to the integral solution process showing Euler's method, there are:

\begin{aligned} ε_{i + 1}^{p} = ε_{i}^{p} + {\dot{ε}}_{i}^{p} Δ t \\ ε_{i + 1}^{t} = ε_{i}^{t} + {\dot{ε}}_{i}^{t} Δ t \\ K_{1, i + 1} = K_{1, i} + {\dot{K}}_{1, i} Δ t \\ K_{2, i + 1} = K_{2, i} + {\dot{K}}_{2, i} Δ t \\ \dots \\ K_{i + 1} = K_{1, i + 1} + K_{2, i + 1} + \dots \\ σ_{i + 1} = E (ε_{i + 1}^{t} - ε_{i + 1}^{p}) \\ σ_{i + 1} = σ_{i} + Δ t \end{aligned}

(14)

In equation (14),

Δ t

is the time increment;

i = 1, 2, 3, \dots, n

, where i is the number of steps of the current integral iteration and

i + 1

is the number of steps of the next integral iteration based on the current cumulative value.

Step 2:

The material constant is calculated by the GA algorithm, which minimizes the logarithmic error between the predicted and experimental values. The objective function is the residual difference of the flow stress and fracture strain:

{\begin{matrix} f_{σ} (X) = \frac{1}{M} \sum_{j = 1}^{M} [\frac{1}{N_{j}} \sum_{i = 1}^{N_{j}} {(\ln \frac{σ_{i j}^{c}}{σ_{i j}^{e}})}^{2}] \\ f_{D} (X) = \frac{1}{M} \sum_{j = 1}^{M} {(\ln \frac{D_{j}^{c}}{D_{j}^{e}})}^{2} \end{matrix}

(15)

where X—Material constant vector set, $X = {X_{1}, X_{2}, X_{3}, \dots, X_{n}}$ ;

$f_{D}$ —Residual difference between experimental value and predicted value of true stress;

$f_{d}$ —Residual difference between experimental and predicted fracture strain;

M—Number of true stress–strain curves, $M = 3$ ;

$N_{j}$ —Article j True stress sampling points taken on the true stress–true strain curve;

$L_{j}$ —Article j Sampling number of fracture strain variable taken on the true stress–strain curve;

$σ_{i j}^{c}$ —The predicted value of true stress at the i sampling point on the j curve;

$σ_{i j}^{e}$ —The experimental value of the true stress at the i sampling point on the j curve;

$D_{j}^{c}$ —Fracture strain predicted value on the j curve;

$D_{j}^{e}$ —Experimental value of fracture strain on the j curve.

The global objective function is the sum of various objective functions:

f (X) = f_{σ} (X) + f_{d} (X)

(16)

The material constants are given by Lin and Yang³⁶ based on their physical meaning. The GA parameters are: population size 200, crossover probability 0.8, mutation probability 0.05, solution accuracy 10⁻⁶, iteration number 3000.

Step 3:

Simplify the process further, modify the GA function for different deformation conditions, obtain 14 material constants, fix some non-temperature sensitive constants and use linear fitting for others. The optimized constants are in Table 2.

Table 2.

Determined values of material constants.

Material constants	Optimal value	Material constants	Optimal value	Material constants	Optimal value
$k_{0}$ (MPa)	6.345 × 10⁻¹	$C_{0}$	2.928 × 10⁴	$γ_{40}$	2.259 × 10⁰
$Q_{k}$ (J/mol)	5.109 × 10²	$Q_{C}$ (J/mol)	5.529 × 10⁴	$Q_{γ 4}$ (J/mol)	1.014 × 10⁴
$K_{0}$	2.912 × 10²	$π_{0}$	1.021 × 10⁻³	$E_{0}$ (MPa)	7.017 × 10³
$Q_{K}$ (J/mol)	1.268 × 10²	$Q_{J}$ (J/mol)	3.851 × 10²	$Q_{E}$ (J/mol)	1.180 × 10³
$n_{0}$	1.353 × 10⁻¹	$β_{0}$	1.001 × 10⁻⁶	A	1.060 × 10¹
$Q_{n}$ (J/mol)	1.486 × 10⁴	$Q_{D}$ (J/mol)	1.993 × 10⁴	$γ_{1}$	9.917 × 10⁻³
$B_{0}$ (MPa)	1.138 × 10¹	$γ_{30}$	1.001 × 10⁻²	$γ_{2}$	1.536 × 10⁰
$Q_{B}$ (J/mol)	5.260 × 10³	$Q_{γ 3}$ (J/mol)	1.255 × 10⁴	$δ$	1.360 × 10⁰

Verification of constitutive model

To evaluate the predictability of the proposed model, the correlation coefficient (R), average absolute relative error (AARE) and root mean square error (RMSE) are specified:

R = \frac{\sum_{i = 1}^{N} (E_{i} - \bar{E}) (P_{i} - \bar{P})}{\sqrt{\sum_{i = 1}^{N} {(E_{i} - \bar{E})}^{2} \sum_{i = 1}^{N} {(P_{i} - \bar{P})}^{2}}}

(17)

A A R E = \frac{1}{N} \sum_{i = 1}^{N} | \frac{E_{i} - P_{i}}{E_{i}} | \times 100 %

(18)

RMSE = \sqrt{\frac{1}{N} \sum_{i = 1}^{N} {(E_{i} - P_{i})}^{2}}

(19)

where N is the total number of data points;

E_{i}

is the experimental value;

\bar{E}

is the mean value of experimental values;

P_{i}

is the predicted value;

\bar{P}

is the mean of the predicted values.

The unified constitutive model was validated by comparing numerical and experimental data (Figure 7(a)-(f)). The model shows high accuracy in predicting stress (R = 0.910, AARE = 4.68%, RMSE = 28.1) and low correlation in predicting fracture strain (Figure 7(f)). However, this does not mean that the constitutive model has a worse prediction ability, and the model fits well with the flow stress curves of ZK60 magnesium alloy at different temperatures and reflects the work hardening and dynamic recovery effects. Meanwhile, errors in some experimental measurement processes may lead to prediction errors.

The fracture strain values of the model are the strain values for D = 0.7. The model showed great correlation (R = 0.968) and accuracy (AARE = 9.63%, RMSE = 0.0264) in predicting fracture strain.

The unified constitutive model has two main sources of error. At first, the model's internal state variables cannot fully reflect the microstructure evolution and its correlation with different features during alloy deformation, and some coefficients and their dependence on temperature and strain rate need more investigation, leading to some deviation in the model's physical interpretation. Secondly, the genetic algorithm for solving model constants requires many iteration parameters, causing local convergence.

Finite element implementation and analysis of the constitutive model

User material subroutine

The constitutive equation of ZK60 magnesium alloy considering damage evolution can be implemented into the commercial finite element solver ABAQUS via the user-defined subroutine VUMAT. The flow chart of the subroutine is shown in Figure 8.

Simulation of uniaxial tensile experiment

Considering the distribution characteristics of loading and deformation in the effective length range of the specimen in the actual tensile process, a characteristic segment is taken in the gauge section for modelling. The characteristic region with a size of 2.5mm–0.25mm–0.25 mm is selected as shown in Figure 9(a); as shown in Figure 9(b), considering the symmetry characteristics of the characteristic segment, symmetric constraints are set about the xz plane and xy plane; a fixed constraint in the x-direction is applied to the left side of the characteristic element, and a positive x-direction displacement is applied to the right side. As shown in Figure 9(c), C3D8 elements are used, and the number of uniformly divided elements is 100. Figure 10(a)–(d) shows the comparison between the prediction values and the experimental values of the unified constitutive model.

Evolution of internal variables

The damage prediction values of the unified constitutive model are presented in Figure 11. Within the temperature range of the experiment, as the deformation temperature and strain rate increase, the nucleation and growth rates as well as the damage factor also increase, which is consistent with the microscopic observation results.

At the initial deformation, the damage factor increases slowly from zero. Then, due to the accumulation of microdefects, the damage factor increases rapidly until the material fails. Because the continuous deformation of ZK60 magnesium alloy at lower temperatures mainly depends on dislocation slip, and the atomic motion is not active, a large number of dislocations accumulate at the grain boundaries, the lattice is severely distorted, and the dislocation movement is difficult to continue, resulting in the rapid expansion of internal cracks and macroscopic fracture.

At higher strain rates, the nucleation rate increases sharply with strain, and then followed by a slow increase (Figure 11a2). In contrast, the growth rate increases significantly at larger strains (Figure 11b2). This indicates that microdefects mainly nucleate at the initial stage of deformation, but as the deformation continues, microdefects mainly grow and connect with each other.

As depicted in Figure 12, the normalized dislocation density exhibits a rapid increase in the early stage of deformation, owing to the multiplication and entanglement of dislocations, which leads to work hardening. Then, it reaches a peak and remains stable, which is attributed to the rearrangement and annihilation of dislocations by static and dynamic recovery. Overall, the normalized dislocation density increases with the increase of deformation temperature and the strain rate.

Conclusions

At the same strain rate, the fracture strain of ZK60 magnesium alloy increases with temperature, while the initial yield stress, peak stress and hardening rate all decrease. At the same temperature, the fracture strain decreases with strain rate, and work hardening is evident. The damage of ZK60 magnesium alloy is caused by the nucleation and growth of the microvoids and cracks. The main factor for damage might be the microvoids and cracks formed around the second-phase particles and/or grain boundaries that coalesce to form macrocracks.

A unified damage constitutive equation was established, which coupled the dislocation density and hardening rate, and modified the nucleation and growth of microvoids and microcracks. The calculated fracture strain agrees well with the experimental data, with the maximum values of AARE and RMSE being 9.63% and 0.0264, respectively, indicating the reliability of the model prediction.

The unified constitutive model was embedded into ABAQUS through the VUMAT subroutine, and the evolution laws of physical internal variables during uniaxial tensile process at different temperatures were predicted. The predicted flow stress agrees well with the experimental data. It is predicted that the normalized dislocation density increased with increasing deformation temperature and strain rate; microdefects mainly nucleated in the early stage of deformation and increased rapidly with deformation.

Footnotes

Acknowledgements

This work was supported by the [National Natural Science Foundation of China] under Grant [number 52175285]; [National Natural Science Foundation of China] under Grant [number 52161145407]; [National Natural Science Foundation of China] under Grant [number U22A20186]; High-performance Fastener Project [number TC220H063].

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This work was supported by the National Natural Science Foundation of China, (grant number 52161145407, 52175285, U22A20186), and High-performance Fastener Project (number TC220H063).

References

Zhou

, et al. Interactions between twins and dislocations in ZK60 alloy under dynamic compression. Mater Sci Technol 2022; 38: 342–351.

Cho

Kim

Han

, et al. Microstructure characterization of ZK60 magnesium alloys using TEM and HR-EBSD. Microsc Microanal 2013; 19: 8–12.

Hadadzadeh

Wells

. Analysis of the hot deformation of ZK60 magnesium alloy. J Magnes Alloy 2017; 5: 369–387.

Liu

Cui

. Characterization of hot deformation behavior of extruded ZK60 magnesium alloy using 3D processing maps. Mater Des 2014; 56: 889–897.

Shang

Zhang

Wang

, et al. The effect of stress state and strain partition mode on the damage behavior of a Mg-Ca alloy. Int J Plast 2021; 144: 103040.

Hadadzadeh

Wells

. Analysis of the hot deformation of ZK60 magnesium alloy. J Magnes Alloy 2017; 5: 369–387.

Qian

Sun

, et al. Investigation into the fracture behavior of ZK60 Mg alloy rolling sheet under different stress triaxiality. J Mater Res Technol 2023; 27: 7368–7379.

Cockcroft

. Ductility and workability of metals. J Metals 1968; 96: 2444.

Johnson

Cook

. Fracture characteristics of three metals subjected to various strains, strain rates, temperatures and pressures. Eng Fract Mech 1985; 21: 31–48.

10.

McClintock

. A criterion for ductile fracture by the growth of holes. J Appl Mech 1968; 35: 363–371.

11.

Oyane

. Criteria of ductile fracture strain. Bull JSME 1972; 15: 1507–1513.

12.

Rice

Tracey

. On the ductile enlargement of voids in triaxial stress fields. J Mech Phys Solids 1969; 17: 201–217.

13.

Gurson

. Continuum theory of ductile rupture by void nucleation and growth: part I-Yield criteria and flow rules for porous ductile media. J Eng Technol 1977; 99: 2–15.

14.

Tvergaard

Needleman

. Analysis of the cup-cone fracture in a round tensile bar. Acta Metall 1984; 32: 157–169.

15.

Hayhurst

Wong

Vakili-Tahami

. The use of CDM analysis techniques in high temperature creep failure of welded structures. JSME Int J Ser A Solid Mech Mater Eng 2002; 45: 90–97.

16.

Lin

Mohamed

Balint

, et al. The development of continuum damage mechanics-based theories for predicting forming limit diagrams for hot stamping applications. Int J Damage Mech 2014; 23: 684–701.

17.

Huo

Lin

Bai

, et al. Prediction of microstructure and ductile damage of a high-speed railway axle steel during cross wedge rolling. J Mater Process Technol 2017; 100: 359–369.

18.

Xiao

Wang

, et al. Constitutive modeling of flow behavior and microstructure evolution of AA7075 in hot tensile deformation. Mater Sci Eng A 2018; 712: 704–713.

19.

Huo

Yang

, et al. A novel unified visco-plastic damage constitutive model considering stress state of TC16 titanium alloy during cold deformation. J Mater Eng Perform 2023; 32: 4522–4540.

20.

Liu

. Formation mechanism and evolution of plastic damage in billet during reduction pretreatment. Metals 2023; 13: 747.

21.

Feng

Huang

Meng

, et al. A constitutive and fracture model for AZ31B magnesium alloy in the tensile state. Mater Sci Eng A 2014; 594: 334–343.

22.

Kim

Lee

, et al. Failure prediction of AZ31B magnesium alloy sheet based on a micro-mechanical void model incorporating the asymmetric plasticity constitutive law. Int J Plast 2017; 94: 98–121.

23.

Zhang

Badreddine

Yue

, et al. Failure prediction of magnesium alloys based on improved CDM model. Int J Solids Struct 2021; 217: 155–177.

24.

Qian

, et al. Investigation into the anisotropic damage behavior of LA103Z Mg-Li alloy rolling sheet using X-ray computed tomography and numerical modeling. Fatigue Fract Eng Mater Struct 2023; 27: 7368–7395.

25.

Kondori

Morgeneyer

Helfen

, et al. Void growth and coalescence in a magnesium alloy studied by synchrotron radiation laminography. Acta Mater 2018; 155: 80–94.

26.

Nie

Shin

Zeng

. Microstructure, deformation, and property of wrought magnesium alloys. Metall Mater Trans A 2020; 51: 6045–6109.

27.

Arul Kumar

Clausen

Capolungo

, et al. Deformation twinning and grain partitioning in a hexagonal close-packed magnesium alloy. Nat Commun 2018; 9: 4761.

28.

Besson

. Continuum models of ductile fracture: a review. Int J Damage Mech 2010; 19: 3–52.

29.

Xiao

Tang

, et al. Flow and fracture study for ZK60 alloy at dynamic strain rates and different loading states. Mater Sci Eng A 2018; 724: 208–219.

30.

Tang

Chen

Zhao

, et al. Experimental study on tensile deformation behaviors under room and elevated temperatures induced by microstructure inhomogeneity of ZK60 Mg with a longitudinal weld. J Magnes Alloy 2023; 11: 4338–4356.

31.

Jordon

Gibson

Horstemeyer

, et al. Effect of twinning, slip, and inclusions on the fatigue anisotropy of extrusion-textured AZ61 magnesium alloy. Mater Sci Eng A 2011; 528: 6860–6871.

32.

Lemaitre

. A continuous damage mechanics model for ductile fracture. J Eng Mater Technol 1985; 107: 83–89.

33.

Lin

Liu

Dean

. A review on damage mechanisms, models and calibration methods under various deformation conditions. Int J Damage Mech 2005; 14: 299–319.

34.

Yang

Wang

Liu

, et al. Behavior and modeling of flow softening and ductile damage evolution in hot forming of TA15 alloy sheets. Mater Des 2015; 85: 135–148.

35.

Lin

Liu

Farrugia

DCJ

, et al. Development of dislocation-based unified material model for simulating microstructure evolution in multipass hot rolling. Phil Mag 2005; 85: 1967–1987.

36.

Lin

Yang

. GA-based multiple objective optimisation for determining viscoplastic constitutive equations for superplastic alloys. Int J Plast 1999; 15: 1181–1196.