# Influence of the Constitutive Flow Law in FEM Simulation of the Radial Forging Process

Authors: Olivier Pantalé; Babacar Gueye
Journal: Journal of Engineering
 (2013)

---

## Abstract

Radial forging is a widely used forming process for manufacturing hollow products in transport industry. As the deformation of the workpiece, during the process, is a consequence of a large number of high-speed strokes, the Johnson-Cook constitutive law (taking into account the strain rate) seems to be well adapted for representing the material behavior even if the process is performed under cold conditions. But numerous contributions concerning radial forging analysis, in the literature, are based on a simple elastic-plastic formulation. As far as we know, this assumption has yet not been validated for the radial forging process. Because of the importance of the flow law in the effectiveness of the model, our purpose in this paper is to analyze the influence of the use of an elastic-viscoplastic formulation instead of an elastic-plastic one for modeling the cold radial forging process. In this paper we have selected two different laws for the simulations: the Johnson-Cook and the Ludwik ones, and we have compared the results in terms of forging force, product's thickness, strains, stresses, and CPU time. For the presented study we use an AISI 4140 steel, and we denote a fairly good agreement between the results obtained using both laws.


---

## Body

## 1. Introduction

Metal forming is a widely used tool in the industry to manufacture a large range of parts in different sectors. Researchers and engineers have always resorted to some techniques in order to improve their understanding of the process on the one hand, and on the other hand to properly control the process parameters and their influence on the produced workpieces. The radial forging process is used to reduce the diameter of tubes and bars. The final shape of the workpiece is obtained thanks to a large number of strokes achieved by four anvils radially arranged around the workpiece as shown in Figure1. The preform is maintained in the proper position and is pushed into the dies using a mandrel which grips it securely. The preform is subjected to a combination of an axial translation and a circumferential rotation.Figure 1
Radial forging process used to reduce the diameter of a bar.For those studies, analytical methods have been originally developed.  Hosford and Caddell  [1] in their contribution have detailed many of such methods. One of these methods  is  the work balance method is used to estimate the amount of force required to achieve a metal forming operation. In the proposed approach, the global work is divided into three parts: (a) an ideal work that would be required for only the shape change, (b) a friction work, and (c) an unwanted redundant work. The force required is estimated from the computation of the global energy involved in the forming operation. Other approaches have also been used such as the one proposed by Ghaei et al. [2], based on the slab analysis method, to compute the deformation versus the geometry of the die in a radial forging process. Choi et al. [3] employed the upper-bound method to forecast the forging force. Their results were compared to experimental data, and a good agreement was denoted for the proposed analytical model.  Pitt-Francis et al.  [4] have proposed a three-dimensional formulation, also based on the upper bound method, to analyze the forging process. In their contribution, Donald and Chen  [5] have developed an upper bound approach to analyze the stability in soils and rocks. Even if those methods are very useful and simple to use, they are limited as soon as we need to get more detailed results concerning the material state. For those reasons the finite element method (FEM) is now widely used to simulate and analyze the forming processes.Finite element modeling was originally used in structural mechanics for managing the parts design cycle or for evaluating the performance of an existing system. The development of computational capabilities during the last three decades allows the extension of its use to many engineering fields including the metal forming processes and in particular the radial forging. Numerous authors have used FEM to study the radial forging process. In 1992, Piela [7] has studied the applicability of the finite element method for simulating radial forging process and has proposed his own model in a latter paper in 1997 [8]. In the same period, Domblesky et al. [9, 10] have proposed their operating strategy to model multiple-pass radial forging using the finite element method. Their approach is based on a stroke by stroke axisymmetric simulation with an automatic update of the die and the workpiece positions between each stroke. Due to some convergence difficulties and an excessive CPU time, they modeled the chuckhead as a surrounded ring.  Ameli and Movahhedy  [11] have also proposed a parametric study in cold radial forging process using a finite element model where they analyzed the influence of three parameters (the axial feed, the preform thickness, and the friction coefficient) on the residual stresses. The effect of the die inlet angle and the die land length on the forging force has also been evaluated. In their contribution,  Ghaei and Movahhedy [12] have used the finite element method for the design of the dies. In a later work [13] they have proposed their modeling results in terms of axial stress distribution within the tube.As presented before, the large deformation of the preform during radial forging process results from numerous strokes of the dies. The dies are hence driven in such a way that they are subjected to a high-frequency sinusoidal motion. In such conditions (large strain and large strain rate), the Johnson-Cook dynamic constitutive law [14] seems to be more adapted than a classical power law like the Hollomon or Ludwik laws usually used in most of the proposed analytical and numerical models described above. Usually, the authors make the assumption that the material presents an elastic-plastic behavior in cold conditions; therefore, they have used a classical power law to describe the hardening. As far as we know, in the literature, there are no papers which deal with the validation of this assumption when the simulation of the radial forging process is concerned. Knowing the importance of the flow law choice on the effectiveness of the model, we propose in this paper to compare the results of the same simulated process with regard to the use of the Johnson-Cook and the Ludwik laws, respectively.In this paper, theoretical bases will be first briefly presented. The mechanical problem in its local form and the bases of the finite element method used to transform the local differential equations to a set of nonlinear equations will be analyzed in the second section of the paper. In the third section the model used in this study will be detailed, and its accuracy will be presented by comparing the predicted forging force with experimental results available in the literature. Finally we will compare the Johnson-Cook and the Ludwik constitutive laws in terms of forging force, product thickness, strains, stresses, and CPU time in the last section.


## 2. Theoretical Bases

In this section, we present the mechanical problem in its local form, the discretization of the problem and the construction of the weak form associated. The constitutive law integration based on the radial return algorithm is presented in a second part of this section with the details concerning the elastic prediction and the plastic correction. As the problem is highly nonlinear, an explicit integration scheme has been retained for the time integration of the proposed equations.


### 2.1. Governing Equations

From the knowledge of the initial configuration and the imposed boundary conditions, the mechanical problem consists in predicting the final configuration depending on the initial configuration as depicted in Figure2. As presented in Figure 2 and according to the continuum mechanics approach [15], we consider an arbitrary domain Ω with boundary Γ leading to the mechanical problem presented in (1) that must be solved at each increment.  Consider,
(1)∇→·σ+ρb→=ργ→inΩu→=u→donΓuσ·n→=t→donΓd,
where ∇→ is the divergence operator, σis Cauchy stress tensor, ρb→ is the body forces vector, ρ is the mass density of the material, γ→ is the acceleration vector, u→ is the displacement vector, n→ is the surface normal vector, · is the contraction of inner indices operator,:  is the double contractor of inner indices operator, and Γu and Γd are a partition of the domain’s boundary Γ, where displacements u→d and external loads t→d are imposed, respectively. In this state, the mechanical problem cannot be solved because of a lack of equations compared to the unknown variables. Therefore the additional equations given below are added to the previous system: (i)
the geometrical compatibility equation which is written in a general way as folows:(2)ε=12[∇→u→+(∇→u→)T+∇→u→·(∇→u→)T],(ii)
the constitutive equation used to represent the material behavior as a relation between different variablesf(σ,ε,ε˙,T)=0,  where εis the strain tensor ε˙is the strain rate tensor,   and T is the temperature.Figure 2
Initial and final configuration during a transformation.


### 2.2. Finite Element Discretization

The finite element method is neither more or less than a mathematical way to resolve differential equations. It is an approximate method based on the discretization of the problem’s equations and the domain in which a solution is looked for. In radial forming, the mechanical problem is given by (1). Before resolution, this equation is turned into a weak form by multiplying (1) with an admissible virtual displacement δu→ and integrating in the hole domain Ω. So we obtain the following form on the whole domain:
(3)∫Ω‍δu→ρb→dΩ+∫Γd‍δu→t→ddΓ-∫Ω‍δε:σdΩ=∫Ω‍δu→ρu→¨dΩ.The finite element method allows writing the approximation of the displacement vectoru→e for each element of the decomposed domain Ωe and the ponderation vector δu→e using the following form: (i)
approximation of the displacement:u→e=Neu→, where Ne is the matrix of the shape functions and u→ is the displacement vector of the nodes involved in the current element e,(ii)
approximation of the admissible displacement:δu→e=Neδu→.Using the above approximation, (1) leads to the following system:
(4)Meu→¨e+F→eext=F→eint,
where, (i)
Me=∫ΩeρNeTNedΩe is the elementary mass matrix.(ii)
F→eint=∫ΩeBeTσdΩe is the elementary internal force vector.(iii)
F→eext=∫ΩeρNeTb→dΩe+∫ΓedNet→ddΓe is the elementary external force vector.The global problem is obtained by assembling the matrices above, and subsequently the system is resolved using an iterative method such as the Newton-Raphson algorithm.


### 2.3. Constitutive Law

The problem to be solved here is a large strain simulation; therefore we must ensure the objectivity of all the terms appearing in the constitutive law in order to have correct responses of the model. It is therefore necessary to maintain correct rotational transformation properties all along a finite time step. The symmetric part of the spatial velocity gradientL, denoted by D, is objective while its skew-symmetric part W, called the spin tensor, is not objective as reported, for example, in [16]. Assuming that the Cauchy stress tensor σis objective, its classical material time derivative σ˙is nonobjective, so one must introduce an objective rate notion σ∇ which is a modified time derivative form of the Cauchy stress tensor as the Jaumann or the Green-Naghdi rates. The incremental formulation of the constitutive law is therefore given by σ∇=f(D,…).One of the solutions to this problem consists  in transporting  the Cauchy stressσin a corotational frame defined using a rotation tensor wbased on the following set of equations:
(5)w˙=ωw,w(t=t0)=I.
Defining any quantity () in this rotating frame as a corotational one denoted by ()c, one may obtain in these axes, when transformed by w, the following expressions for the Cauchy stresses:
(6)σc=wTσw,σ˙c=wTσ˙w.
In fact, the choice of ω=W corresponds to the Jaumann rate. The major consequence of corotational rates is that if we choose the local axis system as the corotational one, constitutive laws integration can be performed as in small deformation, leading to a simplified formulation. According to the decomposition of the Cauchy stress tensor σinto a deviatoric term sand a hydrostatic term p, one may obtain:
(7)s˙c=ℂc:Dc,p.=Ktr[Dc],
where K is the bulk modulus of the material, Dc=wTDw,ℂc=wT[wTℂw]w, and ℂ is the fourth order constitutive tensor. In this application, we use a J2 plasticity model with a nonlinear isotropic hardening law. It is generally assumed that the rate of deformation can be additively decomposed into a elastic and a inelastic parts; therefore the approach consisting in decomposing the stress computation into an elastic predictor and a plastic corrector can be used. The associated von Mises yield criterion allows the use of the radial-return mapping strategy, briefly summarized hereafter, to integrate the constitutive behavior of the material along the time increment Δt=tn+1-tn.


#### 2.3.1. Elastic Prediction

Due to the objectivity and the use of a corotational system, all the terms of the constitutive equation are co-rotational ones, so we can drop the upperscriptc in the following equations for simplicity. The predicted elastic stresses at increment  n+1 are calculated from the current known values at increment n using the Hooke’s law, according to (7), by the following equations:
(8)pn+1trial=pn+Ktr[Δe],sn+1trial=sn+2Gdev[Δe],
where Δe=ln[U] is the co-rotational strain increment tensor between increment n and increment n+1 and G is the Lamé coefficient. Hence, the von Mises criterion f(σ,σv) is then defined by
(9)f=32sn+1trial:sn+1trial-σv,
where σv is the current yield stress depending on the history of the deformation in case of a plastic behavior. Hence, if f≤0, the predicted solution is physically admissible and the whole increment is assumed to be elastic. If f>0, then a plastic correction must be taken into account in order to restore computed stresses to a physically admissible value.


#### 2.3.2. Plastic Correction

If the predicted elastic stresses do not correspond to a physically admissible state, a plastic correction has to be performed. The previously trial stresses serve as the initial condition for the so-called return mapping algorithm. This one is summarized by the following equation:(10)sn+1=sn+1trial-2GΓn,
where n=sn+1trial/sn+1trial:sn+1trial is the unit outward normal to the yield surface and Γ is the consistency parameter defined as the solution of the one scalar parameter nonlinear equation below:
(11)f(Γ)=(sn+1trial-2GΓn):(sn+1trial-2GΓn)-23(σv(Γ))2=0.Equation (11) is effectively solved by a local Newton iterative procedure [17]. Since f(Γ) is a convex function, convergence is guaranteed.


### 2.4. Time Integration

All above equations are integrated by an explicit scheme associated with lumped mass matrices. The flowchart for the explicit time integration of the Lagrangian mesh is given in Algorithm1 as proposed, for example, by Belytschko et al. [15].Algorithm 1:Flowchart for explicit time integration.
(
1
) Initial conditions and initialization: n=0; σ0=σ(t0); x→0=x→(t0); v→0=v→(t0)
(
2
) Update quantities: n←n+1; σn=σn-1; x→n=x→n-1; v→n+1/2=v→n-1/2
(
3
) Compute the time-step and update current time: tn=tn-1+Δt
(
4
) Update nodal displacements: x→n=x→n-1+Δtv→n-1/2
(
5
) Compute internal and external force vector F→nint, F→next
(
6
) Integrate the conservative equations and compute accelerations: v→.n=M-1(F→next-F→nint)
(
7
) Update nodal velocities: v→n+1/2=v→n-1/2+Δtv→.n
(
8
) Enforce essential boundary conditions: if node I on Γv
(
9
) Output; if simulation not complete go to (2).


## 2.1. Governing Equations

From the knowledge of the initial configuration and the imposed boundary conditions, the mechanical problem consists in predicting the final configuration depending on the initial configuration as depicted in Figure2. As presented in Figure 2 and according to the continuum mechanics approach [15], we consider an arbitrary domain Ω with boundary Γ leading to the mechanical problem presented in (1) that must be solved at each increment.  Consider,
(1)∇→·σ+ρb→=ργ→inΩu→=u→donΓuσ·n→=t→donΓd,
where ∇→ is the divergence operator, σis Cauchy stress tensor, ρb→ is the body forces vector, ρ is the mass density of the material, γ→ is the acceleration vector, u→ is the displacement vector, n→ is the surface normal vector, · is the contraction of inner indices operator,:  is the double contractor of inner indices operator, and Γu and Γd are a partition of the domain’s boundary Γ, where displacements u→d and external loads t→d are imposed, respectively. In this state, the mechanical problem cannot be solved because of a lack of equations compared to the unknown variables. Therefore the additional equations given below are added to the previous system: (i)
the geometrical compatibility equation which is written in a general way as folows:(2)ε=12[∇→u→+(∇→u→)T+∇→u→·(∇→u→)T],(ii)
the constitutive equation used to represent the material behavior as a relation between different variablesf(σ,ε,ε˙,T)=0,  where εis the strain tensor ε˙is the strain rate tensor,   and T is the temperature.Figure 2
Initial and final configuration during a transformation.


## 2.2. Finite Element Discretization

The finite element method is neither more or less than a mathematical way to resolve differential equations. It is an approximate method based on the discretization of the problem’s equations and the domain in which a solution is looked for. In radial forming, the mechanical problem is given by (1). Before resolution, this equation is turned into a weak form by multiplying (1) with an admissible virtual displacement δu→ and integrating in the hole domain Ω. So we obtain the following form on the whole domain:
(3)∫Ω‍δu→ρb→dΩ+∫Γd‍δu→t→ddΓ-∫Ω‍δε:σdΩ=∫Ω‍δu→ρu→¨dΩ.The finite element method allows writing the approximation of the displacement vectoru→e for each element of the decomposed domain Ωe and the ponderation vector δu→e using the following form: (i)
approximation of the displacement:u→e=Neu→, where Ne is the matrix of the shape functions and u→ is the displacement vector of the nodes involved in the current element e,(ii)
approximation of the admissible displacement:δu→e=Neδu→.Using the above approximation, (1) leads to the following system:
(4)Meu→¨e+F→eext=F→eint,
where, (i)
Me=∫ΩeρNeTNedΩe is the elementary mass matrix.(ii)
F→eint=∫ΩeBeTσdΩe is the elementary internal force vector.(iii)
F→eext=∫ΩeρNeTb→dΩe+∫ΓedNet→ddΓe is the elementary external force vector.The global problem is obtained by assembling the matrices above, and subsequently the system is resolved using an iterative method such as the Newton-Raphson algorithm.


## 2.3. Constitutive Law

The problem to be solved here is a large strain simulation; therefore we must ensure the objectivity of all the terms appearing in the constitutive law in order to have correct responses of the model. It is therefore necessary to maintain correct rotational transformation properties all along a finite time step. The symmetric part of the spatial velocity gradientL, denoted by D, is objective while its skew-symmetric part W, called the spin tensor, is not objective as reported, for example, in [16]. Assuming that the Cauchy stress tensor σis objective, its classical material time derivative σ˙is nonobjective, so one must introduce an objective rate notion σ∇ which is a modified time derivative form of the Cauchy stress tensor as the Jaumann or the Green-Naghdi rates. The incremental formulation of the constitutive law is therefore given by σ∇=f(D,…).One of the solutions to this problem consists  in transporting  the Cauchy stressσin a corotational frame defined using a rotation tensor wbased on the following set of equations:
(5)w˙=ωw,w(t=t0)=I.
Defining any quantity () in this rotating frame as a corotational one denoted by ()c, one may obtain in these axes, when transformed by w, the following expressions for the Cauchy stresses:
(6)σc=wTσw,σ˙c=wTσ˙w.
In fact, the choice of ω=W corresponds to the Jaumann rate. The major consequence of corotational rates is that if we choose the local axis system as the corotational one, constitutive laws integration can be performed as in small deformation, leading to a simplified formulation. According to the decomposition of the Cauchy stress tensor σinto a deviatoric term sand a hydrostatic term p, one may obtain:
(7)s˙c=ℂc:Dc,p.=Ktr[Dc],
where K is the bulk modulus of the material, Dc=wTDw,ℂc=wT[wTℂw]w, and ℂ is the fourth order constitutive tensor. In this application, we use a J2 plasticity model with a nonlinear isotropic hardening law. It is generally assumed that the rate of deformation can be additively decomposed into a elastic and a inelastic parts; therefore the approach consisting in decomposing the stress computation into an elastic predictor and a plastic corrector can be used. The associated von Mises yield criterion allows the use of the radial-return mapping strategy, briefly summarized hereafter, to integrate the constitutive behavior of the material along the time increment Δt=tn+1-tn.


### 2.3.1. Elastic Prediction

Due to the objectivity and the use of a corotational system, all the terms of the constitutive equation are co-rotational ones, so we can drop the upperscriptc in the following equations for simplicity. The predicted elastic stresses at increment  n+1 are calculated from the current known values at increment n using the Hooke’s law, according to (7), by the following equations:
(8)pn+1trial=pn+Ktr[Δe],sn+1trial=sn+2Gdev[Δe],
where Δe=ln[U] is the co-rotational strain increment tensor between increment n and increment n+1 and G is the Lamé coefficient. Hence, the von Mises criterion f(σ,σv) is then defined by
(9)f=32sn+1trial:sn+1trial-σv,
where σv is the current yield stress depending on the history of the deformation in case of a plastic behavior. Hence, if f≤0, the predicted solution is physically admissible and the whole increment is assumed to be elastic. If f>0, then a plastic correction must be taken into account in order to restore computed stresses to a physically admissible value.


### 2.3.2. Plastic Correction

If the predicted elastic stresses do not correspond to a physically admissible state, a plastic correction has to be performed. The previously trial stresses serve as the initial condition for the so-called return mapping algorithm. This one is summarized by the following equation:(10)sn+1=sn+1trial-2GΓn,
where n=sn+1trial/sn+1trial:sn+1trial is the unit outward normal to the yield surface and Γ is the consistency parameter defined as the solution of the one scalar parameter nonlinear equation below:
(11)f(Γ)=(sn+1trial-2GΓn):(sn+1trial-2GΓn)-23(σv(Γ))2=0.Equation (11) is effectively solved by a local Newton iterative procedure [17]. Since f(Γ) is a convex function, convergence is guaranteed.


## 2.3.1. Elastic Prediction

Due to the objectivity and the use of a corotational system, all the terms of the constitutive equation are co-rotational ones, so we can drop the upperscriptc in the following equations for simplicity. The predicted elastic stresses at increment  n+1 are calculated from the current known values at increment n using the Hooke’s law, according to (7), by the following equations:
(8)pn+1trial=pn+Ktr[Δe],sn+1trial=sn+2Gdev[Δe],
where Δe=ln[U] is the co-rotational strain increment tensor between increment n and increment n+1 and G is the Lamé coefficient. Hence, the von Mises criterion f(σ,σv) is then defined by
(9)f=32sn+1trial:sn+1trial-σv,
where σv is the current yield stress depending on the history of the deformation in case of a plastic behavior. Hence, if f≤0, the predicted solution is physically admissible and the whole increment is assumed to be elastic. If f>0, then a plastic correction must be taken into account in order to restore computed stresses to a physically admissible value.


## 2.3.2. Plastic Correction

If the predicted elastic stresses do not correspond to a physically admissible state, a plastic correction has to be performed. The previously trial stresses serve as the initial condition for the so-called return mapping algorithm. This one is summarized by the following equation:(10)sn+1=sn+1trial-2GΓn,
where n=sn+1trial/sn+1trial:sn+1trial is the unit outward normal to the yield surface and Γ is the consistency parameter defined as the solution of the one scalar parameter nonlinear equation below:
(11)f(Γ)=(sn+1trial-2GΓn):(sn+1trial-2GΓn)-23(σv(Γ))2=0.Equation (11) is effectively solved by a local Newton iterative procedure [17]. Since f(Γ) is a convex function, convergence is guaranteed.


## 2.4. Time Integration

All above equations are integrated by an explicit scheme associated with lumped mass matrices. The flowchart for the explicit time integration of the Lagrangian mesh is given in Algorithm1 as proposed, for example, by Belytschko et al. [15].Algorithm 1:Flowchart for explicit time integration.
(
1
) Initial conditions and initialization: n=0; σ0=σ(t0); x→0=x→(t0); v→0=v→(t0)
(
2
) Update quantities: n←n+1; σn=σn-1; x→n=x→n-1; v→n+1/2=v→n-1/2
(
3
) Compute the time-step and update current time: tn=tn-1+Δt
(
4
) Update nodal displacements: x→n=x→n-1+Δtv→n-1/2
(
5
) Compute internal and external force vector F→nint, F→next
(
6
) Integrate the conservative equations and compute accelerations: v→.n=M-1(F→next-F→nint)
(
7
) Update nodal velocities: v→n+1/2=v→n-1/2+Δtv→.n
(
8
) Enforce essential boundary conditions: if node I on Γv
(
9
) Output; if simulation not complete go to (2).


## 3. Finite Element Modeling and Validation

In a general way, the radial forging process does not exhibit an axial symmetry property. So the modeling should be performed under three-dimensional conditions. But this is excessively time consuming as we need a couple of weeks to complete a single simulation. For this reason we have made some assumptions that allow the use of a 2D model (see Figure3). Therefore the following points are not considered in the proposed model: (i)
the tube’s rotation,(ii)
the clearances between the hammers,(iii)
variation of the dies’ shape along their cross section.Figure 3
Axisymmetric finite element model.These assumptions lead to a significant computational cost saving while providing relatively good results as demonstrated in the literature. In this work, we used the Abaqus Explicit [18] finite element commercial code for the analysis, and the preform is modeled using CAX4R elements, which are 4 node quadrilateral elements with reduced integration and hourglass control. On the other hand, the dies are assumed to be rigid, and only the contacting surface with the tube is taken into account in the model. So going under axisymmetric conditions the die is discretized by 2 nodes axisymmetric rigid elements.In the literature the chuckhead is systematically represented, including thereby a supplementary contact problem. So the solving efficiency, in terms of time consumption, can be collapsed. In our work, the chuckhead is suppressed, and the tube feeding is supported by a coupling constraint between a part of the preform’s outer surface and a control point as illustrated in Figure3. For the die-preform interaction, we use a penalty contact formulation with a Coulomb friction coefficient μ=0.2. Furthermore, contrary to previous studies, our strategy does not consist in a stroke by stroke simulation. In fact we perform one single simulation for the entire process, taking into account the dwell time between each consecutive stroke.The validation of the proposed model is made by comparing our forecasted computed forging load with the experimental values available in the literature and proposed by Uhlig  [6]. This way of evaluating the effectiveness of the model was inspired by  Ghaei et al.  [13] and Ameli and Movahhedy [11].   The initial and the final diameter, are reported in Table 1 as well as the experimental and numerically computed loads. A reasonable good agreement (around 5%) is observed; therefore, our model can then be used to study the influence of the flow law on the results as proposed in the next section of this paper.Table 1
Comparison between experimental [6] and computed loads.
Preform diameter (mm)
Product diameter (mm)
Sizing zone (mm)
Experimental load (KN)
FEM load (KN)
15.97
13.18
18.00
172.00
181.63


## 4. Results and Discussion

We propose here to compare the influence of the use of the Ludwik or the Johnson-Cook laws on the numerical results in a radial forming simulation. For the Ludwik law (12), the equivalent plastic stress σeq depends only on the plastic strain ε¯ whereas the Johnson-Cook law (13) considers in addition the influence of the strain rate ε¯˙ and the effect of the temperature T during the forming process. According to the solicitations, the Johnson-Cook law seems to suit better to represent the material behavior. But, in the case of cold forming, the Ludwik law may be relevant and more practical as one must only identify three parameters (A,K, and  n) versus five (A,B,C,n, and m) for the Johnson-Cook law. In this section we will compare the numerical results of the simulations, such as the forging force, the predicting tube thickness, strains, and stresses, obtained by using the two distinct laws:
(12)σeq=A+Kε¯n,(13)σeq=[A+Bε¯n][1+Clnε¯˙ε¯˙0][1+(T-T0Tm-T0)m].We have employed an AISI 4140 material for this study. The Johnson-Cook and Ludwik parameters [16] are, respectively reported in Tables 2 and 3. The Ludwik parameters have been obtained using a classical traction test whereas the Johnson-Cook ones have been identified by using dynamic tests [16] on the same material. For the numerical simulation we used an initial diameter D0=16 mm, a percentage of reduction in diameter ΔD=17%, and a feeding rate Vr=0.37 mm/stroke. Concerning the numerical simulation, similar CPU times are denoted for both laws. In fact, a CPU time of 13 min 08 s is required when using the Johnson law, whereas a CPU time of 12 min 52 s is necessary for the Ludwik law, that is, a deviation of 2%.Table 2
AISI 4140 Johnson-Cook parameters.
A   (MPa)
B   (MPa)
n
C
m
806
614
0.168
0.0089
1Table 3
AISI 4140 Ludwik parameters.
A   (MPa)
K   (MPa)
n
817
699
0.156


### 4.1. Forging Force Evolution

In Figure4 the evolution of the forging force during forming is plotted for both constitutive laws. We can denote the existence of a transient and a stationary phase. From the beginning to the exit of the die land, the forging force increases because the amount of material being deformed and the frictional work increase continuously. Once the die land exit is reached, deformation and friction energies remain stationary and the load is stabilized around F=425 kN. Predicting the forging force can be very interesting as it allows knowing whether the capabilities required to manufacture a given product are available in a workshop or not. On the other hand, a parametric or optimization study can help to minimize the force and consequently to reduction cost by saving energy. For the proposed simulation, both laws exhibit the same evolution and a good concordance of results is denoted.Numerically predicted forging force and geometry.
(a)
(b)


### 4.2. Geometry and Product Thickness

Knowing that our study focuses on the radial forging without a mandrel, the inner surface of the tube can be radially deformed when the dies stroke. Consequently, the tube thickness is greater after than before deformation. In Figure4, we have plotted the variation of the product thickness versus the distance from the exit. We denote a maximum thickening of 10% with regard to the initial thickness. Furthermore we can see that this greatest thickness concerns the cylindrical extremity of the forged product. At the product end a little bit pinch is depicted. That leads to a thickness decrease of around 3.5% compared to the neighboring points. The conical region represents the transition between the undeformed cylindrical domain and the extremity. Hence, the thickness in this conical zone increases from the initial value up to the maximum one. Johnson-Cook and Ludwik diagrams follow the same trend, and the maximum deviation is about 1.12%. So, once again, we obtain a fairly good agreement between both laws.


### 4.3. Plastic Strain Distribution

In Figure5, the equivalent plastic strain after forming is shown. As we expect, plastic strain is more important at the extremity of the tube as this region has the biggest amount of diameter reduction. So the maximum strain is about 29% for both Johnson-Cook and Ludwik laws. In addition, the non deformed region is highlighted by a zero value of the plastic strain. Within the transition between the two regions, represented by the conical shape, the plastic strain increases from zero up to the maximum value. The two constitutive equations used in this document have the same results on both the contourplot strain and the numerical values​​. The deviations in absolute values denoted are almost equal to zero.Equivalent plastic strain after forming.
(a)
(b)


### 4.4. Stress Distribution

In Figure6 the von Mises stress distribution when the die strokes is depicted. The maximum value is reached under the die land and is about 1300 MPa. We can also see in Figure 6 that, in the die inlet, the change in the material flow direction causes a shift on the stress compared with the back and front neighboring regions.Von Mises stress when the die strokes.
(a)
(b)Even if the Ludwik law does not account for the strain rate, the forecasted stresses are close to those given by the Johnson-Cook law. We can see in Figure6 similar stress distributions. The comparison made up till now is thereby confirmed.


## 4.1. Forging Force Evolution

In Figure4 the evolution of the forging force during forming is plotted for both constitutive laws. We can denote the existence of a transient and a stationary phase. From the beginning to the exit of the die land, the forging force increases because the amount of material being deformed and the frictional work increase continuously. Once the die land exit is reached, deformation and friction energies remain stationary and the load is stabilized around F=425 kN. Predicting the forging force can be very interesting as it allows knowing whether the capabilities required to manufacture a given product are available in a workshop or not. On the other hand, a parametric or optimization study can help to minimize the force and consequently to reduction cost by saving energy. For the proposed simulation, both laws exhibit the same evolution and a good concordance of results is denoted.Numerically predicted forging force and geometry.
(a)
(b)


## 4.2. Geometry and Product Thickness

Knowing that our study focuses on the radial forging without a mandrel, the inner surface of the tube can be radially deformed when the dies stroke. Consequently, the tube thickness is greater after than before deformation. In Figure4, we have plotted the variation of the product thickness versus the distance from the exit. We denote a maximum thickening of 10% with regard to the initial thickness. Furthermore we can see that this greatest thickness concerns the cylindrical extremity of the forged product. At the product end a little bit pinch is depicted. That leads to a thickness decrease of around 3.5% compared to the neighboring points. The conical region represents the transition between the undeformed cylindrical domain and the extremity. Hence, the thickness in this conical zone increases from the initial value up to the maximum one. Johnson-Cook and Ludwik diagrams follow the same trend, and the maximum deviation is about 1.12%. So, once again, we obtain a fairly good agreement between both laws.


## 4.3. Plastic Strain Distribution

In Figure5, the equivalent plastic strain after forming is shown. As we expect, plastic strain is more important at the extremity of the tube as this region has the biggest amount of diameter reduction. So the maximum strain is about 29% for both Johnson-Cook and Ludwik laws. In addition, the non deformed region is highlighted by a zero value of the plastic strain. Within the transition between the two regions, represented by the conical shape, the plastic strain increases from zero up to the maximum value. The two constitutive equations used in this document have the same results on both the contourplot strain and the numerical values​​. The deviations in absolute values denoted are almost equal to zero.Equivalent plastic strain after forming.
(a)
(b)


## 4.4. Stress Distribution

In Figure6 the von Mises stress distribution when the die strokes is depicted. The maximum value is reached under the die land and is about 1300 MPa. We can also see in Figure 6 that, in the die inlet, the change in the material flow direction causes a shift on the stress compared with the back and front neighboring regions.Von Mises stress when the die strokes.
(a)
(b)Even if the Ludwik law does not account for the strain rate, the forecasted stresses are close to those given by the Johnson-Cook law. We can see in Figure6 similar stress distributions. The comparison made up till now is thereby confirmed.


## 5. Conclusion

In this work, a 2D axisymmetric finite element model was presented to study the radial forging process. Unlike some previous works, our strategy was not based on a stroke by stroke simulation, but we perform a single simulation accounting for the dwell time between two successive strokes. The comparison of the predicting forging force with experimental results available in the literature shows the accuracy of our model. The numerical model proposed in this paper has provided some insights that can be summarized as follows:(i)
the forging force evolution during the process exhibits a transient and a stationary phase. The limit between the two phases is reached when the part gets out of the die land;(ii)
when the process is performed without a mandrel, the maximum 10% thickening of the tube is denoted in a case of 17% diameter reduction.Our purpose in this paper is to compare the Johnson-Cook and the Ludwik constitutive laws in modeling the radial forging process. In the proposed approach, we have shown a very good accordance of the results between the two laws in terms geometry, forging forces, strains, and stresses. We can therefore conclude that even if the preform undergoes high-speed stroke, the strain rate is not great enough in the case of cold forming to request the use of an elasto-viscoplastic formulation. Therefore a simple Ludwik or Hollomon law can be used to represent the material behavior and lead to quite good results. Concerning the CPU time, we have denoted a deviation of 2% between the two laws with the advantage to the simplest one.



---

*Source: 231847-2013-05-08.xml*