3. Numerical model
In order to establish a comparison between simulations and experiments, this section presents the numerical model used to simulate actual welding experiments using the Finite Element Method (FEM), and taking into account initial and boundary conditions. Corresponding thermomechanical calculations using an equivalent heat source approach were then performed.
Simulation results obtained using this approach allows determining the evolution of temperature during welding, as well as to evaluate the deformations and angles of the sheets distortion resulting from the welding operations.
The "heat source identification" and "thermomechanical calculations" sections are based on numerical simulation and experimentations. The equivalent heat sources were identified for both selected end-to-end configurations.
The calculation of deformations requires identification of a behavior law for both the base material and the molten zone. These laws were obtained from the literature, [
9,
10,
11].
In the chosen approach, the stages of numerical simulation of welding involve, firstly, using the geometric and energy parameters to identify an equivalent heat source for each configuration. Secondly, simulated temperature fields during heating and cooling are compared to those obtained experimentally. The next step is to determine the displacements/deformations of the metal sheets from the thermomechanical calculations.
The Finite Element (FE) mesh size was chosen on the basis of literature comparable to our study [
12,
13,
14,
15]. In order to reduce calculation times, the mesh was defined including two transition zones 1 and 2 along the transverse direction. To optimize calculation time and accuracy, the size of elements is finer in and around the MZ and HAZ (Heat Affected Zone). On the other hand, elements become coarser away from the weld zone. In the direction of the weld axis, the mesh size is a constant 0.25 mm.
In the weld zone, eight elements are used across the sheet thickness, providing an element size of 0.125 mm. These eight elements become four when separated from the weld in the transverse direction, and the change is operated through a first transition zone (label 1). Similarly, the four elements go through a second transition (label 2) to become two elements.
Thus, the mesh becomes coarser between zone 1 and 2 based on the progressive refinement along the lateral lines outside the weld seam. The minimum mesh size is 0.125 mm and the maximum 0.5 mm.
The FE model, including two transition zones and a progressive refinement along the lateral lines, is shown in
Figure 5, it is inspired by literature [
16] and [
17].
The use of these transition zones allowed reducing the number of elements, nodes and total CPU time for each configuration, as detailed in
Table 3.
An isostatism condition was considered for the mechanical boundary conditions and C3D8RT elements were used as suggested in the literature review [
14,
18,
19]. This type of element (3D bricks with eight nodes) is suitable for thermomechanical calculations.
Within ABAQUS, the steps used in the present study are divided into three groups, the first group of steps is used to apply initial conditions and clamping, and consists of two steps, the second group is used to simulate progressive activation of elements, and the number of steps corresponds to the number of elements, e.g. 25 steps for 25 elements, the last group is used for debriding and cooling, and consists of a single step. Thermomechanical properties (Poisson's ratio, Young's modulus, yield strength, density, specific heat, thermal conductivity, latent heat of fusion, solidus and liquidus temperatures) as a function of temperature, were taken from the literature [
17,
20,
21].
Several values of the heat exchange coefficient h were applied to investigate its effect on the maximum predicted temperature. An adjustment of the heat exchange coefficient was hence done in order to provide the best fit between experimental and calculated temperatures. Finally, the selected coefficient h takes into account both convection and radiative contributions.
For this, values of h, ranging between 10 and 60 W.m-2.°C-1 were tested. The choice of welding on a copper lath was made in view of the very low thickness of the sheets to be welded. It is essential to avoid any risk that the weld pool collapses when assembling this type of sheet. Copper was also chosen because it cools the weld pool very quickly, thus avoiding oxidation problems on the reverse side of the weld seam.
A Bilinear Isotropic Strain Hardening Model was selected to set up the elastoplastic behavior. This type of model allows considering the dependence of properties (Young’s modulus 𝐸, Tangent modulus 𝐸
𝑡 and Yield stress 𝜎
𝑒) on temperature. The Plastic strain 𝜀
𝑝𝑙 was determined from literature data [
19] obtained from experimental nominal stress-nominal strain curves,
Figure 6.
This type of behavior law includes three material parameters that are necessary for the calculation: 𝐸, 𝐸𝑡 and 𝜎𝑒. The strain-hardening parameter µ represents the ratio between the 𝐸𝑡 and the 𝐸.
A sensitivity study was carried out by
Akbari Mousavi [
9], and
Brickstad [
10] in order to determine the value of µ for 304L and 316L stainless steels at different temperatures
, and to assess the influence of varying this value on the obtained final deformation. Finally, the µ ratio was set to 0.0133 for temperatures up to 800°C, and to 0.0001 above 800°C.
In the present work, the elastoplastic behavior model of 316L stainless-steel was assumed to be identical to that of 304L, since Young's moduli as a function of temperature are almost identical [
16,
17,
22].
In ABAQUS software, thermal strains are determined by the following mathematical formula:
With:
: Thermal strain;
: Average coefficient of thermal expansion at a current temperature (°C-1);
: Current temperature (°C);
: Average coefficient of thermal expansion at an initial temperature (°C-1);
: Initial temperature (°C);
: Metal reference temperature (°C).
The reference temperatures for the base metal and seam metal used in our calculations are 25°C (ambient temperature) and 1450°C (melting temperature) respectively.
The ABAQUS "NLGEOM" parameter was used to activate the large deformations option since the small deformation option was found to provide non-realistic results.
Residual deformations induced by the sheet manufacturing process,
e.g. rolling, were not taken into account. It was assumed that, at the melting point of the material, local plastic deformation are suppressed and new plastic deformation is prevented after melting (
i.e. the transition from liquid to solid state), so that the material is reconstructed and the history of the material is assumed to be erased locally under the effect of previous treatments [
23].
It was assumed that the weld seam is divided into several elements of equal length, which were activated progressively during welding.
This was simulated using the "Element Birth and Death" technique available within ABAQUS and well explained in [
24,
25,
26]. The elements are considered inactive at the beginning of welding process,
i.e. "dead" at the start of the calculation. They are then activated as the heat source progresses,
i.e. “births". This technique is adopted to simulate a real welding process, as the seam is formed by a mixture of base metal and filler metal. A quantity of filler metal material is deposited evenly as the filler wire melts. The filler metal is modeled by the progressive activation of elements of a previously meshed seam. This approach assumes that the geometry of the seam is known before the calculation is carried out (geometry of the equivalent heat source,
Figure 4 and
Table 2).
An equivalent heat source approach was adopted, with an imposed cross-section propagated along a length
Lw in the direction of the seam axis. The cross-section is made up of two disc-segments separated by an intermediate trapezoid: it is therefore a prismatic type source. The elements are progressively activated during calculations, and the source is moving linearly in time: an example is shown in
Figure 7. Note that
Lw represents the optimum value for element length, for which calculation results close to reality have been obtained.
This source was used to represent the distribution of absorbed energy across the thickness of the sheets. Temperature fields were calculated taking into account thermal convection effects. The initial temperature of the sheets was 25°C, and an overall heat exchange coefficient was taken into account, as elaborated earlier. Finally, the coupling between the thermal problem, and the mechanical one is weak.
This approach, which enables the shape of the molten zone to be represented, was implemented in a DFLUX subroutine coded in FORTRAN. It takes the history of elements into account, which means changes from solid to liquid state, and liquid to solid. It also takes the actual geometry of the seams into account, which varies strongly with the filler material addition, as suggested on
Figure 3.
A three-dimensional thermomechanical model was developed to calculate deformations for a suitable equivalent volumetric heat source. A welding speed of 0.65 m/min was considered, as well as the average real power input provided for each configuration (homogeneously distributed in the volume of the prismatic source as shown on
Figure 7.