COMPLETE 3D HEAT AND FLUID FLOW MODELING OF KEYHOLE LASER WELDING AND METHODS TO REDUCE CALCULATION TIMES

COMPLETE 3D HEAT AND FLUID FLOW MODELING OF KEYHOLE LASER WELDING AND METHODS TO REDUCE CALCULATION

TIMES

M. COURTOIS*, M. CARIN*, P. LE MASSON* and S. GAIE**

*Univ. Bretagne Sud, UMR CNRS 6027, IRDL, F-56100 Lorient, France

**ArcelorMittal Global R&D Montataire, F-60160 Montataire, France

DOI 10.3217/978-3-85125-615-4-01

ABSTRACT

The fluid flow calculation inside the melt pool in welding processes is a complex challenge. It can be useful for defects prediction in the weld seam or to study the influence of some process parameters. The cost in time of these calculations makes these models not widely used, although they are rich in information.

The aim of this paper is to present a complete model solving the heat and fluid flow equations in all the states (solid – liquid - gas) of the metal during keyhole laser welding. In order to track dynamically the keyhole shape, a level set function is employed and the momentum equations are solved to predict accurately the melt pool behaviour. The particularity of this model is a drastic optimization to reduce computation time (less than 1 day on a commercial workstation without reducing significantly the accuracy of the fluid flow computation. To achieve this goal, a simplified approach of the recoil pressure effect is employed and an original method with three different meshes adapted for each equation and physics problem is proposed and discussed. Finally, to validate the model simplification, a complete experimental validation is added with longitudinal and transverse micrograph cuts for different welding configurations, micro-thermocouples measurements in solid phase [20-1200°C], pyrometer measurements in liquid phase [1500-3000°C].

Keywords: Laser welding; heat and fluid flow; keyhole, numerical simulation

INTRODUCTION

In automotive industry, laser welding is widely used for its high production rate and its capacity to produce clean and discreet welded beams. Moreover, the cost of the laser sources and the increased available power offer always more applications.

Laser welding of steel is currently employed in industry for 20 years. Nevertheless, the use of new materials or the necessity to weld new configurations raises sometimes difficulties to achieve a satisfying welded joint. In most cases, questions come when the weld seam presents some defects like porosity, lake of materials or partial penetration. In

all these cases, the fluid flow in the melt pool during the welding process is the responsible:

gravity with large melt pools can create collapsing, the keyhole opening and behaviour can impact the melt pool stability.

With this observation, and the necessity to understand the liquid and keyhole behaviour in all the complexes configurations, many experimental works have been carried in the last 30 years. More recently (15 years), various numerical works have been proposed to provide a better understanding when experiments were unable to give sufficient information.

Models that propose a self-estimation of the keyhole shape remains limited but numerous methods have been developed to describe the keyhole effect as summarized by Dal in 2016 [1]. The keyhole is the place where the laser energy is introduced and by its strong action, it is responsible for many movements in the molten pool. However, this free surface is difficult to calculate due to the very complicated vaporization effects. The first works that really propose a method to treat dynamically the keyhole are those of Ki in 2002 [2].

Although very accurate, the calculations were performed on a supercomputer not available to an industrial environment. The emergence of similar models with a complete prediction of the keyhole creation with numerical methods like level set or V.O.F. depending on the software begin from 2008 with works of Geiger [3], Pang [4], Cho [5]. Since this time, few (less than 10) teams have proposed this kind of complete models but the experimental validation is often too brief, and the calculation times keep often these models in the academic research field.

The model developed and presented here aims to conciliate keyhole prediction and reasonable calculation times. For that, a method is proposed with 3 meshes and a voluntary simplified method to treat the vaporization and recoil pressure effect. Equations and numerical procedure are presented and next, to discuss on the simplified parts (vaporization and energy deposition) a complete experimental validation is proposed. Finally, the model is used in different configurations and compared to micrograph cuts.

MATHEMATICAL FORMULATION

In order to model all the required phenomena, a certain number of equations must be solved.

Classically, in all phases, the heat transfer equation is solved (Eqn. (1)). Two specific source terms are introduced: Slaser Eqn. (2) (discussed in [6]) introduces the laser energy and Qvap

(Eqn. (3)) is a term to remove the energy of the latent heat of phase change during vaporization. This energy is depending of the mass flow rate of evaporation (Eqn. (4)) and the saturation pressure (Eqn. (5)). Regarding the latent heat of fusion-solidification, the equivalent heat capacity is used (Eqn. (6)).

( ) ( )

eq

p laser vap

c T u T k T S Q



^^t +   ^=    + + ⁽¹⁾

( ) ( ) ( ) (

0

) (

² 0

)

²

2 2

0 0

2.5 cos ^laser exp

laser

x x y y

S P

r r

    



− − − − 

=  

 

  ⁽²⁾

vap v ( )

Q = −L m



 (3)

( )

( ) 1 2

sat

r b

p T m m

k T



=



− (4)

( ) exp ^v 1 ^vap

sat a

b vap

H T

p T p

k T T

   

=   − 

  

  ⁽⁵⁾

( )^²

( )²

exp

( )²

melt

l s

T T T T eq

p p s f

l s

c c L

T T



 − 

− 

 − 

= − +

−

(6)

Where T is the temperature, k the thermal conductivity,  the density, cp the heat capacity, u the velocity vector, θ the surface inclination, r0 the laser beam radius, m the mass flow rate of the evaporation,  the level set variable, kb the Boltzmann constant, psat

the saturation pressure of the vapor phase, m the atomic mass, r the retrodiffusion coefficient, pa the ambient pressure, Hv the creation enthalpy of the vapor, Tvap the vaporization temperature, Tl the liquidus temperature, Ts the solidus temperature, Tmelt the melting temperature average between Ts and Tl , Lf the latent heat of fusion.

In liquid (melt pool) and gas (keyhole), the momentum conservation equation is solved (Eqn. (7)) for incompressible Newtonian fluids. p is the pressure, gthe gravity acceleration,

l the expansion coefficient,  the surface tension coefficient, nthe normal vector,  the curvature, fl the liquid fraction and C and b some numerical coefficients. Some terms have been introduced:grepresents the gravity effect.  l l

(

T−Tmelt

)

g considers the buoyancy effect but only in the molten steel with the term . This effect is neglected in gas. Kuis called the Darcy condition. This term will cancel the velocity if the temperature is below the fusion point. This well-known numerical technique allows a simple treatment of the solid with the fluid flow calculation and a moving fusion/solidification front. Finally, the last term is necessary with the level set method to take into account the surface tension effect.

Eqn. (10) is the mass conservation equation for an incompressible fluid and Eqn. (11) is the transport equation of the level set variable. The coupling with fluid mechanic is operate with the u term. The term on the right of the equal sign is a numerical adjunction to facilitate the convergence. All theses equations are solved in the commercial code Comsol Multiphysics. Some of them are already implemented in the software (Eqn. (1,7,10,11)), the other have to be added (Eqn. (2-6,8,9)).

( ) ( ( ) )

( ) ^{( )}

T

l l melting

u u u pI u u

t

g T T g Ku n

 

    

 +    =   − +  +   

    

 

+ − −  + + 

(7)

( )

²

3

1 _l

l

K C f

f b

 − 

= −  

 + 

  ⁽⁸⁾

0

1

s l

l s l

l s

l

for T T T T

f for T T T

T T

for T T

 

 −

= −  

 



(9)

  =u 0 ⁽¹⁰⁾

(

¹

)

ls ls

t u



^



^

 +  =    −  −   

   

(11)

THERMOPHYSICAL PROPERTIES

The material used is a DP600 steel. This common material in automotive industry and is well known in the solid state [7]. When known, temperature dependent properties are used.

When not known, at high temperature for example, some assumptions are made. The property is either for pure iron or the last known value. The main properties are listed in table 1.

Table 1 Thermophysical properties used

Property [unit] Symbol Value or reference

Thermal conductivity of solid / liquid / gas [W.m^-1.K^-1] ks/l/g

f (T) [internal data] / 3.75.10^-3 T + 23.25 / 0.12

Density of solid / liquid / gas [kg.m^-3] ρs/l/g f (T) [internal data] / 7287 / 1 Specific heat of solid / liquid / gas [J.kg^-1.K^-1] cp s/l/g f (T) [internal data] / 573 / 373 Dynamic viscosity of solid / liquid / gas [Pa.s] μs/l/g 10³ / 5 .10^-3 / 1 .10^-5

Latent heat of fusion [J.kg^-1] Lf 2.54 .10⁵

Latent heat of vaporization [J.kg^-1] Lv 6.1 .10⁶

Liquidus temperature [K] Tl 1808

Solidus temperature [K] Ts 1788

Melting temperature [K] Tmelting (Tl+ Ts )/2

Vaporization temperature [K] Tvap 3134

Surface tension coefficient [N. m^-1] γ 1

Coefficients in Darcy’s law (7) C / b 1 .10⁶ / 1 .10^-3

Coefficient of thermal expansion [K^-1] bl 1 .10^-4

SIMPLIFICATION STRATEGY

All the previous equations are widely used by many authors [2-6]. On the other hand, the way to treat the metal vaporization differs from one author to another because of the numerical complexity. The most rigorous method requires a very thin modelling of the Knudsen layer at the interface. This modelling is potentially very difficult to achieve (discussed in [7]) and can be very time-consuming). For this reason, a simplified modelling with a simple recoil pressure is proposed by applying only the resulting pressure (Eqn. (13)) in momentum conservation equation (Eqn. (12)). This surface force (in Pa or N/m²) is applied on the surface of the keyhole with the ^{ }

( )

term (in m^-1), the spatial derivative of the  function. This term is nonzero only on the free surface of the keyhole. This method, used by many other authors, does not generate vapor plume and takes only into account its effect on the liquid. Although less accurate, this method can be relatively fast and give interesting prediction of keyhole as discussed in the next section.

𝜌 (𝜕𝑢⃗

𝜕𝑡 + 𝑢⃗ ⋅ (𝛻⃗ ⋅ 𝑢⃗ )) = 𝛻⃗ ⋅ [−𝑝𝐼 + 𝜇 (𝛻⃗ 𝑢⃗ + (𝛻⃗ ⋅ 𝑢⃗ )^𝑇)]

+𝜌𝑔 − 𝜌_𝑙𝛽_𝑙(𝑇 − 𝑇_{𝑚𝑒𝑙𝑡𝑖𝑛𝑔})𝑔 𝛷 + 𝐾𝑢⃗ + 𝛾𝑛⃗ 𝜅𝛿(𝛷) + 𝐹_{𝑟𝑒𝑐𝑜𝑖𝑙}

(12)

( )

exp ^v 1 ^vap

recoil a

b vap

H T

F p

k T T



  

=   −  

  

  ⁽¹³⁾

MESH USED

In this model, 3 main equations are solved: heat transfer Eqn. (1), fluid flow Eqn. (7) and level set transport equation Eqn. (11). Each equation does not have the same requirement in size of meshing. For example, the heat transfer problem needs a refined mesh in the laser deposition zone only. On the other hand, the fluid mechanic problem requires a finer mesh than heat transfer and with a refinement in the melt pool. Finally, the transport equation of level set needs a homogeneous mesh everywhere the interface (free surface) can potentially be.

The common way to solve this kind of model (with commercial codes) is by using a unique mesh often built for the more restrictive equation to solve. In this work, we suggest to use a specific and appropriate mesh for each main equation. Figure 1 presents the different meshes for each physics problem.

Fig. 1 View of the 3 mesh used. (a) Level Set, (b) Heat transfer, (c) fluid flow.

(a) (b) (c)

During the solving process presented in figure 2 with an iterative method, the coupling variables must be transferred from a mesh to another. For that, Comsol Multyphysics offers interpolation tools to project a solution on another mesh. If the calculation points do not coincide, an interpolation between the two closers points is realized. Of course, when a solution is computed on a fine mesh and then projected on a coarser mesh to operate the coupling, the solution is degraded. Therefore, the next section will analyze if the model remains predictive by comparing numerical results to experimental data.

The addition of the 3 meshes represents 400 000 degrees of freedom (DOF) to be compared to the 2 000 000 DOF of a classical formulation (more details in [6]). Time steps are set to 100 µs and the model is solved on 300 ms, which represents the physical time necessary to reach an established state. The solver used with the previously described meshes requires 7 GB of RAM. The model is solved on a common commercial workstation with a 2 cores, 4 threads i5-6500 at 3.20 Ghz and 16 GB RAM DDR4. The application of the simplified recoil pressure and the use of multi-meshes make it possible to increase the time step and so the calculation speed. In comparison of a model more complex and classic in the meshing method [6], the computation time is reduced from 3 weeks to 24 hours on the same computer.

Fig. 2 Iterative procedure principle allowing the use of three different meshes.

RESULTS AND DISCUSS

All the previous equations and methods are solved to describe first a fusion line inside a 1.8 mm steel plate. A 4 kW laser with a 600 µm diameter spot is used and the sheet moves with a speed of 6 m/min. As results, figure 3 presents the different stages of the creation of the keyhole and the melt pool growing. The vaporization temperature is reached in less than 2 ms and the keyhole begins to dig under the effect of the recoil pressure. When the keyhole is opened on the rear face, the melt pool finally stabilizes at 6 mm long in approximatively

200 ms. Figure 4 shows temperature and velocity fields. The temperature is around 3100 K inside the keyhole controlled by the evaporation latent heat. The velocity field is plotted only in the melt pool: indeed, although the fluid flow is calculated in it, the vapour plume velocity is not well described because of the choice of the recoil pressure simplification. To evaluate the interest of this model, a complete experimental validation must be carried out.

Fig. 3 Keyhole and melt pool creation at 2, 20, 60,100, 150 and 200 ms.

Fig. 4 Temperature field (max : 3150K) and velocity field (max : 0,6 m/s).

EXPERIMENTAL VALIDATION

To validate the model, a complete experimental set-up is proposed. First and classically, figure 5 compare the melted zone in a transversal cut after chemical etching. If the melted zone is well predicted for the case presented above (4 kW and 6 m/min), it is more interesting to note that a reduced velocity leads to a partial penetration on the experiment and in the model. The model predicts correctly the keyhole depth and can be used to predict

this kind of defects. To obtain an accurate keyhole shape, a thorough calculation of the energy deposition and the recoil pressure is essential.

Fig. 5 Comparison of melted zones. 4 kW – 6 m/min (left) 4 kW – 8 m/min (right)

Another indicator of the correct calculation of the energy deposition and the recoil pressure concerns the keyhole inclination angle. In figure 6, an experiment has been carried out with a sudden stop of the laser. With a small uncertainty thanks to the very thin liquid film beside the front wall of the keyhole, the longitudinal cut gives the inclination angle.

Again, the model gives a consistent result.

Fig. 6 Calculated keyhole inclination (left) compared to experiment (right).

In addition, temperatures in liquid and solid are compared to experiments. At the surface of the melt pool, a pyrometry method developed at PIMM laboratory is used (detailed in [6]). Different areas close and away from the keyhole are compared. At all the chosen points, the difference is lower than 5% showing a good prediction of the melt pool temperature.

Fig. 7 Comparison between calculated temperatures and measurements from pyrometry in different regions.

Time [ms]

In the solid part, at lower temperatures, the figure 8 shows for 11 positions a comparison between micro-thermocouples and model calculation. Again, the prediction of the temperature level is correct but with more discrepancies. The bigger difference is observed during the cooling. Indeed, to reduce calculation times, a small heat transfer domain has been modelled and so the solution is very sensitive on the boundary conditions. The size of the domain is here responsible of the disparity. We can treat this problem either by increasing the domain size (best way but with an increase of time calculation) or by modifying the boundary conditions (more difficult and have to be done for each experiment). Nevertheless, near the melt pool, where the thermal gradient is high, the model gives a good prediction validating also the fluid flow computing.

For this configuration, the model gives interesting results with a possible prediction of the melt-pool shape and the global temperatures reached. The model will be now used for configurations close to industrial needs to see the versatility.

Fig. 8 Comparison of temperatures calculated and measured by micro-thermocouples of 25µm.

(Left: below the sheet. Right: above the sheet. The indicated distance is the spatial shift from the welding axis)

VERSATILITY ON OTHER CONFIGURATIONS

Without any other modifications than the operating parameters (here sheets geometries or laser speed), the model is now compared to more complicated geometries but closer to industrial issues. First, figure 9 shows a result for a 3 mm thick plate. For the same laser power than previously, the welding speed is decreased. Naturally, the keyhole is modified with almost a 90° inclination and the melt pool is shorter. The micrograph cut can be seen in figure. 10.

[time] [time]

Fig. 9 3mm thickness – 4kW 3 m/min

Thanks to the level set method and the free surface method, it is easy to treat gap between sheets or a difference of thickness leading to the collapse of the thicker sheet on the other.

Figure 10 shows various configurations and the micrograph cut in comparison. We can note again that many phenomena are correctly modelled like partial penetration (4 kW 7 m/min), melt pool wider at the surfaces or melt pool globally wider with an initial gap between sheets.

For the results presented in figure 10, the more expansive calculation requires 48h which is close to industrial requirements for analysing defects, or new welding assembly configurations.

Fig. 10 Comparison between model and experiment (PIMM Lab.) for 6 welding configurations.

CONCLUSION

In this paper, a complete and predictive model is proposed to describe dynamically heat and fluid flow during laser welding. The keyhole is self-calculated. In order to propose a fast calculation running on classic computers, some approximations and improvement are proposed.

First, the vaporization – recoil pressure effect is simplified by a simple surface equivalent force. The very strong vapor plume is so neglected to facilitate the convergence. Then, a

new approach of meshing multiphysics models is proposed. The use of 3 meshes allows to refine only in the necessary locations and reduce drastically the number of DOF. This method can degrade the solution and so a complete experimental validation is performed.

The prediction of the thermal levels and the sizes of the melted zones is satisfying. All the calculations are performed on classical computers (4 cores, 16 GB ram) in less than 48h allowing the use of these models in an industrial context.

ACKNOWLEDGMENTS

The authors wish to thank the PIMM laboratory (Arts et Metiers – Paris France) for their support for the experimental validation of the model and ArcelorMittal for the funding of this study and the works associated.

REFERENCES

[1] M.DAL,R.FABBRO: ‘An overview of the state of art in laser welding simulation, Optics &

Laser Technology Vol 78, Part A, 2-14, 2016.

[2] H.KI,P.MOHANTY,J.MAZUMDER: ‘Modeling of laser keyhole welding: Part I. Mathematical modeling, numerical methodology, role of recoil pressure; multiple reflections and free surface evolution’, Metall. Mater. Trans A 33A 1817-1830, 2002.

[3] M.GEIGER,K.H.LEITZ,H.KOCH,A.OTTO: ‘A 3D transient model of keyhole and melt pool dynamics in laser beam welding applied to the joining of zinc coated sheets’, Prod. Eng. Res.

Devel. 3 127-136, 2009.

[4] S.PANG,L.CHEN,J.ZHOU,Y.YIN,T.CHEN: ‘A three-dimensional sharp interface model for self-consistent keyhole and weld pool dynamics in deep penetration laser welding’, J. Phys. D : Appl. Phys. 44 025301, 2011.

[5] W.I.CHO,S.J.NA,C.THOMY,F.VOLLERTSEN: ‘Numerical simulation of molten pool dynamics in high power disk laser welding’, J. Mater. Process. Technol. 212 262-275, 2012.

[6] M. COURTOIS, M. CARIN, P. LE MASSON, S. GAIED, M. BALABANE: Guidelines in the experimental validation of a 3D heat and fluid flow model of keyhole laser welding’. J. Phys : D Appl. Phy 49 155503 (13pp), 2016.

[7] M.COURTOIS,M.CARIN,P.LE MASSON,S.GAIED,M.BALABANE: ‘A new approach to compute multi-reflections of laser beam in a keyhole for heat transfer and fluid flow modeling in laser welding’. J. Phys : D Appl. Phy 46 505305 (14pp), 2013.