An efficient strategy for 3D numerical simulation of friction stir welding process of pure copper plates

Copper and its alloys have a wide spectrum of engineering applications such as heat exchangers, hot water tanks or nuclear pressure vessels. Most of these structures are obtained by welding. Unfortunately, the use of conventional arc welding process is affected by several factors such as the thermal conductivity of the alloy being welded, the shielding gas, the joint design, the welding position, and the surface condition and its cleanliness. Friction stir welding process could be an interesting alternative as it can be performed without melting the material, it involves a non-consumable tool, and provides good mechanical properties. To understand in depth both the physical and the thermal mechanisms involved in this process, numerical modelling is essential. The aim of this paper is to propose an efficient simulation strategy based on the coupled Eulerian Lagrangian finite element method. The mass scaling procedure, which is used to decrease the computation time will be presented, as well as its effect on the temperature field distribution and on the down force. This model will then be used for a parametric study in order to improve the friction welding process’s parameters.


Introduction
Conventional fusion welding processes cause a series of defects in the materials to be joined, due to the excessive rise in temperature. As a result, most assemblies of aeronautical structures are made using riveting methods. This increases the complexity and cost of the processes, the amount of material used and the final weight of the structure. Friction stir welding (FSW) is a welding technique proposed by The Welding Institute in 1991 [1], which makes it possible to assemble alloys (such as copper, aluminum, titanium alloys or hybrid combinations: aluminum with copper) difficult to weld using conventional welding techniques.
The FSW process comprise, generally, four stages: plunging, dwelling, translational movement and retreating. In the first stage, the cylindrical active element (the tool), with a rotational movement and vertical translation, penetrates into the material, in the area of the welding joint, until its shoulder comes in contact with the welding parts. During dwelling stage, the active element is held in rotational position. In the translational movement stage, the tool has a rotational motion and a translational movement along the joint line (figure 1). The FSW process is completed by the vertical retraction of the tool (retreating stage). The rotation of the tool and the contact of its shoulder with the surface of the welding parts leads to the frictional generation of an important amount of heat and, thus, to increase the temperature in the joint area. The increase of the temperature plasticises the material, still remaining in a solid state, at about 70 -90% of the melting temperature and, together with the pressure exerted by the tool, allows the material to be mechanically mixed. The complex flow generated [2] produces a "solid state" joint composed of small equiaxial grains which have good mechanical properties. Figure 1. A scheme of the FSW process [3].
The advantages of the FSW process could be grouped [4] into three main categories: (i) metallurgical: (solid phase process, low part distortion, good metallurgical properties in the welded joint, small grain size), (ii) environmental: (elimination of gases, elimination of the need to clean surfaces, elimination of added materials) and (iii) energy: (reduction in the amount of material used, reduced energy consumption, reduction in weight of the structures obtained).
Despite the many industrial applications of FSW already in place, its development remains largely empirical and based on great experimental knowledge of the process. It is in this direction that numerical simulation can play a decisive role in facilitating the understanding and development of the process. Numerical observation of the phenomena that takes place inside the joint during the passage of the tool will allow a more in-depth study of the sensitivity of the metallurgical characteristics to the process parameters. One of the first attempts to model the FSW process is the analytical model of Rusell et al. [5] where the heat generated is estimated analytically from a simple equation which depends on the kinematic parameters and the properties of the material. Later, the finite element method was widely used by scientists and engineers to deal with a large number of problems associated with friction stir welding [6]. This method offers two possible approaches: (i) the Eulerian approach: mainly adapted to the simulation of flows but they manage with difficulty the conservation of mass at interfaces and extremely long computation times and (ii) the Lagrangian approach which do not suffer from this problem of mass conservation but they must be intensively remeshed in order to adapt it to the flows, which generates numerical errors and is difficult to overcome in order to complete the simulations of the FSW process.
The CEL (Coupled Eulerian-Lagrangian) method includes both benefits and advantages of the Lagrange method and of the Eulerian method, being one of the most used in numerical simulations. The CEL method allows to manage significant deformations, to predict the distribution of the temperature and to track the material flow (including the occurrence of defects during the FSW process) [7][8][9][10][11]. M A Constantin [12] proposed a simulation strategy based on an CEL approach to study the FSW process of two copper sheets. The main criticism that can be made to this modeling strategy is the extremely long computation time. However, friction stir welding simulations for complex configurations lead to extremely small and stable time increments with consequently a very long simulation time (few weeks or even few months). To reduce this computation time, several strategies are possible such as, mass scaling [13] (additions of fictitious masses) or time scaling [14] ModTech 2020 IOP Conf. Series: Materials Science and Engineering 916 (2020) 012021 IOP Publishing doi:10.1088/1757-899X/916/1/012021 3 (reduction of the simulation time). The aim of this work is to propose an efficient simulation strategy with a reasonable computation time, that can be used as parametric studies.

Finite elements Model
The model used in this work was developed in ABAQUS/Explicit environment and has been explained in detail in paper [12]. Therefore, in this paper only the main elements used in the construction of the model are presented.

Geometry of the parts
The geometry of the parts used in the numerical model is shown in figure 2. The finite element model consists of four parts: a tool, two parts (work-pieces) that represent the two welding plates and a Eulerian domain. The work-pieces are 100 mm x 40 mm x 3 mm. The welding tool has a shoulder diameter of 20 mm, the pin is conical and smooth with a length of 2.8 mm. The tools are considered to be rigid elements and have been assigned reference nodes, and the work-pieces are deformable. The Eulerian domain includes the work-pieces and has a thickness of 4 mm. The simulation was carried out in two stages: stage 1 -penetration phase with plunge velocity ‫ݒ‬ ሬሬሬሬ⃗=2.8 mm/s and for a time equal to 1s, the angular velocity was 800 rpm during the whole steps, stage 2 -translational movement with the translational velocity ‫ݒ‬ ௧ ሬሬሬ⃗ =2.5 mm/s for an advanced time equal to 64s. Total time of simulation is set to 65 seconds. The values of the process parameters and those of the dimensions of the parts were the same as the ones used experimentally 2.2. Thermomechanical behaviour of the material The parts to be joined are made of Cu-DHP (phosphorus deoxidized copper: 99.9% copper, 0.015 -0.04% phosphorus). The constitutive equation for linear thermo-elastic behaviour is expressed in the following equation: were is the stress tensor, is the linear elastic strain tensor, and are the Lamé coefficients, T is the actual temperature, is the reference temperature, ࢻ is the thermal expansion coefficient and I is the identity matrix. The thermal field equation is written as: with the Thermal conductivity ݇, Density ߩ and Specific heat ܿ . The first term on the right-hand side considers the influence of the strain rate on the temperature field. The mechanical behavior is formulated by Navier's equations for thermoelectricity [15] is expressed in equation (3) where u is the displacement vector and t is the time.
The main mechanical and thermal properties of the material [16] are presented in table 1. The Inelastic Heat Fraction is 0.9. The plastic flow of the material is modeled by the Johnson-Cook's empirical law that correlates the relationship between the flow stress σ p , the effective plastic strain ߝ , the effective plastic strain rate ߝ̄ and the temperature T: with: ߝ -the normal strain rate; A, B, C, n and m are material constants, listed in table 2; T ref -the temperature at which the parameters of material A, B, n are determined; T melt -the material's melting temperature [17].

Contact between the friction surfaces
The contact between pin and parts is a "general contact", in which tangential behaviour is described by a Coulomb's friction law, with the coefficient of friction μ assumed to be 0.35. The normal behaviour is described by a "hard contact" enabling the parts to be separated after contact.

Boundary conditions
Because in CEL modeling the computation time is very high, only the plunging and the translational stages were simulated. The main process parameters are associated with the tool: for the plunging stage -the rotation movement around its axis and the movement of translation in the Y direction; for the translation stage -the rotational movement around its own axis and the translational movement in the X direction. For the Eulerian domain: -the lower surface was taken all degrees of freedom to avoid any movement.
-a convection type model was used for heat exchange, with a heat transfer coefficient h = 30 W/m 2 -ºC, on all surfaces.

Integration scheme
The direct integration scheme adopted for the resolution is explicit because the problem treated here is dominated by numerous nonlinearities such as nonlinear material behaviour, dynamic phenomena, introduction of the thermal effects or contact.

Strategies for minimizing the computation time
Explicit integration methods are called conditionally stable because they require that the time increment value is less than a critical value ‫ݐ∆‬ ௧ . The critical time increment ‫ݐ∆‬ ௧ is calculated from the properties of mass and stifness of the entire system. ‫ݐ∆‬ ௧ = ݉݅݊ ቀ ቁ, where ‫ܮ‬ is the characteristic length of each element i of the mesh and ‫ܥ‬ ௗ = ට ா ఘ is the wave propagation speed within the material.
Mass scaling and time scaling are two strategies for minimizing the computation time. In both cases, the properties must be adapted in order to keep the thermal time constant. Mass scaling strategy increase artificially the mass of elements and can be applied even when there is rate dependency, or when the parameters are time dependent. Scaling time is only possible when the parameters do not depend on time, and when rate dependency is negligible. In this work only mass scaling strategy is used. The mass in equation (3) is scaled by replacing density ߩ by the fictious density ߩ * = ߢ ߩ, with ߢ > 0. The mass scale factor ߢ has to be chosen such a way that the inertial forces represented by the right-hand side in equation (3) remain small. When substituting the density ߩ by a fictitious density ߩ * the thermal time constant in equation (2) changes. This effect can be compensated by introducing the fictitious specific heat ܿ * = ܿ ߢ ିଵ . Now, we obtained the following scaled thermo-elastic equations. Mass inertia effects could be seen explicitly in the right term of equation (6).

Numerical results
Four values of the ߢ was used: ߢ ∈ [1, 10 ଶ , 10 ଷ ܽ݊݀ 10 ସ ] in order to investigate the effect of the additional mass. In figure 3(a), we had represented the distribution of the temperature field at the end of plunging step, but also during the translational step, for a translational displacement of the tool, equal to 15 mm. The mesh and the whole temperature field distribution are shown in figure 3(b). In order to compare quantitatively the effect of the mass scaling procedure on the temperature field, the path (red line in figure 3(b)) located in the middle of the joint was used.
(a) (b) Figure 3. The mesh and the simulated temperature distribution field during the FSW joining process. (a) the temperature field in the Z symmetry plane of the FSW joint, at the end of plunging step and during translational step for ߢ = 1 (not mass scaling); (b) the temperature field distribution during the translational step. The red line represents the path used to plot the temperature (here ߢ =1).  Figure 4(a) shows the down force evolution during only the plunging step, considered as a parameter of the time and mass scale factor ߢ . An indicator for the appropriateness of the solution is the smoothness of the down force. A smooth down force, as expected, is obtained for ߢ = 1. For values of ߢ ≥ 100 we can see a degradation of the smoothness of the down force. This is explained by the presence of inertial effects due to the addition mass, used in order to decrease the critical time increment. Nevertheless, the force obtained for ߢ ≥ 100 remains, in average, rather closed to that obtained without mass scaling. So, the down force ‫ܨ‬ ௬ is relatively slightly affected by the mass scaling strategy. The distribution of the temperature field is compared in figure 4(b) for different values of ߢ along the path AB (path illustrated in figure 3(b), during the translational step at 11.05 s. Due to the extremely long computation time, this result of the simulation is the most advanced we have to date, and it is for this reason that we compare the other results, to those one. For information, we need 5.5 months to finish the simulation without mass scaling, for a period of 65 seconds with a conventional computer. It can be noted that in figure 4(b), there are differences in the distribution of the temperature, mainly in the vicinity of the pin. The increase in temperature with ߢ is explained by (i) the fact that a higher density of the material could undoubtedly generate a larger local plastic deformation because of the inertia of the material, (ii) but also greater contact pressures, which will increase the quantity of heat generated by friction.
The Temperature error obtained for simulations with values of ߢ > 1 is presented in figure 5. This error is calculated in the vicinity of the pin, where the difference in temperature is the biggest one of all. This error obtained on the temperature field is proportional with the parameter ߢ .
Finally, the simulation time has been analyzed. In figure 5, we have plotted the influence of the mass scale factor ߢ on the CPU-time. Let's notice that a significant diminution of the computation time is observed. Having said this, for a factor ߢ = 100 the simulation time is divided by 8 with an  error in the temperature field equal to 6.9%. If we accept a Temperature error slightly less than 10%, it is possible to reduce the calculation time by 25 times. Concretely for ߢ = 1 (not mass scaling) the calculation time is 4000 hours and for ߢ = 1000 we need only 162 hours.

Conclusions
In this work, an efficient simulation strategy, using finite element method of the FSW process was presented. The mass scaling strategy introduced to reduce the CPU computing time fulfils well its role. It has shown that the computation time of the finite element simulation can be reduced significantly, with a slight degradation of the obtained results. This simulation strategy will then be used mainly for parametric studies which are necessary in numerical computation for the improvement of the FSW process or to simulate complex industrial applications.