Evaluation of Heat Sources for the Simulation of the Temperature Distribution in Gas Metal Arc Welded Joints

: Residual stresses can affect both the static strength and the fatigue endurance of welded joints. Residual stresses can be assessed by numerical simulation; however, the simulation of the welding process is a complex task that requires knowledge of several parameters, many of which can only be estimated with some uncertainty. The reduction in the number of these parameters can lead to a more feasible and efﬁcient study. In this work, the ﬁnite element method is used to assess the capability of different thermal methods used to simulate a single pass of the gas metal arc welding process in reproducing the temperature distribution around the weld. Results of the simulations are compared to experimental measurements of the surface temperature close to the welding area. The thermal techniques analyzed adopt different levels of complexity, from the basic implementation of a constant initial temperature assigned to a given material volume, to the more comprehensive and widespread Goldak’s double-ellipsoid model. The study shows that, close to the weld seam, very similar thermal behaviors can be achieved by employing each one of the analyzed methods. Secondly, considering the constant initial temperature method, the comparison between experimental measurements and numerical simulations showed a fairly good agreement, suggesting that a relatively simple method (i.e., requiring the setting of only one parameter) can be used to efﬁciently reproduce the thermal history of a welding process.

Residual stresses after welding are difficult to detect and are present in most welded structures. Most of the time, residual stresses are detrimental to the integrity of a component; they can cause accidental failure and dimensional defects in the part. Therefore, tensile residual stresses limit the load-carrying capacity and safety margins of mechanical components under operation, and an effective engineering design can only be achieved if they are adequately known.
On the other hand, compressive residual stresses are sought after to delay the initiation and propagation of cracks and thus to increase the fatigue life (e.g., the crack closure effect explained as by [5]). Micro-and macromechanisms are the origin of a stress field inside a part, e.g., inelastic deformations, temperature gradients, and phase transformation to name a few. During welding, all three of these mechanisms are significant, occurring simultaneously all through the production process. Nowadays, the use of the finite element method (FEM) is increasing in different fields of study, being also employed for the residual stresses evaluation [6]. The use of FEM has some advantages with respect to experimental methods, such as the number of physical prototypes that can be analyzed, the associated costs, and, above all, the total time that is required for a complete analysis [7]. However, it is a challenge to deal with all the parameters involved in this kind of process, so simplifications are required [8]. Usually, some input parameters are varied and then the simulation output is studied in order to evaluate their influence. In particular, for the simulation of the welding process, independently from the used software, the focus can be on the process parameters that can be set up for the joining process (e.g., welding speed, welding power, welding technology, etc.) or non-process parameters, such as materials and geometric variables (e.g., material parameters, joint geometry, elements thickness, etc.).
The focus of this work relies on the numerical aspect of the thermal simulation [9], and in which way the material thermal properties affect the results. The process parameters and the joint geometry are kept fixed, as they are both well known, while uncertainties in the material parameters are considered as being the main issue. Other authors already studied the effect of the production process; see, for example, the effect of the heat source distribution, energy input and welding speed on temperature variation as discussed by [10]. Furthermore, a comparison of different simulation techniques of the welding process are discussed to understand to what extent simplified models can be adopted to simulate the temperature distribution close to the weld. Three numerical methods are considered in the following: An initial constant temperature of the welding seam's elements, constant heat flux and the well known Goldak heat source distribution [11]. The Goldak's double-ellipsoid heat source distribution was taken as the reference welding simulation process, as it is widely used in the literature ( [12][13][14]). The constant starting temperature technique [15] and the volumetric heat flux [16] were compared with Goldak's double-ellipsoid model to appreciate the differences in the simulated thermal behavior of the three models. Finally, a sensitivity analysis was performed to rank the parameters in order of importance, thus determining the ones that exhibit the greatest influence on simulation results.

Material
The studied model refers to the tube-to-plate welded joint geometry shown in Figure 1a, with a T-weld seam, as shown in Figure 1b. A technical drawing of the specimen is shown in Figure 1c. The specimen (of steel grade S355JR) was composed of a tube with an internal diameter of 44 mm and a thickness of 10 mm, which was joined by seam welding to a square plate with side length of 190 mm and a thickness of 25 mm. Four holes with a diameter of 21 mm were drilled in the plate corners as well as one 10 mm hole in the centre of the plate to make it possible to carry out fatigue tests ( [17][18][19]). Furthermore, a reinforcement was placed inside the tube at a distance of 15 mm from the plate surface. Both of the welding procedures were performed in a single pass using gas metal arc welding (GMAW) with the process parameters shown in Table 1, with a TransPuls Synergic 4000 (FRONIUS INTERNATIONAL GMBH, Pettenbach, Austria) machine. The welding of the tube to the plate was performed manually, fixing in advance the position of the tube by means of welding spots. The welding torch was kept still, and, using a rotating table, the welding speed was adjusted directly by the operator. The welding process can be considered as steady-state in the middle part, while transient effects are more pronounced at the start and the end of the weld seam.

Thermal Analysis
The numerical simulations were performed using the Ansys c 19.2 software (Ansys Inc., Canonsburg, Pennsylvania, United States). Specifically, a transient nonlinear simulation was used for the thermal analysis. A full 3D model was used to take into account the asymmetry caused by the welding operation, which distributes heat over the specimen progressively in a hoop direction.
A SOLID70 element was adopted for the 3D analysis; it is a brick with eight nodes that uses the temperature as the only degree of freedom. The temperature-dependent thermal material properties were collected from different sources ( [16,[20][21][22][23]). In the following, the properties taken as references are the ones derived from [20]. In particular, referring to Figure 2, a fictitious increase in the coefficient of thermal conductivity was introduced, specifically, for the heat transfer in the melting pool at 2000 • C, as suggested by [24]. Meanwhile, both convection and radiation were considered constant. The heat convection coefficient was set to a value of 20 W m −2 K −1 on all the surfaces except the bottom one in order to simulate the heat transfer between steel and air. In particular, for the bottom surface, the heat sink effect was considered, as a result of contact with the metallic support during welding; the heat convection coefficient was set to a value of 500 W m −2 K −1 , which is comparable to the one used by [25]. In addition, the radiation was taken into account through the emissivity coefficient set at 0.79 at ambient temperature. Generally, the convection is significant for lower temperatures than radiation; instead, for temperatures close to the material melting point, radiation is more relevant. However, no sensitivity was detected for emissivity and convection coefficients in the set of values found in the literature, as shown in Section 4.1.
The main focus of the simulation was the welding process. Three different methods were considered and used to evaluate their capability in reproducing the experimental results. The selected methods are listed below: The thermal analysis follows the equilibrium Equation (1), as described by [26] : where Q int is the internal heat source rate, and k, c p , and ρ are the conductivity, specific heat and density of the material, respectively. Figure 3 qualitatively illustrates the three techniques applied in the FEM; in particular, the parameters of Figure 3a-c are described in the following sections. In Figure 3a, a set of elements is chosen and an initial temperature of T i is assigned to them, after which, the cooling process is simulated by a transient thermal simulation. In Figure 3b, constant volumetric heat flux is applied to the elements belonging to the weld seam. The model activation time is a function of the total welding time and the number of activated elements-while, in Figure 3c, a double-ellipsoid function is used to select the nodes depending on their coordinates with subsequent application of a volumetric heat flux whose intensity depends from the node's coordinates. It can be immediately noticed that one of the main differences among the applied methods is the number of parameters on which they depend, one for CIT, three for CVHF, and ten for GHSD. It comes out that the use of the CIT and CVHF methods can be handled more easily than the other one. Nevertheless, from a thermal point of view, the CIT has some limitations in accurately describing the thermal profile during welding. A transitory is necessary before the simulated temperature profile matches the real temperature profile, as described later in Section 3.2. In addition, the activation of elements at a fixed temperature does not allow for the study of a moving, tilted heat source [27]. Methods like CIT and CVHF are usually based on tuneable parameters, which are necessary for the correct setting of the heat source [28].

Constant Initial Temperature
The first method consisted of a sequential activation of elements with the same temperature; to do this, an "element birth and death" technique was applied. The procedure, which is now implemented in all of the most modern FEM software, deactivates the elements by multiplying their conductivity (or stiffness depending on the simulation type) for a severe reduction factor; as soon as the element is reactivated, the key quantities return to their full original values. The elements to which the CIT is applied were chosen according to the real geometry of the weld seam, as shown in Figure 4. After their activation, the elements exchange heat with the surrounding material; the transient temperature field (T), function of time (t), follows the equilibrium Equation (1). It is worth noting that the use of the "element birth and death" technique approximates the welding process in a discrete way, by activating only a finite set of elements. In this case, Q int = 0 because only the initial temperature condition was set. Another example of the application of this technique can be found in [15]. The application of this method is immediate and easy, even for complex joint geometry; actually, only the elements selection is needed and just one parameter has to be set. Nonetheless, it is worth noting that the real complex process of material deposition with phase transition (e.g., [29]) is not simulated, and, therefore, the temperature that has to be settled on has to be considered as a fictitious parameter that has no real basis.

Constant Volumetric Heat Flux
The volumetric heat flux consists of a constant uniform heat distribution that follows Equation (2): where q vol (W m −3 ) is the volumetric heat flux, U (V) is the welding voltage, I (A) is the current, η is the efficiency of the welding process, and V (m 3 ) is the volume to which the flux is applied. The main used parameters are given in Table 2. Additionally, in this case, the factor c c has no physical interpretation, but its purpose is to fine-tune the volumetric heat flux to reach the desired temperature as observed in the experimental test.
The heat flux was activated sequentially on a set of elements belonging to the welding seam through the "element birth and death" technique. Again, the simulation has a discrete nature due to the activation of specific set of elements. With reference to Equation (1), the internal heat source is given by Q int = ηU I. Different authors already applied this technique (see e.g., [16,30]).

Goldak's Heat Source Distribution
The double-ellipsoid heat source proposed by Goldak is based on the following set of equations: where Q s (W) is the heat distribution of the welding heat source, U (V) is the welding voltage and I (A) is the current of the welding process, while η is the efficiency of the process. In Equations (4) and (5), q f and q r (W m −3 ) are, respectively, the front and rear ellipsoidal heat distribution amount at point (x,y,z), while a, b, C f and C r are geometrical parameters of the ellipse. The geometrical parameters were collected directly from the fusion zone, as shown in Figures 5 and 6; in this way, the ellipsoid geometry is experimentally described, as was already done by [31]. Finally, the terms f f and f r are the strength distributions of the ellipsoid. For the application of the Goldak model, it is necessary to define a coordinate system, which, for instance, has to keep track of the welding torch movement; in this way, the Goldak distribution can replicate the real heat distribution without using the "element birth and death" technique. For the analyzed case, the coordinate system required for this particular joint geometry is represented in  The two systems (X, Y, Z) and (x, y, z) of Figure 5 can be linked by the homogeneous transformation matrix written in the global reference system: In Matrix (8), t is the welding time, ω is the angular rotation speed of the welding torch with respect to the Y axes, and D T and L S are geometrical parameters, representing the external tube diameter and welding width (whose magnitude is given in Figure 1), respectively.

Weld Pool Shape
All three models were validated by comparing the molten pool cross section between the real specimen and simulation results. The reference melting temperature of the material is 1773.15 K.    On the other hand, it could be expected that the CIT method, being an extremely simplified model, would not provide an accurate description of the the temperature distribution. For instance, a comparison between the CIT method and the weld pool geometry is shown in Figure 9. Clearly, the very physics of the CIT model make it impossible to correctly represent the weld pool geometry at the element activation time illustrated in Figure 4. Nonetheless, the geometry of the melted zone is well represented after a few tenths of a second after the element activation, Figure 9b. As for the GHSD and CVHF, and even for the CIT method, the characteristic temperature distribution, the stretched-left shape of the melt area, could be detected. This typical shape is due to the different cooling rates caused by the different pipe and plate thicknesses. This effect is also highlighted in the micrography of Figure 9c. In addition, the same distinction made for CVHF between the molten zone and deactivated elements is also shown here by a dashed line.
The comparison between the methods reasonably shows that the CIT method can be meaningfully applied even for a subsequent mechanical simulation.

Experimental Set-Up
The results of the numerical simulations were discussed and verified through experimental activity. Temperature over time was measured in different points of the plate surface by thermocouples, as shown in Figure 10. For practical experimental reasons, a simplified specimen geometry was selected instead of a real one, in order to measure the temperature over time distribution. The simplified geometry was required because the welding process of the original specimen does not allow for easy positioning of the thermocouples, as already described in Section 2. In this case, the welding process took place by keeping the sample still, with the subsequent movement of the welding torch. As in the case of the real test specimen, the welding process was also steady, except for a starting and for an ending transient phase. The experimental procedure was set up in order to use the same weld speed as in the real tube-to-plate specimen. The specimen geometry shown in Figure 11 was selected, in which the main parameters that influence the thermal behavior, i.e., material, plates and wall thickness, were kept the same as in the pipe-to-plate specimen object of this study.
Eight thermocouples of K type (Nickel/Chrome-Nickel/Aluminium) with a maximum detectable temperature of 1473 K and a diameter of 1 mm were placed at different distances from the weld seam in order to compare the numerical and experimental results at different positions, as summarized in Table 3. The temperature readings were taken simultaneously through a data logger Pico TC-08 (Pico Technology, Cambridgeshire, UK), operating at a frequency of 1 Hz, with cold junction compensation.  . Real specimen's geometry (left), simplified specimen's geometry (right). Table 3. Thermocouples' distance from the welding seam: d = distance, N • = thermocouple's number.

Sensitivity Analysis
A sensitivity analysis was performed to evaluate the influence of different input parameters on the temperature distribution over time. As an example, the temperature-time history on the plate surface, the point at 3 mm from the weld toe is discussed in the following. A different study was carried out on the parameters described below: Concerning thermal conductivity, heat capacity, convection and emissivity, the results are shown in Figure 12. For these quantities, a variation of 30% had been taken as a reference since it turned out to be the maximum difference found by comparing different literature sources: [20][21][22].
The maximum temperature measured in the cooling curve of Figure 13 was chosen as the benchmark to assess the influence of the parameter under consideration. The error was defined as shown in Equation (13): where T Max_Re f erence is the maximum temperature calculated for the reference curve, while T Max_Variation is the maximum temperature calculated on the input variate curve. Approximatively, the maximum temperatures are the ones highlighted from the dashed line in Figure 13; the image reports an example of how the cooling behavior changes with a variation in a generic input parameter. As can be seen from Figure 12, only thermal conductivity and thermal capacity influence the results. A variation of 30% can produce an output variation on the temperature of up to about 5%, while the influence of the other parameters is negligible in this range.

Numerical Results
First of all, a comparison was done between the three numerical models, as explained in Section 3. The comparison is shown in Figure 14; the Goldak's double-ellipsoid result (continuous lines in Figure 14) was taken as the reference curve, while the temperature over time was measured in three points on the surface of the plate.
The legend in Figure 14 refers to three points on the plate surface with the following distance from the weld toe: 3 mm for the point 1, 6 mm for the point 2, and 10 mm for the point 3.
Starting from the Goldak's numerical solution, the other methods were tuned in order to qualitatively have the same thermal behavior. The parameters that were changed were the initial temperature T i for the CIT and the parameter c c for the CVHF method. In particular, the setting parameters were varied in an attempt to reproduce the maximum temperature at point 1, T i = 7900 K and c c = 36 were obtained. Figure 14 shows that, near the weld seam, the same thermal behavior can be achieved with all the implemented methods. Some subtle difference can be observed inside the weld seam as already discussed in Section 3.5. For instance, with reference to the CIT method, it is assumed that the total heat input derives from the molten metal droplets, whereas other heat sources also play an important role (e.g., the electric arc of plasma has a temperature range from 3000 K to 20,000 K, as reported in [32,33]); this leads to the slightly higher initial temperature that has to be set. Similar considerations apply to the CVHF, where only a few parameters were taken into account compared to GHSD. Obviously, it is very complicated to consider all the parameters that can influence the thermal behavior during welding; these numerical models can only take into account some of them. Another example is given by the diverse heat source distribution caused by a different shielding gas composition, as illustrated by [34].

Comparison with Experimental Temperature Measurements
The CIT method was then used in an attempt to reproduce the experimental temperature measurements. The initial temperature was varied until a value of T i = 8500 K, while the number of elements to be activated was determined by the deposition velocity of the weld material. With reference to Figure 15, the calibration on the CIT method with respect to experimental measurements was performed at the point closest to the weld seam (i.e., point 1 at 9 mm). Figure 15 shows the experimental curves in a continuous line and those resulting from the simulation are shown as a dashed line.
Distance from weld seam: A fairly good agreement was observed between the experimental and simulated temperature; the observed maximum difference in temperature during the cooling phase is less than 50 K. With reference to the heating-up phase, the measurements show the typical wave-form of manual welding; for this reason, greater temperature differences were sometimes obtained. However, the temperature gradient during the first phase was captured fairly well for all the curves.

Conclusions
In this study, the finite element method was employed to conduct a thermal numerical analysis for a welded T-joint. Three different methods for the simulation of thermal heat source were applied, i.e., constant initial temperature assigned to a given volume of elements, constant heat flux, and the Goldak's double-ellipsoid model. The methods require a different number of parameters to be set for the analysis. First of all, a sensitivity analysis on several parameters showed that thermal conductivity and heat capacity have a noticeable effect, while the influence of the other parameters on the output is negligible. Thereafter, it was shown that even the most simple numerical model (CIT) can be used to reasonably reproduce the thermal history close to the weld; with a proper choice of the initial temperature, the results are very similar to those obtained by the most sophisticated Goldak's double-ellipsoid model, whose accuracy has been discussed in different studies.
In general, the results show that a negligible discrepancy in the temperature history can be obtained using the different models if the governing parameters are properly set. Consequently, the CIT model was adopted as the main procedure to simulate the welding process.
The use of the CIT method reduces the number of parameters that have to be established, leading to a time-saving thermal study that describes a fairly similar thermal behavior as compared to more complex methods. It must be noted that the initial temperature described in this paper, as well as the other methods, only applies for the specific case analyzed in this work; different setting-parameters have to be determined for different welded joints.
Even if this study focused on the thermal solution of the welding process, it can be argued that the method can also be applied for the simulation of residual stresses.