Proceedings of COBEM 2005 Copyright © 2005 by ABCM 18th International Congress of Mechanical Engineering November 6-11, 2005, Ouro Preto, MG FINITE ELEMENT SIMULATION OF TIG WELDING: THERMAL ANALYSIS Alexandre Campos Bezerra acbezerra@mecanica.ufu.br Domingos Alves Rade domingos@ufu.br Américo Scotti ascotti@mecanica.ufu.br Federal University of Uberlândia – School of Mechanical Engineering, Av. João Naves de Ávila, 2121, 38400-902, Uberlândia/MG - Brazil Abstract. In welding processes, the thermal history involved is responsible for, amongst other effects, the generation of residual stresses and their detrimental consequences. Thus, several studies have been carried out in order to evaluate and minimize the occurrence of welding residual stresses. In this preliminary study, which is a part of a broader project that aims at the development of finite element-based procedures for the prediction of welding residual stresses, a thermal modeling is performed to simulate TIG welding of a rectangular plate without filler metal. The commercial code ANSYS® is used. The plate is meshed using solid elements. The modeling of the heat source is discussed. For a better match to real welding conditions, temperature dependent thermo-physical material properties are considered, which leads to a non-linear analysis. Convection and radiation heat losses are also taken into account. A transient analysis is carried out to obtain the temperature field as a function of time. The approach is applied to a stainless steel AISI 316L plate. The results obtained are shown to be in good agreement with experimental results found in the literature, which demonstrates the viability and efficiency of the finite element method for simulating welding thermal cycles. Keywords: TIG welding, numerical simulation, thermal modeling, thermal cycle, heat distribution 1. Introduction Welding is the most important joining process in the world. TIG welding (also known as GTAW – Gas Tungsten Arc Welding) is a widely used welding process and, unlike processes such as laser and electron beam, it merely causes metal fusion, leading to a very localized vaporization (Depradeux, 2004). The process is based upon the creation and maintenance of a voltaic arc between a non-consumable tungsten electrode and a metal workpiece. The arc is shielded from the environment by a continuous flow of an inert gas. The heat transferred from the arc produces a thermal cycle in the workpiece. In general, during this process, a strong heating of certain regions of the workpiece occurs, while the remaining parts of it are subjected to much lower temperatures, leading to a non-uniform and transient temperature distribution. As a result, and in view of the natural trend of expansion of the heated regions, which are constrained by the colder adjacent areas, non-uniform elastic and plastic strains are generated. These plastic strains are responsible for the major part of welding residual stresses. However, material phase transformations, which lead to volume changes (expansions and contractions), also generate residual stresses. Based on the discussion above, obtaining the transient temperature field becomes the first step for determining the residual stresses field. Several studies address this problem, in which different methodologies and assumptions are proposed to solve it. Initially, analytical methods were considered, such as the one proposed by Rosenthal (1941), which assumes a concentrated heat source. According to Depradeux (2004), these analytical models are appropriate when the melted zone is considered small with respect to the workpiece dimensions. Nevertheless, new analytical solutions have been proposed, which take into account a distributed heat source (Nguyen et al., 1999; Fassani and Trevisan, 2003). The solution of the welding thermal problem using numerical methods has been widely used lately, specially the finite element method (FEM). FEM allows taking into account temperature dependence of thermo-physical material properties (like thermal conductivity, enthalpy, specific heat and density), different modes of heat transfer with the environment (convection and radiation) and the possibility of modeling structures with complex geometries (Depradeux, 2004). Since the first simulations, such as the one carried out by Hibbitt and Marcal (1973), complexity and performance of the models have evolved, mainly due to the technological advances. During this evolution, models have progressed from two-dimensional representation of structure cross-sections (Papazoglou and Masubuchi, 1982; Free and Goff, 1989; Hong et al., 1998) and further to three-dimensional transient models, including non-linear and temperature dependent material properties (Fricke et al., 2001; Francis, 2002; Depradeux, 2004). The present work seeks to obtain the transient temperature field due to a TIG welding process without filler metal on a stainless steel AISI 316L plate. For this purpose, the finite element method is employed, using a solid 3D element. Despite the fact that the use of FEM to solve this problem is well developed, its application still needs validations and a larger diffusion. 2. Thermal analysis During TIG welding, an electric power source provides a difference of potential, U, between a tungsten electrode and the workpiece. This voltage supports the formation of an electric arc, through which passes a current, I. The total power generated is calculated as the product of the voltage U and the current I. However, due to losses that occur via different mechanisms, like convection and radiation in both the electrode and the arc, only part of this power is actually used to heat the material. The ratio between the useful power and the total power is known as arc efficiency, denoted by η. Thus, the actual heat input per unit of time can be expressed as: Q = η .U .I (1) Heat flow in welding is a non-linear phenomenon due to the dependence of material thermo-physical properties on temperature. Consequently, the heat diffusion equation is written as (Incropera and de Witt, 1990): ρ (T )c(T ) ∂T ∂T  ∂  ∂  ∂T  ∂  ∂T  = q +  K x (T )  +  K y (T )  +  K z (T )  ∂t ∂y  ∂z  ∂x  ∂x  ∂y  ∂z  (2) where ρ(T) is the mass density, c(T) is the specific heat, q is the heat generated per unit volume, Kx(T), Ky(T) e Kz(T) are the thermal conductivity coefficients in each direction, T is temperature and t is time. Taking into account that enthalpy is given by: H = ρ (T )c(T )dT ∫ (3) one obtains the enthalpic formulation of heat diffusion: ∂H ∂T  ∂  ∂  ∂T  ∂  ∂T  = q +  K x (T )  +  K y (T )  +  K z (T )  ∂t ∂y  ∂z  ∂x  ∂x  ∂y  ∂z  (4) In order to consider the fusion latent heat (phase change), the value of the heat capacity suffers an abrupt increase followed by an abrupt decrease at fusion temperature (Fig. 1a). In the enthalpic formulation, an abrupt increase on enthalpy value also occurs (Fig. 1b). In the case of pure materials, fusion occurs at a constant temperature. However, in non-pure materials (alloys, for example), fusion occurs between solidus (Ts) and liquidus (Tl) temperatures. (a) (b) Figure 1. Influence of the fusion latent heat on (a) heat capacity and (b) enthalpy. Heat losses on the workpiece surfaces by convection, qc, and radiation, qr, are introduced as boundary conditions, using the following equations (Incropera and de Witt, 1990): qc = h(T − T∞ ) ( qr = εσ T 4 − T∞4 (5) ) (6) where h is the convection coefficient, T∞ is the room temperature, ε is the emissivity of the body’s surface and σ is the Stefan-Boltzmann constant. An important topic in welding simulation is the modeling of the heat source (distribution of heat input). In general, the distribution of heat input can be classified as superficial (the contribution of the plasma is essentially considered) and volumetric (the contribution of the welding pool is also included). According to Depradeux (2004), among the most frequently used heat source models for superficial distribution, one can mention the constant distribution on a disc of radius Rd, infinite Gaussian distribution and finite Gaussian distribution on a radius Rg, which are illustrated in Fig. 2. For a volumetric distribution, finite 3D Gaussian distribution on an ellipsoid and on a double ellipsoid are the most common models (Fig. 3). The choice of a model and its features (dimensions and intensities) depends on the welding process and parameters used. Therefore, it would be useful to have some information obtained from experiments, like the size of fusion zone and/or temperature as a function of time in a few points, to enable comparison with numerical results and thus, to adjust or even change the chosen model. Figure 2. Examples of superficial heat input models. Figure 3. Examples of volumetric heat input models. 3. Numerical Modeling A finite element modeling was performed using the commercial code ANSYS®. The modeled piece is a stainless steel AISI 316L plate with dimensions of 250 × 160 × 10 mm. Welding is carried-out along the mid-line of the plate parallel to its larger dimension. The material properties as functions of temperature are presented on Tab. 1 in the range [293 K - 1473 K]. For temperature values above 1473 K, the properties are extrapolated linearly up to the fusion temperature. The values of the fusion temperature and the latent heat for the AISI 316L steel reported in the literature varies significantly from a bibliographic source to another. The values used in this study are the ones from the work of Dupas and Waeckel (1994), namely fusion temperature between 1643 K and 1698 K, and fusion latent heat of 1.7982 × 109 J/m3. According to Depradeux (2004), above fusion temperature, enthalpy is considered to vary linearly, with the same slope of the enthalpy curve in the pre-fusion region, while thermal conductivity has its value doubled to simulate convection losses and to homogenize the temperature in the welding pool. Convection and radiation losses are considered in the model, where the convection coefficient is h=10 W/m2K and the emissivity is ε=0.75. Room temperature is considered to be 301 K. The element used in model meshing is the SOLID70, which is a 3D solid element with 3D thermal conduction capability. It has eight nodes with a single degree of freedom per node (temperature). The meshed model (containing 17,574 nodes) is shown in Fig. 4. The mesh was refined in the welding zone with hexahedral elements of dimensions 2.5 × 2.5 × 2.0 mm, getting gradually coarser up to elements of dimensions 2.5 × 10.0 × 2.0 mm at the two edges parallel to the weld bead. The element SURF152 was used to include the thermal radiation effect. This element needs an extra node away from the plate to simulate radiation. Therefore, a node was created on an axis normal to the surface of the plate, passing through its midpoint, and at 1.0 m away from the plate. To this node it is ascribed the value of the room temperature. Table 1. Temperature dependent thermo-physical properties of the AISI 316L steel (Depradeux, 2004). Temperature (K) 293 373 473 573 673 773 873 973 1073 1173 1273 1473 Thermal Conductivity (W/m/K) 14.0 15.2 16.6 17.9 19.0 20.6 21.8 23.1 24.3 26.0 27.3 29.9 Enthalpy (*109 J/m3) 0.000 0.300 0.704 1.130 1.560 2.000 2.450 2.910 3.370 3.860 4.360 5.360 (a) (b) Figure 4. Finite element model of the studied plate (a) and its cross-section view (b). As suggested by Depradeux (2004), the welding parameters used were: current of 150 A, voltage of 10 V, travel speed of 1 mm/s and an efficiency of 68 %. The features of the heat distribution (such as depth and area of the region where the heat is applied, and its intensity) were adjusted in such a way to provide satisfactory match to the experimental thermal results presented by Depradeux (2004). Based on this, the heat introduced was distributed in two planes of the plate: the upper face and a parallel plane located 2.0 mm below it, as shown in Fig. 5. In this figure, it can be seen that the distribution of the heat input looks like a Gaussian distribution. (a) (b) Figure 5. Heat input distribution on the upper face (a) and a parallel plane located 2.0 mm below it (b). Simulation of the movement of the heat source was performed in the following way: as the torch moves, the heat flow is applied successively to the subsequent set of nodes. Since, at first, ANSYS® code does not perform this action automatically, it was necessary the implementation of two routines in MATLAB® to write an input file which accomplishes this task. The first routine builds the list of nodes to be passed by the heat source, and the second one sequentially applies the heat flow to nodes as function of time, according to the adopted welding velocity. Welding starts at 10 mm and finishes at 240 mm from the edge x=0. Therefore, the length of the weld bead is 230 mm. 4. Results The numerical results obtained were compared to experimental results presented by Depradeux (2004). The crosssection and the points of the plate used for comparison are illustrated in Figure 6. The cross-section is located at 95 mm from the edge x=0. This way, the heat source, which moves at 1 mm/s, passes through this section after 85 s of welding. y Cross-section Cross-section Comparison Points: 0, 20 e 35 mm from y=0 Comparison Points: 10, 20 e 30 mm from y=0 x Weld bead 95 mm 160 mm 160 mm y x Weld bead 95 mm 250 mm 250 mm (a) (b) Figure 6. Cross-section used for comparing experimental and numerical results: (a) bottom and (b) upper face. Figures 7 and 8 show numerical and experimental results related to the temperature evolution as a function of time on the bottom and upper faces of the plate, respectively. Figure 8 also shows numerical results for a point on y=0. From the results of the bottom face (Fig. 7), it is possible to verify a good agreement between numerical and experimental results, showing a maximum deviation of 5.2 % for the point on y=0, 7.2 % for the point at 20 mm from y=0 and 3.2 % for the point at 35 mm from y=0. For the upper face, a good agreement between results is confirmed as well, with deviations of 6.9 % for the point at 10 mm, 4.2 % for the point at 20 mm and 3.5 % for the point at 30 mm, all the distances relative to y=0. For both faces, it is observed that the highest deviations occurred at the beginning of the heating at each point, where temperature rises earlier in the experiment than in the simulation. This could be explained by the fact that the torch movement is discrete in the model and continuous in actual welding, that is, in the model, the torch “jumps” from a set of nodes to another at each instant of time, instead of moving continuously. Figure 7. Temperature evolution as a function of time for the bottom face. Figure 8. Temperature evolution as a function of time for upper face. Figure 9 illustrates the temperature profile along the cross-section at 95 mm from the edge x=0 as a function of time. It can be seen that far from y=0, temperature evolution becomes smoother. Note that, on the farthest point from y=0 (y=80 mm), temperature rises 15 s after the torch passes through the section. One also verifies that there is a trend of temperature stabilization over the section. Figure 9. Temperature profile throughout the cross-section x=95 mm as a function of time. Figure 10 shows the contour plots of temperature field over the plate at different instants of time. The evolution of the temperature field is illustrated up to a time instant just before welding termination. It should be noticed the border effect, mainly through the isothermal line of 400 K, which, in the case of an infinite plate, should have an elliptic shape (like other isothermals). Nonetheless, as the plate has finite dimensions, in the vicinity of the edges, heat concentration occurs due to the fact that heat transfer by convection on the edges is lesser than the transfer by conduction in the bulk of the plate. (20 s) (100 s) (200 s) (230 s) Figure 10. Contour plots of the plate temperature field at different instants of time (expressed in K). Figure 11 illustrates the fusion zones obtained numerically in this work and experimentally by Depradeux (2004). The fusion zone size is compared for a cross-section at x=205 mm. It can be verified that the fusion zone is well represented by the numerical result, showing an error of 0.79 % for the width of the welding pool and 8.78 % for the depth. Figure 11. Comparison of the fusion zones obtained numerically (contour values are the same of Fig. 10) and experimentally by Depradeux (2004). 5. Conclusions This work presented a FEM based computational simulation of a welding thermal problem. A thermal transient analysis was carried out taking into account the temperature dependence of the material properties as well heat losses by convection and radiation. In the modeling, some hypotheses and simplifications have been adopted, among which the most relevant is the geometry of the heat input, for which a Gaussian-like distribution has been chosen. Also, a trade-off had to be established between the model accuracy and computational effort, leading to the choice of an adequate mesh. One of the major difficulties in this type of simulation is the lack of information about the values of material properties as functions of temperature. The numerical results, presented in terms of the temperature distributions and geometry of the fusion zone, showed good agreement with experimental results reported in the literature. This fact enables to validate the modeling procedure and the hypotheses assumed and confirms the efficiency and feasibility of the FEM for simulation of welding thermal problems. In the context of the research project currently carried out by the authors, the main interest is the calculation of welding-induced residual stresses via structural analysis. This topic is addressed in the second part of this work, presented as a companion paper. 6. Acknowledgements The authors wish to thank the agencies CNPq, of the Brazilian Ministry of Science and Technology, and CAPES, of the Brazilian Ministry of Education, for the grant of doctorate and research scholarships. The support of agency FAPEMIG through Project TEC1634/04 is greatly appreciated and acknowledged. 7. References Depradeux, L., 2004, “Simulation Numerique du Soudage – Acier 316L – Validation sur Cas Tests de Complexite Croissante”, Ph.D Thesis, Ecole Doctorale des Sciences de L’Ingenieur de Lyon, L’Institut National des Sciences Appliquees (INSA) – Lyon. Dupas, P. and Waeckel, F., 1994, “Recueil bibliographique de caractéristiques thermomécaniques pour l’acier de cuve, les revêtements inoxydables et les alliages 182 et 600”, Rapport EDF/DER HI-74/93/097, HT-26/93/058A, 46p. Fassani, R.N.S and Trevisan, O.V., 2003, “Analytical Modeling of Multipass Welding Process with Distributed Heat Source”, Journal of the Brazilian Society of Mechanical Sciences and Engineering, Vol. XXV, No. 3, pp. 302-305. Francis, J.D., 2002, “Welding Simulations of Aluminum Alloy Joints by Finite Element Analysis”, MS Thesis, Virginia Polytechnic Institute and State University. Free, J.A. and Goff, R.F.D.P., 1989, “Predicting Residual Stresses in Multi-Pass Weldments with the Finite Element Method”, Computers and Structures, vol. 32, n. 2, pp. 365-378. Fricke, S., Keim, E. and Schmidt, J., 2001, “Numerical Weld Modeling – a Method for Calculating Weld-Induced Residual Stresses”, Nuclear Engineering and Design, vol. 206, pp. 139-150. Hibbitt, H.D. and Marcal, P.V., 1973, “A numerical, thermo-mechanical model for the welding and subsequent loading of a fabricated structure”, Computers and Structures, Vol. 3, 5, pp. 1145-1174. Hong, J.K., Tsai, C.-L. and Dong, P., 1998, “Assessment of Numerical Procedures for Residual Stress Analysis of Multipass Welds”, Welding Journal, vol. 77, n. 9, pp. 372-382. Incropera, F.P. and de Witt, D.P., 1990, “Fundamentals of Heat and Mass Transfer”, 3rd edition, John Wiley. Nguyen, N.T., Ohta, A., Matsuoka, K., Suzuki, N. and Maeda, Y., 1999, “Analytical Solution of Double-Ellipsoidal Moving Heat Source and Its Use for Evaluation of Residual Stresses in Bead-on-Plate”, International Workshop on Fracture Mechanics & Advanced Materials, Sydney University, Dec. 8-10. Papazoglou, V.J. and Masubuchi, K, 1982, “Numerical Analysis of Thermal Stresses during Welding Including Phase Transformation Effects”, Journal of Pressure Vessel Technology, vol. 104, pp. 198-203. Rosenthal, D., 1941, “Mathematical theory of heat distribution during welding and cutting”, Welding Journal, N. 20, pp. 220-234. 8. Responsibility notice The authors are the only responsible for the printed material included in this paper.