A Generalized Non-Linear Flow Law Based on Modified Zerilli-Armstrong Model and Its Implementation into Abaqus/Explicit FEM Code

Pierre Tize Mha, Amèvi Tongne, Olivier Pantalé^*

Production Engineering Laboratory, INPT/ENIT, Toulouse University, Tarbes, France.

DOI: 10.4236/wjet.2022.102021                

Abstract

Non-linear numerical modeling, widely used in research and development to understand many complex processes such as forming or machining, does not guarantee the success of a study to be performed. Indeed, the numerical simulation uses finite element codes where the models already integrated are not based on shapes adjustable to any type of study. In this study, a new form of non-linear constitutive flow law based on the Modified Zerilli-Armstrong model, which can answer the above problem, has been developed to apply it to the numerical simulation of two different tests (a quasi-static compression test, the necking of a circular bar). This flow law is based on the modified Zerilli-Armstrong model, which, together with the new modified Johnson-Cook model, has been compared to appreciate the relevance of the proposal. For that, an implementation of this new law via the VUHARD subroutine into the Abaqus/Explicit finite element code was made to model the two tests. The comparison of the results obtained (from identification) by our proposed law with those obtained using the NMJC shows that this new law better approaches the experiments than the other one. This is also shown through the numerical results using the Abaqus software. It can be said that this way of formulating a flow law allows highlighting the great performance of the proposed approach. Although this law has been only used for quasi-static tests, we can say that it can also be used in dynamic tests.

Keywords

Zerilli-Armstrong Flow Law, Constitutive Behavior, Finite Element Method, Numerical Implementation, Johnson-Cook Flow Law

1. Introduction

Due to technological advances, industries are constantly looking for innovations to maintain their competitiveness. For the metallurgical industry, these innovations involve the development of new materials. In conjunction with the implementation of these new materials, it becomes necessary to characterize their thermomechanical behavior, in particular in order to feed the databases used by the numerical simulation models more and more widespread nowadays [1] [2].

Among the different approaches used, the finite element method (FEM) is one of the numerical methods that have become the most popular nowadays. Like any numerical method, the accuracy of the solution depends very strongly on the constitutive equation describing the behavior of the material and implemented in the computational code [3] [4]. This constitutive equation is linked to the mathematical model that describes it. This choice is not entirely free, insofar as the choice of the behavior law (speaking of an analytically formulated model) used to describe the behavior of a material is guided by its availability in the finite element codes [5]. This forces an adaptation between the studied material (which is real) and the behavior law (representative mathematical model) available in the FEM code which is not necessarily the one able to describe the behavior of the studied material. The problem of the non-flexibility of the models available in the finite element codes appears immediately and limits the adequacy of the model to the reality of the behavior since the behavior laws established and available in the FEM codes are for the most part specific to a particular material or class of materials. In other words, as soon as they are applied to materials different from the one for which they were originally developed, they find their limits because they do not manage to represent the reality described by the experiment as one would wish. It is therefore interesting to think about developing a method that can be applied to any material of any class and implemented in a finite element code. This approach involves a number of parameters depending on the material studied, but with a more global mathematical form adaptable to a wider range of materials. In this study, we propose an approach allowing better adaptability of the behavior law to various materials. To justify the relevance of this proposal, let us first examine the models developed so far and try to see if some of them meet this need.

According to the literature, analytical models describing the behavior of materials can be grouped into three categories: empirical models, semi-empirical models, and physical models. Empirical models describe the behavior of the material by considering macroscopic quantities without considering its microstructure, whereas physical models consider the microstructure. As for the semi-empirical models, they try to be at the border between the empirical and physical models. Among a large number of models, the Johnson-Cook (JC) [6] [7] and the Zerilli-Armstrong (ZA) [8] models are the best known and available in finite element codes. The JC model is widely used because it is simple to use and has fewer constants to be determined [9] [10]. However, this model, in its original state, contains some problems related to the phenomena describing the behavior of the material. We refer here to strain hardening, strain rate hardening, and thermal effects that are described separately in the JC model as analyzed by Jia et al. [11]. To overcome this difficulty, several modified forms of the Johnson-Cook model have been developed in the past [12] - [18].

Due to its more physical formulation, the ZA model is preferred over the JC model because it allows for the coupling between strain rate and temperature effects [19] [20] [21]. However, the drawbacks of this model are related to the temperature limit value of the test. Indeed, this model applies only to low strain rates and for temperatures limited to $0.6T_m$ (where $T_m$ is the melting temperature of the material) [22] [23] [24] [25]. Moreover, this model, in its original form, does not take into account the coupling between deformation and temperature effects [26]. In order to reduce these shortcomings (applying this model over a wider temperature range) of the ZA model, Modified Zerilli-Armstrong forms (MZA) have been proposed by several authors over the last past years [27] [28] [29] [30] [31].

These two approaches are among the most popular and the choice depends on the effectiveness of the method for a given problem. So we have, on the one hand, the JC formulation and its variants and on the other hand the ZA approach with its variants too. Another approach proposed by Lin et al. [14] is to combine the JC model with the MZA model (in its modified form proposed by Samantaray et al. [26] ). This approach allows combining on the one hand the factor taking into account the strain hardening effects described by the JC model and on the other hand the factor describing the strain rate effects and the thermal softening described by the MZA model in a unified formulation. The major difficulty of all the variants mentioned above is the restriction of the parameters on which the model depends, thus limiting the applicability domain of the developed model. Therefore, it seems necessary to look for a method that does not limit the number of parameters, unlike the approaches mentioned above, and that can describe more accurately the reality observed by the experiment.

Artificial neural networks (ANNs) offer an alternative way to solve these problems. ANNs are currently used to solve problems that are difficult to conceptualize using traditional computational methods. Thus, unlike a classical approach based on a regression method, an artificial neural network does not need to know the mathematical form of the model it seeks to reproduce and can therefore learn directly from experimental data without any prior assumption on the nature of the input parameters and their interactions. Thus, artificial neural networks can enable new approaches for modeling material behavior and have been successfully applied for the prediction of constitutive relationships of some metals and alloys [32] [33] [34] [35]. Although this approach is very effective, it does not allow the interpretation of the obtained results, because it does not take into account the physics of the phenomena for the considered material.

From the above, we can retain that there are two main families of models to describe the behavior of a material. On the one hand, analytical models present a problem of flexibility, and on the other hand, implicit models, based on neural networks, are more adaptable but present a problem of physical interpretation. It is therefore interesting to develop a physical model that is able to adapt to any material and can be efficiently implemented in a finite element code. In this perspective, we propose here a new formulation based on the MZA model, whose interest is to answer the two problems mentioned above. To appreciate the relevance of our approach, we used the Abaqus/Explicit software which offers the possibility to implement a new behavior law using a FORTRAN VUMAT or VUHARD subroutine following the example of the works proposed by [36] [37] [38] [39] [40].

The main paragraphs of this work are as follows: Section 2 presents two known constitutive laws which will serve as the basis of the proposed one presented in Section 3. Since this model is new and will be used to perform a numerical model using Abaqus, the three derivatives of flow stress ${\sigma }^y$ with respect to the strain ${\epsilon }^p$ , the strain rate ${\stackrel{\dot{}}{\epsilon }}^p$ and temperature are required. So, determination of these derivatives is done in this Section. Section 4 is devoted to the comparison of the proposed model to the other one where identification is done on a 100Cr6 steel. Section 5 details the numerical results of the models identified in Section 4 using the cylindrical compression test and the necking of circular benchmark tests. Conclusion and future works are detailed in Section 6.

2. Some Constitutive Laws

The complete characterization of the thermomechanical behavior of a material is very important before its use in numerical simulations. This characterization usually involves experimental tests that serve as a basis for the development of constitutive laws. As mentioned in the introduction, many constitutive laws exist today, and selecting one of them based on explicit criteria hinders the development of computational models for new materials. In this section, we present a few of the best-known behavior laws. These are the Zerilli-Armstrong (ZA) model and one of its derived forms, the Modified Zerilli-Armstrong (MZA). The MZA model is the one that was used as the basis for the development of the PTM model presented in Section 3.

The Zerilli-Armstrong model

As mentioned earlier, Johnson-Cook’s law, widely used in the scientific community and implemented as a base in the Abaqus FEM code, is an empirical model for which the effects of strain, strain rate, and temperature on material behavior are considered separately. Since a coupling exists in reality between these three variables, it is preferable to opt for a Zerilli-Armstrong (ZA) model. In fact, the Zerilli-Armstrong model is based on a physical consideration of heat-activated dislocation motion and the interaction of individual dislocations with each other as shown by Hull et al. [41]. The ZA model is therefore a constitutive law based on a study of the material microstructure. For a FCC material, the original equation describing the Zerilli-Armstrong model [8] is given by:

${\sigma }^y=A_0+A_1{\epsilon }^{p^n}\exp \left[-A_{2}T+A_{3}T\ln \left(\frac{{\stackrel{\dot{}}{\epsilon }}^p}{{\stackrel{\dot{}}{\epsilon }}_0}\right)\right]$ (1)

where ${\sigma }^y$ is the flow stress of the material, ${\epsilon }^p$ is the equivalent plastic strain, ${\stackrel{\dot{}}{\epsilon }}^p$ is the equivalent plastic strain rate, ${\stackrel{\dot{}}{\epsilon }}_0$ is the reference strain rate and T is the absolute temperature within the material. The coefficients $A_0,A_1,A_2,A_3$ and n are the parameters of the model to be identified for a given material and ${\stackrel{\dot{}}{\epsilon }}_0$ is the reference strain rate as proposed by many models in the literature. It is important to note that this model is globally split in two terms. The first term is a constant stress threshold representing the yield strength of the material and the second term is a thermo-viscous component considering the variation of the work hardening threshold of the material with temperature and strain rate. With face-centered cubic (FCC) materials, the strain hardening shows a greater dependence on both temperature and strain rate. The thermal component of strain hardening in FCC materials is physically interpreted as the resistance to dislocation motion by the mutual interaction of dislocations [42]. The flow stress threshold of an FCC material is independent of temperature and strain rate, so that the flow stress threshold is constant and equal to the initial stress [27] [42]. The limitation of the ZA model is that it cannot be used for flow stress prediction at high strain rates and temperatures above $0.6T_m$ . To circumvent this problem, Samantaray et al. [26] proposed in 2009 a modified form of the Zerilli-Armstrong law (MZA), which will be presented in the next subsection.

The Modified Zerilli-Armstrong model

Formulation of the Modified Zerilli-Armstrong model

As mentioned at the end of Subsection 2, the Modified Zerilli-Armstrong model (MZA) model is an extension of the Zerilli-Armstrong model. Since the new approach proposed in this article is based on the MZA model, it is necessary to recall, on the one hand, the motivations for the formulation of the MZA model and, on the other hand, the methodology for identifying its parameters. Equation that describes the MZA model, as proposed by Samantaray et al. [26], is given by:

${\sigma }^y=\left(C_1+C_2{\epsilon }^{p^n}\right)\exp \left[-\left(C_3+C_4{\epsilon }^p\right)\left(T-T_0\right)+\left(C_5+C_6\left(T-T_0\right)\right)\ln \left(\frac{{\stackrel{\dot{}}{\epsilon }}^p}{{\stackrel{\dot{}}{\epsilon }}_0}\right)\right]$ (2)

where the 7 coefficients $C_i$ and n are the parameters of the model to be identified for a given material. This formulation was initially developed from a compression test of 15Cr-15Ni-2.2Mo Ti-modified austenitic stainless steel (the D9 alloy). This test was performed over a range of plastic strains ( ${\epsilon }^p=0.1\text{-}0.5$ ), strain rates ( ${\stackrel{\dot{}}{\epsilon }}^p={10}^{-3}\text{-}\mathrm{1\; <mtext>\; \hspace{0.17em}\; </mtext>}{\text{s}}^{-1}$ ) and temperatures ( $T=800\text{-}{1200}^{\circ }\text{C}$ ). In attempting to apply the original Johnson-Cook model, the authors of this work found that it was unable to describe the behavior of this type of material over the ranges of strain, strain rate, and temperature mentioned above. The same conclusion was made regarding the ZA model which was also unable to correctly describe the behavior of this material. Based on the ZA model, Samantaray et al. [26] have therefore proposed a modified form of the ZA model to better describe the behavior of this material. This modified model takes into account the effect of isotropic strain hardening, strain rate hardening, thermal softening and the coupled effects of strain, strain rate and temperature on the flow stress. It predicts the evolution of the high-temperature flow threshold of the D9 alloy over the above ranges of strain, strain rate and temperature with good correlation to the experimental data and good generalization. The correlation coefficient obtained is 0.995 while the average error for the whole data set is only 5.3%. The proposed model can predict the flow threshold of the D9 alloy in the hot working range. In the first phase, we will examine the different steps in the development of the model proposed by Samantaray and their justifications.

Equation (1) of the original Zerilli-Armstrong model defines the flow stress in two terms: thermal stress and athermal stress that depends on the grain size of the material and is independent of temperature. However, several studies have shown that this independence decreases when the temperature increases beyond a certain value. Based on these observations, Samantaray questioned the athermal consideration for the D9 material. This approach then allows for the consideration of temperature effects on the yield strength of D9 steel and even other materials.

According to the ZA model (1), the term $A_1{\epsilon }^{p^n}$ is calculated from the y-intercept of the line $\ln \left({\sigma }^y\left(T\right)-A_0\right)$ giving the flow stress at 0˚K. These y-intercepts are calculated for various values of strain to obtain separately the values of $A_1$ and n. Since experimental data are not always available at 0˚K, Samantaray then introduced a reference temperature to incorporate the effect of thermal softening. In view of the above, the term $A_1{\epsilon }^{p^n}$ has been replaced by $C_1+C_2{\epsilon }^{p^n}$ in the MZA model (2).

According to the ZA model (1), in the absence of a strain rate effect, $A_2$ is calculated from the slope of the function ${\sigma }^y\left(T\right)-A_0$ which is independent of the plastic strain in the original formulation. Samantaray having shown on a graphical representation a linear dependence to the plastic strain, he replaced the term $A_2$ by $\left(C_3+C_4{\epsilon }^p\right)$ in his model (2).

Finally, the term $A_{3}T$ is determined from the slope of $\ln {\sigma }^y$ as a function of $\ln \left({\stackrel{\dot{}}{\epsilon }}^p/{\stackrel{\dot{}}{\epsilon }}_0\right)$ , and then the constant $A_3$ is calculated from the slope of the line $A_{3}T$ . In its original formulation, the ZA model predicts the same flow stress for all values of the strain rate. In reality, the slope of $\ln {\sigma }^y$ as a function of $\ln \left({\stackrel{\dot{}}{\epsilon }}^p/{\stackrel{\dot{}}{\epsilon }}_0\right)$ , denoted $S_2$ , shows a linear dependence with ( $T-T_0$ ). This reflects the effect of the strain rate on the flow stress at the reference temperature. Taking into account these remarks, Samantaray introduces the term $C_5+C_6\left(T-T_0\right)$ in his model. For more details, see [26]. With the main steps of the model development assimilated, the identification of the parameters is the subject of the next subsection.

Determination of the Modified Zerrilli-Armstrong model’s constants

This section is devoted to the calculation of the parameters of the MZA model formulated from the set of remarks described previously. Considering the expression of the model given by Equation (2), parameter $C_1$ is deduced direcly from ${\sigma }^y\left({\epsilon }^p=0,T=T_0,{\stackrel{\dot{}}{\epsilon }}^p={\stackrel{\dot{}}{\epsilon }}_0\right)$ . When the current strain rate ${\stackrel{\dot{}}{\epsilon }}^p$ is identical to the reference strain rate ${\stackrel{\dot{}}{\epsilon }}_0$ , Equation (2) becomes:

${\sigma }^y=\left(C_1+C_2{\epsilon }^{p^n}\right)\exp \left[-\left(C_3+C_4{\epsilon }^p\right)\left(T-T_0\right)\right]$ (3)

Taking the logarithm of the previous relationship we have:

$\ln {\sigma }^y=\ln \left(C_1+C_2{\epsilon }^{p^n}\right)-\left(C_3+C_4{\epsilon }^p\right)\left(T-T_0\right)$ (4)

By plotting the function $\ln {\sigma }^y\left(T-T_0\right)$ the y-intercept and the slope of this line can be deduced as $I_1=\ln \left(C_1+C_2{\epsilon }^{p^n}\right)$ and $S_1=-\left(C_3+C_4{\epsilon }^p\right)$ . From the expression of $I_1$ we can write:

$\ln \left(\exp \left(I_1\right)-C_1\right)=\ln C_2+n\ln {\epsilon }^p$ (5)

By plotting the function $\ln \left(\exp \left(I_1\right)-C_1\right)$ versus $\ln {\epsilon }^p$ we can deduce n and $C_2$ . The plot of $S_1$ versus ${\epsilon }^p$ allows deducing $C_3$ and $C_4$ from the y-intercept and the slope, respectively. Once the first four terms of the model are calculated, we can now take the logarithm of the entire MZA model given by Equation (2).

$\ln {\sigma }^y=\ln \left(C_1+C_2{\epsilon }^{p^n}\right)-\left(C_3+C_4{\epsilon }^p\right)\left(T-T_0\right)+S_2\ln \left(\frac{{\stackrel{\dot{}}{\epsilon }}^p}{{\stackrel{\dot{}}{\epsilon }}_0}\right)$ (6)

with $S_2=C_5+C_6\left(T-T_0\right)$ . Finally, the graph $\ln {\sigma }^y$ versus $\ln \left({\stackrel{\dot{}}{\epsilon }}^p/{\stackrel{\dot{}}{\epsilon }}_0\right)$ allows deducing $C_5$ and $C_6$ from the slope $S_2$ .

3. Proposed Model (PTM Model)

The model proposed in this study is an extension based on the formulation of the MZA model proposed by Samantaray. In the MZA model, the functions $I_1\mathrm{,}{S}_1$ and $S_2$ are defined as linear functions. However, there are materials where these functions cannot be approximated by linear functions as illustrated in Figure 2(c) showing the evolution of $S_1$ as a function of plastic strain ${\epsilon }^p$ . Therefore, the MZA model cannot reproduce the experimental results as expected [43] - [49]. It should also be remembered that this proposal is not only to extend the MZA model for better results, but also to make it suitable for both low and high strain rate tests.

Formulation of the PTM model

The PTM model proposed in the present study is an extension of the MZA model for which the linear terms $I_1\mathrm{,}{S}_1$ and $S_2$ described by Equation (2) are now polynomials with an unknown degree depending on the experimental data. In this way, our approach allows us to make the model flexible and, therefore, applicable to several different types of materials. Moreover, in order to further improve the generalization of the proposed model, the coefficients of the polynomial representing $S_2$ are still other polynomials depending on the deformation. In the work of Samantaray et al. [26], these coefficients were considered as constants. The equation describing the PTM model is therefore:

$\begin{array}{c}{\sigma }^y\left({\epsilon }^p,{\stackrel{\dot{}}{\epsilon }}^p,T\right)=\left({\displaystyle \underset{i=0}{\overset{q}{\sum }}}{A}_i{\epsilon }^{p^i}\right)\exp [\left({\displaystyle \underset{j=0}{\overset{r}{\sum }}}{B}_j{\epsilon }^{p^j}\right)\left(T-T_0\right)\\ +\left({\displaystyle \underset{k=0}{\overset{s}{\sum }}}\left({\displaystyle \underset{l=0}{\overset{t}{\sum }}}{C}_k^l{\epsilon }^{p^l}\right){\left(T-T_0\right)}^k\right)\ln \left(\frac{{\stackrel{\dot{}}{\epsilon }}^p}{{\stackrel{\dot{}}{\epsilon }}_0}\right)]\end{array}$ (7)

where $A_i\mathrm{,}{B}_j\mathrm{,}{C}_k^l$ are the parameters of the model to be identified using the procedure proposed in the following subsection. Quantities q, r, s and t define the degree of the polynomials used to describe the behavior of the material. The larger these quantities are, the more parameters need to be identified for the PTM model.

Determination of the PTM model’s parameters

In accordance with what has been presented in the previous sections and in view of the modifications made by the formulation of the PTM model, the determination of the parameters of this model is calculated thanks to the LMFIT Python library [50] and according to the following three steps:

Step 1: Determination of the two functions $I_1\left({\epsilon }^p\right)$ and $S_1(\; \epsilon \; p\; )$

At the reference strain rate, when ${\stackrel{\dot{}}{\epsilon }}^p={\stackrel{\dot{}}{\epsilon }}_0$ , Equation (7) is written as:

${\sigma }^y\left({\epsilon }^p,T\right)=\left({\displaystyle \underset{i=0}{\overset{q}{\sum }}}{A}_i{\epsilon }^{p^i}\right)\exp \left[\left({\displaystyle \underset{j=0}{\overset{r}{\sum }}}{B}_j{\epsilon }^{p^j}\right)\left(T-T_0\right)\right]$ (8)

Taking the logarithm of each member of this equation we obtain:

$\ln {\sigma }^y\left({\epsilon }^p,T\right)=I_1\left({\epsilon }^p\right)+S_1\left({\epsilon }^p\right)\left(T-T_0\right)$ (9)

where:

$I_1\left({\epsilon }^p\right)=\ln \left({\displaystyle \underset{i=0}{\overset{q}{\sum }}}{A}_i{\epsilon }^{p^i}\right)\text{and}{S}_1\left({\epsilon }^p\right)={\displaystyle \underset{j=0}{\overset{r}{\sum }}}{B}_j{\epsilon }^{p^j}$ (10)

Plotting the function $\ln {\sigma }^y\left(T\right)$ for a given value of ${\epsilon }^p$ allows us to deduce $I_1\left({\epsilon }^p\right)$ and $S_1\left({\epsilon }^p\right)$ from the y-intercept and the slope, respectively. These functions are calculated for all plastic strain values ${\epsilon }^p$ available in the experimental database. Analysis of those results allows us to propose some adequate values for the q and r quantities. Once those values have been selected, an identification procedure allows us to compute the q parameters $A_i$ and r parameters $B_j$ of the PTM model.

• Step 2: Determination of the function $S_2\left({\epsilon }^p\mathrm{,}T\right)$ as polynomials depending on the temperature T

Once Step 1 is done, it is possible to take the logarithm of the whole Equation (7) to obtain:

$\ln {\sigma }^y\left({\epsilon }^p\mathrm{,}{\stackrel{\dot{}}{\epsilon }}^p\mathrm{,}T\right)-I_1\left({\epsilon }^p\right)-S_1\left({\epsilon }^p\right)\left(T-T_0\right)=S_2\left({\epsilon }^p\mathrm{,}T\right)\ln \left(\frac{{\stackrel{\dot{}}{\epsilon }}^p}{{\stackrel{\dot{}}{\epsilon }}_0}\right)$ (11)

where:

$S_2\left({\epsilon }^p,T\right)={\displaystyle \underset{k=0}{\overset{s}{\sum }}}\left({\displaystyle \underset{l=0}{\overset{t}{\sum }}}{C}_k^l{\epsilon }^{p^l}\right){\left(T-T_0\right)}^k$ (12)

From this latter the plot of $\ln {\sigma }^y-I_1\left({\epsilon }^p\right)-S_1\left({\epsilon }^p\right)\left(T-T_0\right)$ versus $\ln \left({\stackrel{\dot{}}{\epsilon }}^p/{\stackrel{\dot{}}{\epsilon }}_0\right)$ allows us to deduce $S_2\left({\epsilon }^p\mathrm{,}T\right)$ as its slope for a given plastic strain value ${\epsilon }^p$ . The evolution of $S_2$ versus $\left(T-T_0\right)$ allows to find the same function for all plastic strains which is described in the polynomial forms of degree k in Equation (7). The work does not stop at this stage; it is necessary to evaluate the dependence of each coefficient of $S_2\left({\epsilon }^p\mathrm{,}T\right)$ to the plastic strain ${\epsilon }^p$ . This has never been considered in the previous works [43] - [49] [51].

• Step 3: Determination of the coefficients of the function $S_2\left({\epsilon }^p\mathrm{,}T\right)$ as polynomials depending on the plastic strain ${\epsilon }^p$

This step is a direct consequence of the previous one, which consists in determining each coefficient of the function $S_2\left({\epsilon }^p\mathrm{,}T\right)$ as polynomial functions. Once we know the order of the polynomial that describes the function $S_2\left({\epsilon }^p\mathrm{,}T\right)$ as described in the previous step, for each strain value the corresponding coefficients of $S_2\left({\epsilon }^p\mathrm{,}T\right)$ should be recorded in the empty lists initialized before. Thus, with the empty lists initialized, for a given strain value, it is necessary to go through all the temperatures and deduce the values of each coefficient of the function $S_2\left({\epsilon }^p\mathrm{,}T\right)$ for this strain which will be recorded in the empty lists initialized at the beginning. Thus, an algorithm must be written which goes through all the values of the strains and adds to each list of coefficients of the polynomial $S_2\left({\epsilon }^p\mathrm{,}T\right)$ the value deduced from the corresponding strain. It is therefore necessary to evaluate the dependence of each parameter (described by one list) on the strain in order to finally deduce its polynomial from which we also deduce its coefficients. This is how the operation is performed for all the lists. Finally, we calculate all the parameters of the proposed model. This part will be detailed more deeply in Section 4.

Implementation of PTM model into Abaqus

The usual method when implementing a new behavior law in the Abaqus finite element code is to program a VUMAT or VUHARD in FORTRAN subroutine. To introduce a new behavior law using the VUHARD (in our case), Abaqus only required the definition of the function that describes the behavior law ${\sigma }^y\left({\epsilon }^p\mathrm{,}{\stackrel{\dot{}}{\epsilon }}^p\mathrm{,}T\right)$ and its derivatives with respect to the plastic strain ${\epsilon }^p$ , the plastic strain rate ${\stackrel{\dot{}}{\epsilon }}^p$ , and the temperature T. Thus, in accordance with the proposed approach, the three derivatives of the PTM model are given by Equation (13) obtained directly by analytical derivation from Equation (7):

$(\begin{array}{l}\frac{\partial {\sigma }^y}{\partial {\epsilon }^p}=\{\left({\displaystyle \underset{i=1}{\overset{q}{\sum }}}i{A}_i{\epsilon }^{p^{i-1}}\right)+\left({\displaystyle \underset{i=0}{\overset{q}{\sum }}}{A}_i{\epsilon }^{p^i}\right)[\left({\displaystyle \underset{j=1}{\overset{r}{\sum }}}j{B}_j{\epsilon }^{p^{j-1}}\right)\left(T-T_0\right)\\ +\left({\displaystyle \underset{k=0}{\overset{s}{\sum }}}\left({\displaystyle \underset{l=1}{\overset{t}{\sum }}}l{C}_k^l{\epsilon }^{p^{l-1}}\right){\left(T-T_0\right)}^k\right)\ln \left(\frac{{\stackrel{\dot{}}{\epsilon }}^p}{{\stackrel{\dot{}}{\epsilon }}_0}\right)]\}\exp \left(C_{mf}\right)\\ \frac{\partial {\sigma }^y}{\partial {\stackrel{\dot{}}{\epsilon }}^p}=\left({\displaystyle \underset{i=0}{\overset{q}{\sum }}}{A}_i{\epsilon }^{p^i}\right)\left[{\displaystyle \underset{k=0}{\overset{s}{\sum }}}\left({\displaystyle \underset{l=0}{\overset{t}{\sum }}}{C}_k^l{\epsilon }^{p^l}\right){\left(T-T_0\right)}^k\right]\frac{\exp \left(C_{mf}\right)}{{\stackrel{\dot{}}{\epsilon }}^p}\\ \frac{\partial {\sigma }^y}{\partial T}=\left({\displaystyle \underset{i=0}{\overset{q}{\sum }}}{A}_i{\epsilon }^{p^i}\right)\left\{\left({\displaystyle \underset{j=0}{\overset{r}{\sum }}}{B}_j{\epsilon }^{p^j}\right)+\left[{\displaystyle \underset{k=1}{\overset{s}{\sum }}}k\left({\displaystyle \underset{l=0}{\overset{t}{\sum }}}{C}_k^l{\epsilon }^{p^l}\right){\left(T-T_0\right)}^{k-1}\right]\ln \left(\frac{{\stackrel{\dot{}}{\epsilon }}^p}{{\stackrel{\dot{}}{\epsilon }}_0}\right)\right\}\exp \left(C_{mf}\right)\end{array}$ (13)

where the $C_{mf}$ term is given by:

$C_{mf}=\left({\displaystyle \underset{j=0}{\overset{r}{\sum }}}{B}_j{\epsilon }^{p^j}\right)\left(T-T_0\right)+\left({\displaystyle \underset{k=0}{\overset{s}{\sum }}}\left({\displaystyle \underset{l=0}{\overset{t}{\sum }}}{C}_k^l{\epsilon }^{p^l}\right){\left(T-T_0\right)}^k\right)\ln \left(\frac{{\stackrel{\dot{}}{\epsilon }}^p}{{\stackrel{\dot{}}{\epsilon }}_0}\right)$ (14)

4. Application

After presenting the models that will be used in this study, setting up the formulation of the proposed model and the theory of identifying its parameters, we will now show in detail identifying the parameters of the PTM model. In the numerical simulations and validation tests, only one material is used, and the identification results are compared each time.

Identification of the PTM parameters for a 100Cr6 material

Data collection from literature

As mentioned above, the compression test is performed on 100Cr6 steel for a temperature range of 750˚C - 1300˚C and a strain rate range of 0.001 - 0.1 s⁻¹. For more details, see Zhou et al. [17]. Figure 1 shows the plot of flow stresses versus plastic strain values for different strain rates and temperatures as extracted from the above publication.

Since we did not have the opportunity to contact the author of this paper to obtain the raw data, we used WebPlotDigitizer-4.4-linux-x64 to extract these data from the graphs of the publication, which explains, for example, the non-conformity of the stress values at ${\epsilon }^p=0$ with those of Zhou et al. [17]. It should also be noted that in Zhou’s publication, the values at this level of plastic strain are not considered in the parameter identification because they are unreadable. Only 14 strain values were used for the identification of the model

Figure 1. True stress-strain curves for 100Cr6 steel at different strain rates [17].

parameters. As presented by Zhou et al. [17], the results will be displayed in two temperature ranges, namely 750˚C - 850˚C on the one hand and 900˚C - 1300˚C on the other. This implies that the parameters of our model depend on these two temperature ranges. For both temperature ranges, the reference strain rate is the same while the reference temperatures are different: 750˚C as reference temperature for the first range and 900˚C for the second.

Identification of the parameters of the PTM model

To have a clearer overview of the method developed in this study for the determination of the PTM model parameters as presented in Section 3, we detail here after, only for the second temperature range (900˚C - 1300˚C), how to identify some parameters of the flow law. The other parameters are obtained using the same method, but for the sake of brevity, only the generic part is presented. For the identification procedure, in this temperature range, we have extracted 378 data points ${\sigma }^y\left({\epsilon }^p\mathrm{,}{\stackrel{\dot{}}{\epsilon }}^p\mathrm{,}T\right)$ with 14 different plastic strains ${\epsilon }^p$ , 3 different plastic strain rates ${\stackrel{\dot{}}{\epsilon }}^p$ and 9 different temperatures T.

Referring to the explanation of the method, the first step concerns the determination of the functions $I_1\left({\epsilon }^p\right)$ and $S_1\left({\epsilon }^p\right)$ that appears in Equation (9) when ${\stackrel{\dot{}}{\epsilon }}^p={\stackrel{\dot{}}{\epsilon }}_0=0.1$ . We first plot 14 curves $\ln {\sigma }^y\left(T\right)$ (1 for each plastic strain ${\epsilon }^p$ ) and extract the y-intercept and the slope of those 14 curves as $I_1\left({\epsilon }^p\right)$ and $S_1\left({\epsilon }^p\right)$ respectively. This process is illustrated in Figure 2(a) where 3 $\ln {\sigma }^y\left(T\right)$ plots for different values of the plastic strain ${\epsilon }^p$ are reported. Thus, from the left part of Equation (10), and after plotting $\exp \left(I_1\left({\epsilon }^p\right)\right)$ as illustrated in Figure 2(b) and choosing an appropriate value for the parameter q (here $q=4$ ) defining the order of the polynomial function used to approximate

Figure 2. Determination of the PTM model’s parameters $A_i$ and $B_j$ for a range of temperature (900˚C - 1300˚C).

$I_1\left({\epsilon }^p\right)$ , the five parameters $A_i$ , reported in the right column of Table 1, can be identified. Parameters $B_i$ of function $S_1\left({\epsilon }^p\right)$ are computed using the same approach as illustrated in Figure 2(c) after selecting an adequate order r of the polynomial function (here again $r=4$ ). The obtained values of parameters $B_j$ are also reported in the right column of Table 1.

Once the parameters of those two functions have been identified, the final step then concerns the identification of the parameters of the $S_2\left({\epsilon }^p\mathrm{,}T\right)$ function. The identification of that function is divided into three steps:

1) From Equation (11), we first plot $\ln {\sigma }^y-I_1\left({\epsilon }^p\right)-S_1\left({\epsilon }^p\right)\left(T-T_0\right)$ versus $\ln \left({\stackrel{\dot{}}{\epsilon }}^p/{\stackrel{\dot{}}{\epsilon }}_0\right)$ for all plastic strains ${\epsilon }^p$ and temperatures T and gather the $126=14\times 9$ slopes defining the set of $S_2$ parameters depending on ${\epsilon }^p$ and T respectively.

2) Those 126 data points are now grouped on their plastic strain values and are used to identify 14 $S_2\left(T\right)$ curves. Figure 3(a) illustrates this process for 3 values of the plastic strain. From this later, and after selecting a correct value $s=4$ of the polynomial order we fit each curve $S_2\left(T\right)$ and obtain $5\times 14=70$ new data points.

3) From those later, we now plot them the plastic strain ${\epsilon }^p$ for each value k of the parameter s, as reported in Figure 3(b) for $k=0$ . Then, after choosing the adequate value of $t=4$ , we finally obtain the 25 $C_k^l$ parameters reported in the right column of Table 1.

The exact same procedure has been done to obtain the same parameters of the PTM model, for the temperature range (750˚C - 850˚C), reported in the left column of Table 1.

Comparison and accuracy of the PTM model

To evaluate the relevance of the model proposed in this work, we need to compare its results to those provided by the existing well-known models. This comparison is based on the MZA and NMJC (New Modified Johnson-Cook) models presented by Zhou et al. [17]. Figure 4 (first temperature range) and Figure 5 (second temperature range) show the comparison of the MZA and NMJC models with the PTM model. The comparison is performed using the parameters of the MZA and NMJC models reported in Table 2 and Table 3

Figure 3. Determination of the PTM model’s parameters $C_k^l$ for a range of temperature (900˚C - 1300˚C).

Table 1. Parameters of the proposed model.

Table 2. Parameters of modified Zerilli-Armstrong model [17].

Figure 4. Comparison between experimental data and predicted data of MZA, NMJC and PTM models for a range of temperature (750˚C - 850˚C).

Table 3. Parameters of new modified Johnson-Cook model [17].

Figure 5. Comparison between experimental data and predicted data of MZA, NMJC and PTM models for a range of temperature (900˚C - 1300˚C).

respectively. Figure 4 and Figure 5 show the performance of each model in approximating the experimental data. To evaluate the fidelity of each law with respect to the experimental data, we compute the correlation coefficient (R) and the average relative error (E_AAR) whose relations are given by the following Equations (15) and (16):

$R\left(\%\right)=\frac{{\displaystyle \underset{i=1}{\overset{N}{\sum }}}\left(E_i-\bar{E}\right)\left(P_i-\bar{P}\right)}{\sqrt{{\displaystyle \underset{i=1}{\overset{N}{\sum }}}{\left(E_i-\bar{E}\right)}^2.{\displaystyle \underset{i=1}{\overset{N}{\sum }}}{\left(P_i-\bar{P}\right)}^2}}\times 100$ (15)

$E_{\text{AAR}}\left(\mathrm{<mi>\; \%\; </mi>}\right)=\frac{1}{N}{\displaystyle \underset{i=1}{\overset{N}{\sum }}}\left|\frac{E_i-P_i}{E_i}\right|\times 100$ (16)

where $E_i$ is the experimental value and $P_i$ is the predicted value of the constitutive law, $\bar{E}$ and $\bar{P}$ are the mean values of the experimental and predicted values, respectively. The values of R and E_AAR for the two temperature ranges are reported in Table 4 and Table 5. These tables clearly show that the MZA model is not able to approximate the experimental data well for both temperature ranges simultaneously, while the NMJC model and the PTM model approximate these experimental data much better. This is also clearly visible in the correlation curve shown for the PTM model. Since the NMJC model and the PTM model better represent the experimental data, we will perform a more refined analysis of their results.

For the first temperature range (750˚C - 850˚C) the E_AAR and R values are 3.25% and 99.34% for the NMJC, and 4.07% and 99.18% for PTM model. We observe that the two models differ slightly. However, for the second temperature range,

Table 4. Comparison of E_AAR for different models.

Table 5. Comparison of R for different models.

these coefficients have values of 6.77% and 99.34% for the NMJC and 4.35% and 99.54% for PTM model. In this temperature range, we observe that the PTM model performs better than the NMJC model. From these results, we can say that the proposed model is efficient. Indeed, the development of the current law aims at making the analytical models flexible by giving them the possibility to behave according to the data they have to predict. It should be noted that this model is developed from the MZA model which, according to the results we have presented, is not able to describe the behavior of this material and which, finally, has described the behavior of the material studied very well. To better understand the meaning of this model, it is better to call it an approach because it gives the overall shape of all the models that can be derived from it. We can conclude by saying that the model developed in this work applies very well to low strain rate tests. To show the relevance of this model, it was also applied in dynamics and proved to be very effective. This will be the main objective of the next subsection. The interested reader can download a Jupyter Notebook version of the identification software used for this paper in [52].

To further show the effectiveness of the model proposed in this work, we have decided to implement the PTM model into the Abaqus Explicit code. This involves carrying out a compression test on the 100Cr6 material. For this purpose, the simulation will be carried out considering only the second temperature range (900˚C - 1300˚C). This simulation will be compared with the results provided by the NMJC model for the compression test.

5. Numerical Simuations

The performance of the PTM model has now been demonstrated in terms of its ability to reproduce the experimental behavior, so it is important to evaluate its suitability for numerical simulation. According to the identification results, the PTM model reproduces the experimental behavior better than the NMJC model in the second temperature range while in the first one they are almost identical. We therefore wish to verify these results by numerical simulations. For these simulations, we will focus on a quasi-static compression test and necking of a circular bar using the 100Cr6 material. Running the circular bar bonding benchmark will allow us to better see the main difference between the two models, even though it was shown in the identification step. The PTM and NMJC behavior models are not natively available in the Abaqus Explicit FEM code, so they have been implemented by means of a VUHARD user subroutine in FORTRAN. The VUHARD subroutine is compiled using the GNU gfortran 9.3.0 and linked to the main Abaqus Explicit executable. All benchmarks tests have been solved using Abaqus Explicit 2021 on a ZBook laptop running Ubuntu 20.04 64-bits with 8 GiB of Ram and one 4 core i7-10510U Intel Processor. All computations have been done using the double precision option of Abaqus, with only one core (no parallel execution).

Cylindrical compression test

The purpose of the compression test we perform here is to compare the numerical results from using the NMJC model flow law and the PTM flow law proposed in this study. The material used in this test is the 100Cr6 steel presented in Section 4. This numerical simulation considers an axisymmetric half model of the specimen, with dimensions of 12 mm for the height and 8 mm for the diameter, as shown in Figure 6. The mesh consists of 400 CAX4R elements (4-node axisymmetric bi-linear quadrilateral elements with reduced integration and hourglass control). The FEM model is an explicit dynamic analysis, and the total simulation time is set to $t=\mathrm{17\; <mtext>\; \hspace{0.17em}\; </mtext>}\text{s}$ . We imposed a total displacement of 1.7 mm of the top crossbeam during compression, along the vertical axis of the specimen. The radial displacement of the top face of the specimen $u_r$ is kept free while its axial displacement $u_z$ is imposed and for the bottom face, the radial displacement $u_r$ is free while its axial displacement $u_z$ is kept at zero. The initial temperature of the whole specimen is 900˚C. The mass-scaling factor has been used, in order to decrease the CPU time, in this simulation and its value is $m_s=100$ which leads to a speed-up of a factor 10.

Figure 7 shows the contour plot of the von Mises stress $\sigma$ of the deformed cylinder for two different models: the NMJC (left side of the specimen) and the PTM model (right side of the specimen), while Figure 8 shows the distribution of the equivalent plastic strain ${\epsilon }^p$ for both models. We see from these figures that the spatial distribution of the equivalent stress $\sigma$ and the equivalent plastic strain ${\epsilon }^p$ are almost the same for both models. We can note that the behavior of the two models, regarding their similarity during the identification phase, was proven in the numerical simulation. An important thing to note is the ability of the PTM model to have a behavior close to that of the NMJC model. Indeed, the NMJC model is formulated from the Johnson-Cook model (which is good for dynamic tests) while the PTM model is formulated from the modified Zerilli-Armstrong model (which is good for quasi-static tests and initially was

Figure 6. Finite element meshing if the initial shape of the specimen during the compression test.

Figure 7. Von Mises stress $\sigma$ contour-plot for the cylindrical compression test.

Figure 8. Equivalent plastic strain ${\epsilon }^p$ contour-plot for the cylindrical compression test.

not able to reproduce the experiment provided by Zhou et al. [17] ). Furthermore, the Johnson-Cook model, from which the NMJC model is derived, is an empirical model that considers the macroscopic nature of the material, whereas the MZA model considers the micro-structure of the material. Therefore, one cannot expect an exact behavior in terms of plastic deformation. Hence, the interest of this numerical simulation, which allows a more effective comparison of these two formulations, could only be observed during the identification phase. We can also note on these figures that the maximum equivalent plastic deformation is located at the center of the cylinder. Thus, we decided to plot on Figure 9 the evolution of the equivalent plastic strain ${\epsilon }^p$ (for 900˚C) with the reduction of the length of the cylinder during compression for the central element (the red element in the bottom left corner in Figure 6). In this figure, we can see that the evolution of the two models is almost similar for a reduction of 1.5 mm and beyond this distance, the difference is observed but minimally. To better appreciate this simulation, we decided to reproduce all the experimental tests performed for 0.1 s⁻¹ for a temperature range $T=\left[900\text{-}1300\right]˚\text{C}$ . Figure 10 shows the comparison of the numerical results of the two models and the experimental

Figure 9. Equivalent plastic strain ${\epsilon }^p$ displacement of the top of the specimen for the cylindrical compression test.

Figure 10. Numerical and experimental von Mises stress $\sigma$ equivalent plastic strain ${\epsilon }^p$ for the cylindrical compression test for ${\stackrel{\dot{}}{\epsilon }}^p=\mathrm{0.1\; <mtext>\; \hspace{0.17em}\; </mtext>}{\text{s}}^{-1}$ .

tests. In this figure, we can observe that the numerical models do not reproduce faithfully the experiment, but the evolution is reproduced. This can be explained by the identification that was done only for 14 plastic strain values. Even if these models seem similar, there is a difference between them, as shown in Table 6. In this table, we notice that the PTM model has a larger increment and CPU time than the NMJC model. This can be explained by the form of the model which is formulated as a combination of polynomial and exponential functions. An

Table 6. Comparison of results for the cylindrical compression test for a displacement of 1.7 mm.

important thing is the equivalent plastic strain which is close to the value that was used for the identifications ( ${\epsilon }^p=0.7$ ). These values allow us to say that our numerical results are in good agreement with the experiment.

Necking of circular bar benchmark test

The necking of a circular bar is useful to evaluate the performance of the nonlinear constitutive law. An axisymmetric half model of the specimen, previously presented in Ming et al. [40], is used as shown in Figure 11. The material used for this test case is again the 100Cr6 steel presented in Section 4. An axial displacement of $u_z=\mathrm{7\; <mtext>\; \hspace{0.17em}\; </mtext>}\text{mm}$ is imposed along the z axis on the left side of the specimen while the radial displacement $u_r$ of the same edge is assumed to remain zero. On the opposite side, the axial displacement $u_z$ is constrained while the radial displacement $u_r$ is free. The mesh consists of 1600 CAX4RT elements (4-node axisymmetric bilinear thermomechanical element with reduced integration and hourglass control) with a refined region of 800 elements on the right side over 1/3 of the total length. The FEM model is an explicit coupled temperature-displacement model (it is a coupled thermo-stress analysis where the heat and mechanical transfer solutions are obtained simultaneously by explicit coupling), and the total simulation time is fixed at $t=\mathrm{3.325\; <mtext>\; \hspace{0.17em}\; </mtext>}\text{s}$ . Since there are many elements in the mesh for this analysis, we used a mass scaling factor equal to $m_s=1000$ . Recall that the main objective of this type of test is to evaluate the stress state within the specimen.

Figure 12 shows the von Mises stress contour-plot $\sigma$ of the deformed bar for two different models: the PTM model (top side) and the NMJC model (bottom side). There is a slight difference between the two models in terms of the spatial distribution of the stress. It can be seen that on the right side and for a certain distance, the distribution of the equivalent plastic stress is the same for both models and beyond this distance, it is different. The maximum values are close. This is also observed for the contour plot of the equivalent plastic strain plotted on Figure 13. The spatial distribution of the plastic strain for both models is almost the same, and the maximum values are very close. This allows us to validate our model against the results from the identification step. This is also confirmed in Table 7 reporting a comparison of the two models for a displacement of 7 mm. From this table, we can see that using the PTM model leads to fewer increments but longer computation time. The increase in computational time of the PTM model compared to the NMJC model is related to several causes, among which, the fact that the PTM model has many parameters and a polynomial form. Since the maximum equivalent plastic strain is located at the

Figure 11. Finite element meshing of the numerical model for the necking of circular bar.

Figure 12. Von Mises stress $\sigma$ contour-plot for the necking of a circular bar test.

Figure 13. Equivalent plastic strain ${\epsilon }^p$ contour-plot for the necking of a circular bar test.

Table 7. Comparison of the results for the necking of circular bar benchmark for a displacement of 7 mm.

center of the bar on the necking side (the red element in the bottom right corner in Figure 11), we have plotted the equivalent plastic strain and equivalent plastic stress for both models for this central element as presented in Figure 14 and Figure 15 respectively. From these two figures, it can be observed that the equivalent plastic strain is the same for both models for an elongation less than 5 mm and beyond this distance, a difference can be observed. Another notable effect to mention is the maximum value of the equivalent plastic strain ${\epsilon }^p$ for both models. The maximum value of the plastic strain for the PTM model is about ${\epsilon }^p=0.7$ while it is about ${\epsilon }^p=0.8$ for the NMJC model. Therefore, the fact that the PTM model reproduces the experiment better than the NMJC model, as presented in the identification part, was also validated in the numerical simulation. More details concerning the numerical results’ comparison are presented in Table 7.

Figure 14. Equivalent plastic strain ${\epsilon }^p$ displacement for the necking of a circular bar.

Figure 15. Von mises $\sigma$ equivalent plastic strain ${\epsilon }^p$ for the necking of circular bar.

6. Conclusions

The development of more complex constitutive laws considering the non-linear behavior of materials is helpful for describing several complex phenomena (such as forming and machining) in the industrial field or other ones. In this perspective, a more general formulation of the Modified Zerilli-Armstrong model (MZA), which can therefore be adjusted according to the study to be carried out, has been proposed in this work. Once the model was formulated, its parameters were identified on 100Cr6 material for a quasi-static compression test. The validation of this identification is done by two other models (NMJC and MZA models). This comparison, linked to the identification, shows that the PTM model gives better results in terms of correlation with the experimental data than the other ones. The objective was to implement this constitutive law in the numerical code, and implementation of this model was also made. Therefore, the PTM model and the NMJC model were implemented into the FEM Abaqus/Explicit code. These implementations were performed using VUHARD user subroutines for two types of tests: a cylindrical compression test and the necking of a circular bar. The numerical results obtained from those simulations show that the PTM model gives better results. For example, the plastic deformation of the PTM model for the cylindrical test is near to the one used for identification, while it is far from that value for the NMJC model. These results allow validating the relevance of the proposed constitutive law.

The specificity of the model proposed in this work is linked to the fact that it is flexible; it can then be envisaged to evaluate its performances in cases of great deformations. From this perspective, it could be possible to perform dynamic tests. It would also be very interesting to apply the model proposed in this study to other materials such as titanium, aluminum alloys, and many others that are widely used in industrial applications such as linear friction welding (LFW).

Nomenclature and Glossary

$A_i$ PTM polynomial terms of $I_1$

$B_i$ PTM polynomial terms of $S_1$

$C_i^k$ PTM polynomial terms of $S_2$

E_ARR Average Relative Error (16)

$\bar{E}$ Mean of Experimental values

$E_i$ Experimental value

$\bar{P}$ Mean of Predicted values

$P_i$ Predicted value

R Correlation coefficient (15)

T Temperature

$T_0$ Reference temperature

$T_m$ Melting temperature

$m_s$ Mass scaling factor

n Hardening exponent of ZA model

q Number of $A_i$ terms

r Number of $B_i$ terms

s Number of $C_i$ terms (T)

t Number of $C^k$ terms ( ε p )

$I_1\left({\epsilon }^p\right)$ A-thermal term of PTM

$S_1\left({\epsilon }^p\right)$ A-thermal dynamic term of PTM

$S_2\left({\epsilon }^p\mathrm{,}T\right)$ Thermal dynamic term of PTM

${\epsilon }^p$ Equivalent plastic strain

${\stackrel{\dot{}}{\epsilon }}^p$ Equivalent plastic strain rate

${\stackrel{\dot{}}{\epsilon }}_0$ Reference strain rate

$\sigma$ Von-Mises equivalent stress

${\sigma }^y\left({\epsilon }^p\mathrm{,}{\stackrel{\dot{}}{\epsilon }}^p\mathrm{,}T\right)$ Yield stress

ANNs Artificial Neural Networks

FCC Face Centered Cubic Structure

FEM Finite Element Method

JC Johnson-Cook model

MZA Modified Zerilli-Armstrong model

NMJC New Modified Johnson-Cook

PTM Proposed model

VUHARD User hardening Abaqus subroutine

VUMAT User material Abaqus subroutine

ZA Zerilli-Armstrong model

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Trimble, D., Agarwal, A., McDonnell, D., Barron, S., Ahearne, E. and O’Donnell, G.E. (2020) Finite Element Simulation of Orthogonal Machining of Biomedical Grade CoCrMo Alloy. CIRP Journal of Manufacturing Science and Technology, 28, 8-14. https://doi.org/10.1016/j.cirpj.2020.01.008
[2]	Yuan, W.-H., Wang, H.-C., Zhang, W., Dai, B.-B., Liu, K. and Wang, Y. (2021) Particle Finite Element Method Implementation for Large Deformation Analysis Using Abaqus. Acta Geotechnica, 16, 2449-2462. https://doi.org/10.1007/s11440-020-01124-2
[3]	Li, H. and Song, Z. (2021) A Reduced-Order Finite Element Method Based on Proper Orthogonal Decomposition for the Allen-Cahn Model. Journal of Mathematical Analysis and Applications, 500, Article ID: 125103. https://doi.org/10.1016/j.jmaa.2021.125103
[4]	Liu, P., Quan, Y. and Wan, J. (2021) Finite Element Simulations of Rail Milling Based on the Modified Johnson-Cook Constitutive Model. Journal of Physics: Conference Series, 1759, Article ID: 012025. https://doi.org/10.1088/1742-6596/1759/1/012025
[5]	Miller, L., Zhou, K., Tang, J., Frame, L.D., Hebert, R.J., Narayan, L.R., Alpay, S.P., Merkouriou, A. and Kim, J. (2021) Thermomechanical Finite Element Simulation and Correlation Analysis for Orthogonal Cutting of Normalized AISI 9310 Steels. The International Journal of Advanced Manufacturing Technology, 114, 3337-3356. https://doi.org/10.1007/s00170-021-07130-2
[6]	Johnson, G.R. and Cook, W.H. (1983) A Constitutive Model and Data for Materials Subjected to Large Strains, High Strain Rates, and High Temperatures. Proceedings 7th International Symposium on Ballistics, The Hague, 19-21 April 1983, 541-547.
[7]	Johnson, G.R. and Cook, W.H. (1985) Fracture Characteristics of Three Metals Subjected to Various Strains, Strain Rates, Temperatures and Pressures. Engineering Fracture Mechanics, 21, 31-48. https://doi.org/10.1016/0013-7944(85)90052-9
[8]	Zerilli, F.J. and Armstrong, R.W. (1987) Dislocation Mechanics Based Constitutive Relations for Material Dynamics Calculations. Journal of Applied Physics, 61, 1816-1825. https://doi.org/10.1063/1.338024
[9]	Khan, A.S., Sung Suh, Y. and Kazmi, R. (2004) Quasi-Static and Dynamic Loading Responses and Constitutive Modeling of Titanium Alloys. International Journal of Plasticity, 20, 2233-2248. https://doi.org/10.1016/j.ijplas.2003.06.005
[10]	Nemat-Nasser, S. and Guo, W.-G. (2003) Thermomechanical Response of DH-36 Structural Steel over a Wide Range of Strain Rates and Temperatures. Mechanics of Materials, 35, 1023-1047. https://doi.org/10.1016/S0167-6636(02)00323-X
[11]	Jia, Z., Guan, B., Zang, Y., Wang, Y. and Lei, M. (2021) Modified Johnson-Cook Model of Aluminum Alloy 6016-T6 Sheets at Low Dynamic Strain Rates. Materials Science and Engineering: A, 820, Article ID: 141565. https://doi.org/10.1016/j.msea.2021.141565
[12]	Rule, W.K. and Jones, S.E. (1998) A Revised form for the Johnson-Cook Strengh Model. International Journal of Impact Engineering, 21, 609-624. https://doi.org/10.1016/S0734-743X(97)00081-X
[13]	Vural, M., Ravichandran, G. and Rittel, D. (2003) Large Strain Mechanical Behavior of 1018 Cold-Rolled Steel over a Wide Range of Strain Rates. Metallurgical and Materials Transactions A, 34, 2873-2885. https://doi.org/10.1007/s11661-003-0188-8
[14]	Lin, Y., Chen, X.-M. and Liu, G. (2010) A Modified Johnson-Cook Model for Tensile Behaviors of Typical High-Strength Alloy Steel. Materials Science and Engineering: A, 527, 6980-6986. https://doi.org/10.1016/j.msea.2010.07.061
[15]	Lin, Y., Li, Q.-F., Xia, Y.-C. and Li, L.-T. (2012) A Phenomenological Constitutive Model for High Temperature Flow Stress Prediction of AlCuMg Alloy. Materials Science and Engineering: A, 534, 654-662. https://doi.org/10.1016/j.msea.2011.12.023
[16]	Li, H.-Y., Wang, X.-F., Duan, J.-Y. and Liu, J.-J. (2013) A Modified Johnson-Cook Model for Elevated Temperature Flow Behavior of T24 Steel. Materials Science and Engineering: A, 577, 138-146. https://doi.org/10.1016/j.msea.2013.04.041
[17]	Zhou, Q., Ji, C. and Zhu, M.-Y. (2019) Research on Several Constitutive Models to Predict the Flow Behaviour of GCr15 Continuous Casting Bloom with Heavy Reduction. Materials Research Express, 6, 1265f2. https://doi.org/10.1088/2053-1591/ab52c2
[18]	Zhang, D.-N., Shangguan, Q.-Q., Xie, C.-J. and Liu, F. (2015) A Modified Johnson-Cook Model of Dynamic Tensile Behaviors for 7075-T6 Aluminum Alloy. Journal of Alloys and Compounds, 619, 186-194. https://doi.org/10.1016/j.jallcom.2014.09.002
[19]	Johnson, G.R. and Holmquist, T.J. (1988) Evaluation of Cylinder-Impact Test Data for Constitutive Model Constants. Journal of Applied Physics, 64, 3901-3910. https://doi.org/10.1063/1.341344
[20]	Voyiadjis, G.Z. and Abed, F.H. (2005) Microstructural Based Models for Bcc and Fcc Metals with Temperature and Strain Rate Dependency. Mechanics of Materials, 37, 355-378. https://doi.org/10.1016/j.mechmat.2004.02.003
[21]	Dey, S., Borvik, T., Hopperstad, O.S. and Langseth, M. (2007) On the Influence of Constitutive Relation in Projectile Impact of Steel Plates. International Journal of Impact Engineering, 34, 464-486. https://doi.org/10.1016/j.ijimpeng.2005.10.003
[22]	Chiou, S.-T., Cheng, W.-C. and Lee, W.-S. (2005) Strain Rate Effects on the Mechanical Properties of a FeMnAl Alloy under Dynamic Impact Deformations. Materials Science and Engineering: A, 392, 156-162. https://doi.org/10.1016/j.msea.2004.09.055
[23]	Chiou, W.S. and Liu, C.Y. (2005) Comparison of Dynamic Compressive Flow Behavior of Mild and Medium Steels over Wide Temperature Range. Metallurgical and Materials Transactions A, Physical Metallurgy and Materials Science, 36, 3175-3186. https://doi.org/10.1007/s11661-005-0088-1
[24]	Chen, C., Yin, H., Humail, I.S., Wang, Y. and Qu, X. (2007) A Comparative Study of a Back Propagation Artificial Neural Network and a Zerilli-Armstrong Model for Pure Molybdenum during Hot Deformation. International Journal of Refractory Metals and Hard Materials, 25, 411-416. https://doi.org/10.1016/j.ijrmhm.2006.11.004
[25]	Lee, W.-S. and Liu, C.-Y. (2006) The Effects of Temperature and Strain Rate on the Dynamic Flow Behaviour of Different Steels, Materials Science and Engineering: A, 426, 101-113. https://doi.org/10.1016/j.msea.2006.03.087
[26]	Samantaray, D., Mandal, S., Borah, U., Bhaduri, A.K. and Sivaprasad, P.V. (2009) A Thermo-Viscoplastic Constitutive Model to Predict Elevated-Temperature Flow Behaviour in a Titanium-Modified Austenitic Stainless Steel. Materials Science and Engineering: A, 526, 1-6. https://doi.org/10.1016/j.msea.2009.08.009
[27]	Nemat-Nasser, S. (2004) Plasticity: A Treatise on Finite Deformation of Heterogeneous Inelastic Materials.
[28]	Lennon, A.M. and Ramesh, K.T. (2004) The Influence of Crystal Structure on the Dynamic Behavior of Materials at High Temperatures. International Journal of Plasticity, 20, 269-290. https://doi.org/10.1016/S0749-6419(03)00037-8
[29]	Lennon, A.M. and Ramesh, K.T. (2017) On the Performance of Modified Zerilli-Armstrong Constitutive Model in Simulating the Metal-Cutting Process. Journal of Manufacturing Processes, 28, 253-265. https://doi.org/10.1016/j.jmapro.2017.06.011
[30]	Cheng, C. and Mahnken, R. (2021) A Modified Zerilli-Armstrong Model as the Asymmetric Visco-Plastic Part of a Multimechanism Model for Cutting Simulations. Archive of Applied Mechanics, 91, 3869-3888. https://doi.org/10.1007/s00419-021-01982-6
[31]	Gurusamy, M., Palaniappan, K., Murthy, H. and Rao, B.C. (2021) A Finite Element Study of Large Strain Extrusion Machining Using Modified Zerilli-Armstrong Constitutive Relation. Journal of Manufacturing Science and Engineering, 143, Article ID: 101004. https://doi.org/10.1115/1.4050652
[32]	Lin, Y., Zhang, J. and Zhong, J. (2008) Application of Neural Networks to Predict the Elevated Temperature Flow Behavior of a Low Alloy Steel. Computational Materials Science, 43, 752-758. https://doi.org/10.1016/j.commatsci.2008.01.039 https://www.sciencedirect.com/science/article/pii/S0927025608000761
[33]	Lu, Z., Pan, Q., Liu, X., Qin, Y., He, Y. and Cao, S. (2011) Artificial Neural Network Prediction to the Hot Compressive Deformation Behavior of AlCuMgAg Heat-Resistant Aluminum Alloy. Mechanics Research Communications, 38, 192-197. https://www.sciencedirect.com/science/article/pii/S0093641311000504 https://doi.org/10.1016/j.mechrescom.2011.02.015
[34]	Ashtiani, H.R. and Shahsavari, P. (2016) A Comparative Study on the Phenomenological and Artificial Neural Network Models to Predict Hot Deformation Behavior of AlCuMgPb Alloy. Journal of Alloys and Compounds, 687, 263-273. https://www.sciencedirect.com/science/article/pii/S0925838816312877 https://doi.org/10.1016/j.jallcom.2016.04.300
[35]	Huang, X., Zang, Y. and Guan, B. (2021) Constitutive Models and Microstructure Evolution of Ti-6Al-4 Alloy during the Hot Compressive Process. Materials Research Express, 8, Article ID: 016534. https://doi.org/10.1088/2053-1591/abdaf0
[36]	Pantalé, O., Tize Mha, P. and Tongne, A. (2022) Efficient Implementation of Non-Linear Flow Law Using Neural Network into the Abaqus Explicit FEM Code. Finite Elements in Analysis and Design, 198, Article ID: 103647. https://doi.org/10.1016/j.finel.2021.103647
[37]	Gao, C. (2007) Fe Realization of a Thermo-Visco-Plastic Constitutive Model Using VUMAT in Abaqus/Explicit Program. In: Yao, Z.H. and Yuan, M.W., Eds., Computational Mechanics, Springer, Berlin, 301. https://doi.org/10.1007/978-3-540-75999-7_101
[38]	Jansen van Rensburg, G. and Kok, S. (2012) Tutorial on State Variable Based Plasticity: An Abaqus UHARD Subroutine. Eighth South African Conference on Computational and Applied Mechanics, Johannesburg, 3-5 September 2012, 158-165. https://researchspace.csir.co.za/dspace/handle/10204/6612
[39]	Duc-Toan, N., Tien-Long, B., Dong-Won, J., Seung-Han, Y. and Young-Suk, K. (2012) A Modified Johnson-Cook Model to Predict Stress-Strain Curves of Boron Steel Sheets at Elevated and Cooling Temperatures. High Temperature Materials and Processes, 31, 37-45. https://doi.org/10.1515/htmp.2011.127
[40]	Ming, L. and Pantalé, O. (2018) An Efficient and Robust VUMAT Implementation of Elastoplastic Constitutive Laws in Abaqus/Explicit Finite Element Code. Mechanics & Industry, 19, Article No. 308. https://doi.org/10.1051/meca/2018021 https://www.mechanics-industry.org/articles/meca/abs/2018/03/mi170114/mi170114.html
[41]	Hull, D. and Bacon, D.J. (2011) Introduction to Dislocations. Elsevier, Amsterdam. https://doi.org/10.1016/B978-0-08-096672-4.00002-5
[42]	Voyiadjis, G.Z. and Abed, F.H. (2005) Microstructural Based Models for Bcc and Fcc Metals with Temperature and Strain Rate Dependency. Mechanics of Materials, 37, 1816-1825. https://doi.org/10.1016/j.mechmat.2004.02.003 https://www.sciencedirect.com/science/article/pii/S0167663604000894
[43]	Gao, F., Liu, W., Zhu, Q., Gao, Z., Misra, R.D.K., Liu, Z. and Yu, F. (2020) Flow Behaviour and Constitutive Modeling for Hot Deformation of Austenitic Stainless Steel. Materials Research Express, 7, Article ID: 116512. https://doi.org/10.1088/2053-1591/abb151
[44]	Yu, H., Shao, Z. and Ren, F. (2020) Flow Stress Modeling Using the Material Constitutive Model of Modified Zerilli-Armstrong. Journal of Physics: Conference Series, 1676, Article ID: 012079. https://doi.org/10.1088/1742-6596/1676/1/012079
[45]	Gurusamy, M.M. and Rao, B.C. (2017) On the Performance of Modified Zerilli-Armstrong Constitutive Model in Simulating the Metal-Cutting Process. Journal of Manufacturing Processes, 28, 253-265. https://doi.org/10.1016/j.jmapro.2017.06.011
[46]	Samantaray, D., Mandal, S. and Bhaduri, A.K. (2009) A Comparative Study on Johnson-Cook, Modified Zerilli-Armstrong and Arrhenius-Type Constitutive Models to Predict Elevated Temperature Flow Behaviour in Modified 9Cr1Mo Steel. Computational Materials Science, 47, 568-576. https://doi.org/10.1016/j.commatsci.2009.09.025
[47]	He, A., Xie, G., Zhang, H. and Wang, X. (2014) A Modified Zerilli-Armstrong Constitutive Model to Predict Hot Deformation Behavior of 20CrMo Alloy Steel. Materials & Design (1980-2015), 56, 122-127. https://doi.org/10.1016/j.matdes.2013.10.080
[48]	Zhan, H., Wang, G., Kent, D. and Dargusch, M. (2014) Constitutive Modelling of the Flow Behaviour of a β Titanium Alloy at High Strain Rates and Elevated Temperatures Using the Johnson-Cook and Modified Zerilli-Armstrong Models. Materials Science and Engineering: A, 612, 71-79. https://doi.org/10.1016/j.msea.2014.06.030
[49]	Gupta, A.K., Anirudh, V.K. and Singh, S.K. (2013) Constitutive Models to Predict Flow Stress in Austenitic Stainless Steel 316 at Elevated Temperatures. Materials & Design, 43, 410-418. https://doi.org/10.1016/j.matdes.2012.07.008
[50]	Newville, M., Stensitzki, T., Allen, D.B., Rawlik, M., Ingargiola, A. and Nelson, A. (2016) Lmfit: Non-Linear Least-Square Minimization and Curve-Fitting for Python, Astrophysics Source Code Library.
[51]	Samantaray, D., Mandal, S., Bhaduri, A.K., Venugopal, S. and Sivaprasad, P.V. (2011) Analysis and Mathematical Modelling of Elevated Temperature Flow Behaviour of Austenitic Stainless Steels. Materials Science and Engineering: A, 528, 1937-1943. https://doi.org/10.1016/j.msea.2010.11.011
[52]	Tize Mha, P., Pantalé, O. and Tongne, A. (2022) PTM Identification Jupyter Notebook. https://doi.org/10.5281/zenodo.6078423