Journal o[ Materials Processing Technology, 27 (1991) 119-133 119 Elsevier 3-D finite element modelling of the forging process with automatic remeshing Thierry Coupez, Nathalie Soyris and Jean-Loup Chenot Ecole Nationale Supdrieure des Mines de Paris, Centre de Mise en Forme des Matdriaux, Sophia Antipolis, 06560 Valbonne, France Industrial Summary The purpose of this paper is to illustrate the interest of using a finite element program to simulate the forging process of industrial parts. The FORGE3code, which can simulate hot forging of industrial parts, is presented: thermo-mechanical formulation, numerical resolution. It is well known that, in an updated lagrangian approach using a convective mesh, the degeneracy of the mesh occurs very rapidly and stops the simulation. An automatic mesh generation procedure for 3-D complex geometries has been developed with which it is possible to create the initial mesh of the billet as well as to remesh it after its degeneracy. This technique enables to simulate the whole forging process of complex industrial parts using quadratic tetrahedra. In order to show the effectiveness of the method, the example of the forging of a tripod has been computed. The simulation results show that the computation can be carried out using the described remeshing procedure and that it can be applied with success to even more complex geometries. 1. Introduction During the last twenty years, more and more studies have been done to model hot forming processes. The forging industry is slowly introducing numerical simulation codes that have been written to calculate the deformation of axisymmetrical or 2-D parts. But for complex industrial parts, a 3-D code has to be used and this is under development. A lot of work has already been published on viscoplastic or rigid-plastic finite element computation [ 1-4 ]. Most of the applications concern the field of metal forming in 2-D approximations or axisymmetrical configurations [5-10]. Three-dimensional approaches were reported for the forging of very simple parts [ 11-13]. More complex industrial shapes have been presented [14-16 ]. But one of the most limitating problems seems to be the evolution of the mesh if an updated lagrangian description is used: the mesh can become so distorted that remeshing is absolutely necessary to prevent the degeneracy of the elements. A lot of work has already been devoted to the construction of meshes with optimum geometric properties, or with some degree of adaptivity to the solution [17-21]. A continuous remesh0924-0136/91/$03.50 © 1991--Elsevier Science Publishers B.V. 120 ing technique has been suggested [22] which allows to keep a smooth and adaptive mesh during the whole process. This method has been illustrated in 2-D examples with four-node linear elements [23] and in 3-D examples with cubic eight-node linear elements [24]. For 2-D cases, forging computations with fully automatic remeshing procedures have been developed [25] involving triangular elements. The problem is much more difficult for complex 3-D geometries. A modular remeshing has been presented that allows the complete forging simulation of gears [26], but this technique is not fully automatic and does not seem easy and possible to apply for complex industrial geometries. To overcome this difficulty, and automatic mesh generation and remeshing procedure has been developed, it allows to simulate the complete forging of complex geometries using ten-node quadratic elements. In the present paper the forging code FORGE3 and its mesh generation procedure are described in detail. The complete forging simulation of a tripod is given as an example to illustrate the effectiveness of the code. 2. T h e r m o - m e c h a n i e a l formulation 2.1 The mechanical and constitutive equations Mechanics The equilibrium equation, neglecting inertia and gravity forces, is written, using the Cauchy stress tensor: V'o=O Material behavior The material is assumed to be isotropic and incompressible with the condition on the velocity field v (or the strain rate tensor ~): V . v = t r &=O In hot forging conditions, the elastic effects have a small contribution in the deformation, they are neglected, and the material is assumed to obey a viscoplastic law derived from a viscoplastic potential ¢t. The deviatoric stress tensor is then related to the strain rate tensor by: 0~(v) The Norton-Hoff law has been chosen: so that 121 S_-- 2K(v/-3~-)m- 1~ The strain rate sensitivity m, characterizing a Newtonian fluid if equal to 1 or the Von Mises rigid-plastic behavior if equal to 0, is usually around 0.1-0.3 for hot metals. The temperature and the isotropic hardening are easily taken into account by introducing a relationship between the consistency K, the equivalent strain and the temperature T: K = Ko (Co + e) nex4 • 2.2 The mechanical boundary conditions A friction law, consistent with the viscoplastic constitutive equation is built from a friction potential ¢~f(Av), where Av is the difference between the velocities of the part and the die on the interface and At~t its tangential component: AU~V--U d The friction shear stress is then: 0~f(hv) 0v A friction potential similar to the Norton-Hoff one has been chosen: ~f(Av) -- a K IIA r t IIp + 1 -p+l so that ~'= -- o~KAut ]1A r t IIp - 1 where a is the friction coefficient a n d p the friction rate sensitivity parameter. There is no force on the free surfaces: f=O The contact is unilateral, a point of the part comes in contact or leaves but cannot penetrate the die; the non-penetrability condition (of the material in the die) is written: (U--Ut)'n~0 where n is the outside normal unit vector on the part. 2.3 The associated functional Let us define the part as the domain ~, Sf and S T the portions of the boundary in contact and where a surface traction T may be imposed. The variational 1'2'2 method leads to the functional solution u: J(v) which is minimum for the velocity field J(v)= fO(v)d~+ fOf(Av)dSf- f T*vdST ~2 St" ST with the constraint V. v = 0 to prescribe the incompressibility on t~. In the special case of the N o r t o n - H o f f law and the conditions previously given we write: d(v)= ~ (v/36) d~+ HA r t IlP+ldSf Q with V- v = 0 on ~2. To enforce the incompressibility requirement, the penalty method is used and the functional is rewritten: Jp(t~): 2 ~ ( ~ f 3 ( ) d~-~- A~)tHP+ldSf+½p K(V'o)2d~Q where p is a large positive constant. 2.4 The heat transfer equation Using an isotropic Fourier law for the heat flux, the differential heat transfer equation is dT pc~= V. (kVT) + I)V where Pl/is the heat production due to the viscoplastic deformation i)d=r6: ~ if r is the rate of heat conversion energy ( r ~ 1 ). For the N o r t o n - H o f f law, this gives: W= rK(~3~-)~+' 2.5 Thermal boundary conditions The boundary conditions are very general: (i) On the free surface there is convection and radiation -kVT'n=h( T - Text) with h=hcv+er(~r(T.+.Text)(T2+ Text) 2 123 where h~v is the convection coefficient, ~r the metal emissivity, ar the Stefan constant and Text the outside temperature. (ii) On the surface in contact with the dies we have both conduction with the tool and surface energy due to the friction: b -kVT'n=hcd( T - Td) -b---~aaKLI AV t IIp+ I where h ~ i s the conductive heat transfer coefficient, Ta the die temperature, b (b=x/kpc) and bd the effusivity of the part and the die respectively. 3. Numerical resolution Using an updated lagrangian method, the whole process to simulate is cut in small time steps on which the mechanical problem and the thermal one are solved separately one after another. 3.1 Resolution of the mechanical equation Using the finite element method, the unknown function V can be approximated throughout ~ at any time t, by the relationship: q V-~ 2 N i V i i=1 where q is the number of nodes of the element, Ni the shape function related to the ith nodal point and V i the nodal velocity. The minimization of the discretized functional gives a set of non linear algebraic equations which is solved by the Newton-Raphson method. An explicite eulerian integration scheme is then used: ~ ( t + At2) = ~ ( t ) + V(t)At2 3.2 Resolution of the thermal equation As it has been done for the velocity, the temperature is written: q T= ~ N i Ti i=l where T / is the temperature at the ith nodal point. The numerical discretization of the equation and its associated boundary conditions is done using a weighted residual method and the Galerkin method. This leads to the differential non linear equations set which can be written in a matrix form as: C?+ KT+ Q = 0 124 Three-level integration schemes, consistent to the second order have been chosen: T=aTt_ ~t, + (~- 2a-g) Tt + (a-- ½+ g) Tt+ At2 dT dt (1 - g ) Tt -- Tt_ at, +gTt+ At1 At2 -- Tt At2 A linearized technique has been tested, so that each non-linear parameter A (C,K,Q) is written: A*= (½ -g)At-/~tl + (½+g)At The differential equation then leads to a set of linear equations where the nodal temperatures T i are the unknowns. 3.3 Resolution algorithm The mechanical problem and the thermal one are solved independently at each time step. At time t the temperatures Tt and Tt_ ~,,, the consistency Kt and the domain ~2, are known; then the viscoplastic resolution is computed: the velocity Vt, the friction flux, the plastic work rate and the domain ~2t+At2 are obtained; at last the thermal calculation is computed: the temperature Tt+ at2 and the consistancy Kt+/,t2 are found. 4. Automatic remeshing A lot of work has been presented in the literature on the mesh generation and on the mesh quality and adaptivity improvement even in 3-D. The different applications concerned rectangular (in 2-D ) and cubic (in 3-D) elements. To make the mesh generation and the remeshing procedure fully automatic with these elements critical geometric difficulties must be overcome. The same operations with triangles and tetrahedra are easier; this allows to create a mesh whatever the complexity of the geometry. Also to maintain a good precision in the incompressibility condition, and using a convective mesh, quadratic elements are well suited. Moreover these elements give a better approximation of the curved surface all along the deformation of the part. The mesh generation procedure that is developed in the following pages gives a 3-D quadratic mesh. It is automatic in the sense that with few interactive informations from the user, it creates a mesh fitting perfectly the geometry of the part and well adapted to a forging code. 125 4.1 Mesh generation The procedure developed to create ten-node quadratic tetrahedra meshes is very general and is applicable to remesh a deformed mesh as well as to generate the initial mesh for the beginning of the simulation. In both cases, the mesh generation requires the description of the closed outside surface of the part and especially its discretization in small linear triangles. T h e n the mesh will be made automatically through three main procedures: (i) creation of a linear surface mesh, (ii) creation of a quadratic surface mesh, (iii) creation of a quadratic volume mesh. Discretization of the outside surface of the part To mesh the initial part, the user has to give the geometrical information: using a commercial CAD code, he designs the surface and represents it with linear triangles (representation no. 0 of the surface). This is an imposed step. For a remeshing, this information is deduced from the previous deformed mesh (volume-mesh no. 1 and surface-mesh no. 1 ): each quadratic triangle of the deformed surface-mesh no. 1 is cut into several smaller linear triangles as close as possible to the real surface (representation no. 1 of the surface). During this operation new nodes have been created, some of them may be inside the dies, so this representation is then modified to avoid the nodes penetrating into the dies: it gives the representation no. 2 of the surface from which the mesh will be created. Creation of a linear surface mesh The representation no. 2 of the surface is not used as a triangular linear mesh of the surface because it has too many nodes. So a new linear triangular mesh for the surface is deduced from the representation no. 2. In this step the number of nodes and triangles on the surface is reduced while maintaining a good precision on the geometry, this gives the surface-mesh no. 2. Creation of a quadratic surface mesh In this step, the linear triangles of the previous surface mesh are curved to match as precisely as possible the real surface of the part at this stage of remeshing given by the representation no. 2. This gives a mesh of the surface with six-node quadratic triangles: surface-mesh no. 3. Creation of a quadratic volume mesh Four-node linear tetrahedra are first created, fitting exactly the linear surface mesh. T h e n this volumic linear mesh is refined by adding nodes inside the volume. Finally nodes are added on the lines of each tetrahedron. Combined with the quadratic surface mesh, this gives a ten-node quadratic tetrahedral mesh: volume-mesh no. 2. 126 4.2 Remeshing procedure • I . For the complete mmulatlon of a forging pass, different stages are necessary independently of how the meshes are created: ( 1 ) the geometrical definition of the part is needed: here, it is given by a linear triangles representation (no. 0); (2) the mesh is created, see Section 4.1 (surface-mesh nos. 2 and 3, volumemesh no. 2); (3) the forging simulation is computed and stopped either due to mesh degeneracy or other criteria such as the penetration of the nodes into the dies: in our procedure the simulation does not stop on other criteria than the degeneracy of the mesh, a manual step is necessary, by visualizing the results it can be decided to remesh before; the user decides himself when he thinks he should remesh (surface-mesh no. i and volumic-mesh no. 1 ); (4) the geometrical definition of the deformed part is deduced from the simulation results (representation no. 1 ); (5) the new mesh is created, see Section 4.1 (representation no. 2, surfacemesh nos. 2 and 3, volume-mesh no. 2); (6) the state variables are interpolated from the old mesh to the new one: here, a very simple interpolation is done for each node of the new mesh using its neighbouring nodes in the old mesh. Then it is possible to carry on the forging simulation, and, by repeating the stages three to six, to finish it. 5. R e m e s h i n g o f a t r i p o d Because of the symmetry, only one sixth of the whole part has been calculated in isothermal conditions• Figure 1 gives an idea of the different steps of the calculation. Figure 1 (a) is a view of the initial mesh. This first mesh has been made using the geometrical data of the outside surface defined by 696 linear triangles which has been obtained by the CAD code CATIA (representation no. 0); this first mesh of the part has 144 boundary quadratic triangles, 195 quadratic tetrahedra and 432 nodes• After 30% reduction, some sides of elements penetrate the dies (Fig. 1 (b)) the mesh is not well adapted for further calculation, it has been decided to make a new one. The geometrical data of the outside surface (representation no. 1 and no. 2) has 2304 linear triangles: the second mesh of the part automatically generated has 110 quadratic triangles on the surface, 154 quadratic tetrahedra and 340 nodes. Note that this second mesh of the part has less nodes than the first one but gives a better precision. Figure l ( c ) is a view of this second mesh, and Fig. 1 (d) shows its deformation after several increments: the axis begins to appear. At 80% reduction a second remeshing has been made using the surface dis- 127 .... :~ i¸ (a) (c) Fig. 1. (a) Initial mesh; (b) first mesh at 30% reduction; (c) first remeshing at 30% reduction. 128 Fig. 1. (d) Intermediate deformation of the second mesh; (e) second remeshing at 80% reduction. cretization (representation no. 1 ) in 1760 linear triangles given by the second mesh of the part (upper part of Fig. 1 (e)). The third mesh and last one with 200 quadratic triangles, 294 quadratic tetrahedra and 633 nodes is shown on 129 Fig. 2.At the end of the forging simulation. (a) final geometry: thickness of the flash; (b) whole tripod. 130 Fig. 3. Evolution of the deformation of the part with the dies during forging: (a) initial position and first mesh; (b) first remeshing at 30% reduction; (c) second mesh deformation; (d) second remeshing. the lower part of Fig. 1 (e). Note on this view that the geometrical precision is excellent: the upper part is the given geometry while the meshed lower part represents the geometry resulting from the mesh generation. Figure 2 (a) shows the part at the end of the computation: the thickness of the flash between the axes is only a few percent of the part height. Only quadratic tetrahedra allow 131 Fig. 3. Evolutionof the deformationof the part with the dies duringforging: (e) almostfinished. such achievements. Figure 2 (b) is a view of the whole tripod at the end of the simulation and forging. Finally Fig. 3 gives the evolution of the part deformation between the dies during the forging: (a) initial position and first mesh; (b) first remeshing at 30% reduction; (c) second mesh deformation; (d) second remeshing; and (e) almost finished! 6. Conclusions A complete forging program F O R G E 3 has been presented here. A new, fully automatic remeshing procedure, that generates quadratic tetrahedral elements, has been described. The complete forging simulation of a tripod has been computed and analyzed. The remeshing scheme, illustrated here, enables the program to simulate the forging of any complex industrial part. Nevertheless several improvements have to be developed: introduction of remeshing criteria in the forging code, integration of the remeshing procedure in the forging code itself so that the whole simulation will be fully automatic, as it is done for example in the 2-D code FORGE2. Other industrial cases that have not been presented here have already been computed in isothermal conditions, thermomechanical resolution has now to be tested on these cases. And of course comparisons between calculation and reality will check the results validity. Acknowledgements We would like to thank the companies Peugeot, Pechiney, Forges de Courcelles, Teksid and the French Ministate de la Recherche et de la Technologie for their financial support. 132 References 1 G.Y. Goon, P.I. Poluckin, Z.P. Poluckin and B.A. Prodcowsky, The plastic deformation of metals, Metallurgica, 1968. 2 G.C. Cornfield and R.H. Johnson, Theoretical prediction of plastic flow in the hot rolling including the effect of various temperature distribution, J. Iron Steel Inst., 211 (1973) 567. 3 C.H. Lee and S. Kobayashi, New solutions to rigid-plastic deformation problems using a matrix method, Trans. ASME, J. Eng. Ind., 95 (1973) 865. 4 O.C. Zienkiewicz and P.N. Godbole, Flow of plastic and visco-plastic solids with special reference to extrusion and forming processes, Int. J. Num. Meth. Eng., 8 (1974) 3. 5 S.I. Oh, N. Rebelo and S. Kobayashi, Finite element formulation for the analysis of plastic deformation of rate sensitive materials on metal forming, Metal forming plasticity, IUTAM symposium, H. Lippmann (Ed.), 1978, pp. 273. 6 J. Perrot, J.L. Chenot, Calcul par 6ldments finis de l'dvolution d'une inclusion clans une matrice d~form~e uniform~ment, Physical non-linearities in structural analysis, IUTAM symposium, Senlis, 1980. 7 J.W.H. Price and M. Alexander, Specimen geometries predicted by computer model of high deformation forming, Int. J. Mech. Sci., 21 (1979) 417. 8 N. Rebelo, H. Rydstad and G. Schroder, Simulation by material flow in closed-die forging by model techniques and rigid-plastic FEM, J.F.T. Pittman, R.D. Wood, J.M. Alexander and O.C. Zienkiewicz (Eds.), Numerical methods in industrial forming processes, 1982, pp. 237. 9 S.I. Oh, J.J. Park, S. Kobayashi and T. Altan, Applications of FEM modeling to simulate metal flow in forging a titanium alloy engine disc, J. Eng. Ind., 105 (1983) 251. 10 Y. Germain, J.L. Chenote, P.E. Mosser, Finite element analysis of shaped lead-tin disk forgings, Proc. Numiform '86 Conf., Gothenburg, K. Mattiasson, A. Samuelson, R.D. Wood, O.C. Zienkiewicz, (Eds.), 1986, pp. 271-276. 11 K. Mori, K. Osakada, K. Nakadoi and M. Fukuda, Simulation of three-dimensional deformation in metal forming by the rigid-plastic finite element method, Adv. Tech. of plasticity, 2 (1984) 1009. 12 J.J. Park and S. Kobayashi, Three-dimensional finite element analysis of block compression, Int. J. Mech. Sci., 26 (1984) 165. 13 G. Surdon and J.L. Chenot, Finite element calculation of three dimensional hot forging, Int. J. Num. Meth. Eng., 24 (1987) 2107-2117. 14 G. Surdon and J.L. Chenot, Finite element calculation of three-dimensional hot forging, Proc. Numiform '86 Conf., Gothenburg, K. Mattiasson, A. Samuelson, R.D. Wood, O.C. Zienkiewicz, (Eds.), 1986, pp. 287-292. 15 N. Soyris, J.P. Cescutti, T. Coupez, G. Brachotte and J.L. Chenot, Three-dimensional finite element calculation of the forging of a connecting rod, J.L. Chenot and E. Ofiate, (Eds.), Modelling of metal forming processes, proceedings of the Euromech 233 colloquium, SophiaAntipolis, 1988, pp. 227-236. 16 N. Soyris, J.J. Brioist, G. Brachotte and J.L. Chenot, Three dimensional finite element calculation of the whole forging process of an automotive part, Proc. adv. technol, plasticity, Kyoto, 1990, pp. 165-170. 17 J.U. Brackbill and J.S. Saltzmann, Adaptive zoning for singular problems in two dimensions, J. Comput. Physics, 46 (1982) 342-368. 18 A.R. Diaz, N. Kikuchi and J.E. Taylor, A method of grid optimization for finite element methods, Comp. Meth. in Appl. Mech. Eng., 41 (1983) 29-45. 19 N. Kikuchi, Adaptive grid-design methods for finite elements analysis, Comp. Meth. Appl. Mech. Eng., 55 (1986) 129-160. 20 L. Demkowicz and J.T. Oden, On a moving mesh strategy based on an interpolation error estimate technique, Int. J. Eng. Sci., 24 ( 1 ) (1986) 55-68. 133 21 22 23 24 25 26 O.-C.Zienkiewicz, Y.C. Liu and G.C. Huang, Error estimation and adaptativity in flow formulation for forming problems, Int. J. Num. Meth. Eng., 25 (1988) 23. J.P. Cescutti and J.L. Chenot, A geometrical continuous remeshing procedure for application to finite element calculation of non-steady state forming process, NUMETA 87, 2nd Int. Conf. Adv. Num. Meth. Eng.: Theory and Applications, Swansea, 1987, pp. 517. J.P. Cescutti, E. Wey, J.L. Chenot and P.E. Mosser, Finite element calculation of hot forging with continuous remeshing, J.L. Chenot and E. Ofiate (Eds.), Modelling of metal forming processes, proceedings of the Euromech 233 colloquium, Sophia Antipolis, 1988, pp. 207-216. J.P. Cescutti, N. Soyris, G. Surdon and J.L. Chenot, Thermo-mechanical finite element calculation of three-dimensional hot forging with remeshing, K. Lange (Ed.), Proc. Adv. Technol. Plasticity, Stuttgart, 1987, pp. 1051-1058. C. Corsini, E. Wey, T. Coupez and J.L. Chenot, Numerical simulation of the forging process with automatic remeshing, European conference on new advances in computational structural mechanics, Giens, France, 1991. D.Y.Yang, J.H. Yoon and N.K. Lee, Modular remeshing: a practical method of 3-D remeshing in forging of complicated parts, Proc. Adv. Technol. Plasticity, Kyoto, 1990, pp. 171-178.