Numerical Investigation of Laser Surface Hardening of AISI 4340 Using 3D FEM Model for Thermal Analysis of Different Laser Scanning Patterns

Laser surface hardening is becoming one of the most successful heat treatment processes for improving wear and fatigue properties of steel parts. In this process, the heating system parameters and the material properties have important effects on the achieved hardened surface characteristics. The control of these variables using predictive modeling strategies leads to the desired surface properties without following the fastidious trial and error method. However, when the dimensions of the surface to be treated are larger than the cross section of the laser beam, various laser scanning patterns can be used. Due to their effects on the hardened surface properties, the attributes of the selected scanning patterns become significant variables in the process. This paper presents numerical and experimental investigations of four scanning patterns for laser surface hardening of AISI 4340 steel. The investigations are based on exhaustive modelling and simulation efforts carried out using a 3D finite element thermal analysis and structured experimental study according to Taguchi method. The temperature distribution and the hardness profile attributes are used to evaluate the effects of heating parameters and patterns design parameters on the hardened surface characteristics. This is very useful for integrating the scanning patterns’ features in an efficient predictive modeling approach. A structured experimental design combined to improved statistical analysis tools is used to assess the 3D model performance. The experiments are performed on a 3 kW Nd:Yag laser system. The modeling results exhibit a great agreement between the predicted and measured values for the hardened surface characteristics. The model evaluation reveals also its ability to provide not only accurate and robust predictions of the temperature distribution and the hardness profile as well an in-depth analysis of the effects of the process parameters. How to cite this paper: Tarchoun, B., El Ouafi, A. and Chebak, A. (2020) Numerical Investigation of Laser Surface Hardening of AISI 4340 Using 3D FEM Model for Thermal Analysis of Different Laser Scanning Patterns. Modeling and Numerical Simulation of Material Science, 10, 31-54. https://doi.org/10.4236/mnsms.2020.103003 Received: February 21, 2020 Accepted: July 14, 2020 Published: July 17, 2020 Copyright © 2020 by author(s) and Scientific Research Publishing Inc. This work is licensed under the Creative Commons Attribution International License (CC BY 4.0). http://creativecommons.org/licenses/by/4.0/ Open Access


Introduction
Laser surface transformation hardening (LSTH) is a heat treatment process applied to mechanical part surfaces in order to reach specific mechanical properties needed to achieve specific tasks [1]. LSTH increases the wear and the fatigue resistance by improving the properties of critical surface areas through the application of concentrated, short and localized thermal energy. Compared to other surface heat treatment technologies such as induction treatment and conventional furnace treatment, LSTH is a very advantageous technique integrating good precision, high productivity, good reproducibility and highly regarded quality. Generating much localized heating zone and a fast hardening speed limiting deformation and producing a fine microstructure, LSTH is an essential process in many applications [1] [2].
In some applications of LSTH where the objective is to heat-treat very large surfaces, the most commonly used technique is the multi-track processing that involves multiple laser beam passes to cover the entire surface to be heat-treated. In addition, it is recommended to use a laser beam with the widest possible focal spot and maintain a uniform distribution of thermal intensity to produce a constant hardness profile throughout the hardened surface. However, this option requires a powerful laser source to produce a very wide focal spot in order to keep the number of passes as well as the processing time to the minimum possible. Unfortunately, despite the advantages of the multi-track laser hardening process, the overlapped area between different laser tracks remains a major problem. Different zones in the heat-treated part after multi-track laser hardening process are illustrated in Figure 1. In the overlapped area located in the middle of two overlapped laser tracks, an over-tempered zone with a lower hardness and a tempered zone with hardness that may not reach the minimum hardness requirements are observed. To overcome this problem, is important to investigate other economical and efficient strategies than the traditional multi linear tracks processing approach [3] [4]. Figure 1. Illustration of the multi-track processing approach: (1) Base material, (2) Hardened zone, (3) Partial hardened zone, (4) Tempered-zone, (5) Over-tempered and partial re-hardened zone, and (6) Re-hardened zone.
Few studies have addressed this problem in a structured and comprehensive way. The proposed studies are based on the idea which intend to produce the most uniform hardness profile for specific areas instead of aiming to heat-treat the entire surface, either by matrix of points, matrix of lines or by more complex patterns such as rotating, sinusoidal, zigzags, etc. [3] [4] [5] [6] [7]. What is interesting with the proposed approaches is that the input energy as well as the processing time can be gradually optimized while ensuring cost minimization. Still, it is necessary to confirm that these approaches can ensure uniform hardness distribution in order to obtain mechanical performances similar to that achieved using other surface hardening techniques. Farrahi et al. and Zhou are the more recent researchers who have worked on complex patterns as alternative to multi-track approach for laser surface transformation hardening [3] [4] [5] [6] [7]. They focused on analytical modeling of laser surface hardening according to various scanning patterns and showed that the maximum temperature reached using rotating patterns was lower than that achieved using zigzag patterns. The relationship between the pattern frequencies and the temperature variation is also investigated [3] [4] [5] [6] [7].
The achieved results are very useful to understand some facts about the complex scanning patterns implementation. However, there is a lot of information missing. The effects of the geometrical attributes of the scanning patterns on the hardness profiles are to evaluate. The interaction between these attributes and the laser system heating parameters are also to assess and analyze. The control of process variables using predictive modeling strategies for producing the expected surface quality without being forced to use the traditional test and error method is an attractive approach to investigate. Numerical and experimental investigations, metallographic and microscopic analysis are essential to characterize hardness profile (HP), hardened depth (HD) and the hardened width (HW) produced using each particular scanning pattern.
This research is devoted to numerical and experimental investigations of laser surface transformation hardening of AISI 4340 steel using a 3D thermal analysis of different laser scanning patterns. The investigations are based on exhaustive modelling and simulation efforts carried out using a 3D finite element thermal analysis and structured experimental study according to Taguchi method. The temperature distribution and the hardness profile are used to evaluate the effects of heating parameters and laser scanning patterns parameters on the hardened surface characteristics. The performance of the 3D FEM model is assessed using structured experimental data produced using a 3 kW Nd:Yag laser system.

Thermal Conduction
Several researchers developed numerical models for predicting the temperature distribution during the laser heat treatment process [8] [9] [10] [11] [12]. Billaud et al. established simple models to predict the hardness profile in the case of mechanical parts treatment based on 3D numerical simulations models using the finite element method [8]. Caron-Guillemette et al. focused their research on modelling and simulating of heat flows resulting from the interaction of a laser with steel parts, while developing a numerical function that automatically and accurately models a laser beam traversing complex surfaces, such as the teeth of a gearwheel [9]. Several other researchers approached the 3D numerical modeling based practically on the same formulations. The differences are related to the type of laser source, the nature of displacement path, and the geometry of the part to be treated. The Fourier-Kirchhoff heat flux equation is used to evaluate the temperature distribution according to time and specific spatial coordinates.
Based on the temperature distribution and using the metallurgical transformation equations, the hardness profile and eventually the hardened depth can be easily estimated. Using the same mathematical formulation, Ashby et al. described the metallurgical aspect of hypo eutectoid laser surface transformation of steel hypo eutectoid [12]. The heat transfer in a solid can then be described by Equation (1): where ρ, C p and λ are respectively the density, the specific heat and the thermal conductivity of the material.
where W represents the Gaussian beam radius, P the laser power, R c the reflection coefficient of the surface of the material to be treated, A C the absorption coefficient of the laser beam power, SS the scanning speed, L z the thickness of the part and finally x 0 and y 0 the laser beam center coordinates at the time t = 0 s. Therefore, E represents a Gaussian heat source that moves along the X-axis at a velocity SS. The emissivity of the surfaces as well as the thermal convection by the surfaces to the ambient air are expressed respectively by the Equations (3) and (4) where ε represents the emissivity of the air, σ the Boltzmann constant, h the thermal transfer coefficient and n the complex refractive index. T represents the temperature distribution, T amb is the ambient temperature and T ext is the external temperature [12].
( ) ( ) According to Komanduri et al., the increase in total temperature at any point M in the surface and at any time t is given by Equation (5) [13]. Where E represents the heat source expressed using Equation (2), SS the scanning speed, x 0 the radius of the moving focal spot, x i the radius of a segmental focal spot, a the thermal diffusivity of the medium (cm 2 /s), λ the thermal conductivity of the material, I 0 the modified Bessel function of first kind order zero, ω a dimensionless number varying from 0 to SS 2 t/4a, θ the tilt angle and X, Y, Z represent the coordinates of any point M in a moving coordinate system where the temperature rise is concerned [13]. This equation is very useful to evaluate the surface temperature distribution of each trajectory performed during the laser hardening process.

Metallurgical Transformation
In equilibrium, austenite formation generally begins at the temperature delimited by the Ac 1 curve found in the iron-carbon phase diagram. From this point, the carbon contained in the perlite grains becomes in solid solution, homogenizes in the perlite before migrating to the ferrite grains, which contain a small percentage of carbon. The diffusivity of carbon is approximately 1 × 10 −5 e −9.0/T m 2 /s in austenite and 6 × 10 −5 e −5.3/T m 2 /s in ferrite [14] [15]. In rapid cooling, only austenite regions containing sufficient carbon will be transformed into martensite. Generally, it is well known that the heating rate greatly affects the formation of austenite and its homogenization. Therefore, it is necessary to allow sufficient time for the carbon to migrate homogeneously into the austenite, otherwise the transformation will not be uniform and may contain mixed zones in the hardened area.
The time required for carbon diffusion is related to the distance it must travel, i.e. the initial microstructure of the steel. The formation of martensite requires that the cooling of the austenite must respect a well-defined rate. Unlike austenite, which absorbs carbon from neighboring regions, martensite is produced by atoms displacements to form crystal structure that traps carbon easily (intersti-tial sites). This crystalline arrangement with a lower density and higher volume than austenite is at the source of the geometric deformations produced by the transformation hardening. The increase in volume with about 4% according to Steen et al. [14], produce some desired residual compressive stresses for certain mechanical loading situations [16]. During the heating cycle, temperature peaks affect the level of the structural change in the material. This phenomenon is expressed by the equation of kinetic force (6).
where Q represents the activation energy of the structural transformation, R the gas constant, T p the temperature peak at a considered temperature, α the thermal diffusivity and τ the time constant [12]. The terms α and τ are determined by Equations (7) and (8).
where T 0 is the initial temperature, R c is the reflection coefficient and λ is the thermal conductivity of the material. The spacing of the perlite layers has been described by Ashby and Easterling using Equation (9).
In heating cycles where the temperature is time dependent, the spacing will be sufficient for the total transformation of the perlite into austenite [12] [16]. D 0 represents the constant of carbon diffusion in ferrite. Indeed, the austenite is produced at the same carbon content as perlite (C e = 0.8%). At temperature Ac 1 , the austenite volume fraction becomes the volume fraction f i previously occupied by the perlite colonies. This represents the minimum volume fraction of martensite.
where C represents the carbon content in steel and C f the carbon content in ferrite which is negligible in this case. According to the time diagram (TTT), the maximum allowable martensitic fraction is expressed in Equation (11).
Furthermore, the fraction of the martensite is then given by Equation (12) with C c is the critical carbon content, D 0 is the pre-exponential for diffusion of carbon and g is the average of the austenite grain size.
Finally, the hardness can be determined by Equation (13)

Laser Beam Displacement
To perform each desired scanning pattern, the laser beam displacement is controlled by specific parametric equations reflecting the desired laser beam path. Figure 2 presents an illustration of the selected scanning patterns. The patterns parametric equations are presented in Table 1. In these equations, a mp represents the amplitude, d the distance between the circle center and the point of the trochoid pattern, SS the scanning speed, t the time and f r the frequency.

Simulation Parameters
Since the material properties change significantly with temperature, it is recognized that the numerical models describe the material transformation process better than analytical models. Generally, analytical models give a preliminary approximation of the phenomenon, then the numerical approach come often to confirm the results. Based on material properties and the previously represented equations, extensive simulations are performed using various scanning patterns applied on 4340 steel prismatic parts (60 mm × 6 mm × 5 mm). The thickness of the part to be treated must be large enough in order to not affect the temperature profile and to provide cooling cycle while allowing metallurgical transformation [5]. Table 2 presents the chemical composition of AISI 4340 steel and Table 3 illustrates different physicochemical properties [17].
In this study, the numerical investigation is based on a 3D finite element thermal analysis. This modelling and simulation approach provides an understanding of the behaviour of material properties as a function of processing temperature. To simplify the model and reduce the calculation time, the material is assumed homogeneous and isotropic, the physicochemical properties are defined at thermodynamic equilibrium, the ambient temperature is set at 20˚C (about 293 K) and the maximum temperature must not exceed 1400˚C. In addition, the simulations take into account the emissivity of the surfaces as well as the phenomenon of thermal convection by the surfaces to the ambient air.  The main parameters used in this numerical investigations are laser power (P), laser source displacement speed (SS), frequency (f r ) and amplitude (a mp ). For each scanning pattern, the temperature distribution is analyzed and its effects on the hardened surface characteristics are evaluated. Table 4 presents the specific simulation parameters for the selected scanning patterns. In this table, X 0 and Y 0 are the coordinates of the start point, x 0 and y 0 are the laser beam radius, T step is the time step for discretizing of model, T end is the time end of the simulation and b is the correction coefficient used from the model calibration. Figure 3 shows the evolution of the thermal conductivity (λ) and the specific heat (C p ) as a function of the temperature evolution. These properties are expressed by Equations (14) and (15) respectively. The variables constituting the parametric equation of the heat source are defined in Table 4 and expressed using Equation (16). The coefficient b is a model correction coefficient obtained by multiple linear regression as function of laser power and scanning speed.     T  T  T  T  C T  T  T  T

Mesh Stability Study
The convergence study is necessary in all simulation problems because it increases the results consistency, ensures the model validity, and saves a lot of computing time. The convergence study consists of creating an initial mesh using a reasonable mesh size, computing and analysing the first results, then decreasing the mesh size step by step and evaluating the current results by comparing them to previous iteration results. The points used for the evaluation of the surface temperature are uniformly and symmetrically distributed on both sides of the laser beam path. Three points (P 1 , P 2 and P 3 ) are fixed at x = 30 mm and z = 5 mm, and y = 2, 3 and 4 mm. Respectively at its points, three temperatures T 1 , T 2 , T 3 are measured for different mesh sizes. Given the symmetrical aspect of the distribution of points, the temperature T 2 must be higher than T 1 and T 3 since it is a centered position, while T 1 and T 3 must be equal. After several simulations, the temperatures become stable between the 0.016 and 0.32 mm mesh size. Outside this range, the temperature fluctuations are due to the accuracy and truncation errors. The selected mesh size of 0.016 mm is a reasonable compromise to achieve realistic precision without requiring excessive computing resources. This choice is also effective and well adapted to the three scanning patterns to be simulated. However, to improve the results accuracy smaller mesh size and higher mesh density are used in the part central areas. Figure 4 illustrates the selected final mesh size and Figure 5 shows the effect of mesh size on the temperatures variation.

Experimental Validation
A series of validation tests are performed according to the adopted experimental design to evaluate the hardness profile obtained for each scanning pattern. The experiments are done on AISI 4340 steel plates (60 mm × 6 mm × 5 mm) with an initial hardness of 25 to 30 HRC. The surface hardening is performed using an IPG YLS-3000-ST2 fiber laser with a maximum power of 3 kW powered by an ND:Yag laser source mounted on a six-axis Fanuc M-710IC robot. The used experimental setup is shown in Figure 6. For each scanning pattern, the appropriate laser beam displacement is ensured by specific path program implemented in the robot control.   speed). Table 5 presents the approved experimental design for each scanning pattern [18].   Figure 7. Typical hardened zones using a power of 800 W and a scanning speed of 20 mm/s for linear, sinusoidal, triangular and trochoid scanning patterns respectively. Figure 8 and Figure 9 represent the hardness curves for typical tests representing the variation effects of laser power and scanning. These hardness curves can be divided into three regions. The first region is characterized by high hardness compared to the part core. This region is totally austenitized and transformed to martensite upon rapid cooling. The second region presents a radical decreasing in hardness before reaching the core hardness. This region is characterized by the hardness loss caused by sharp drop in hardness to reach initial hardness. This zone is composed of a mixture of hard and over-tempered martensite since the temperature is between Ac 1 and Ac 3 . Finally, the third region corresponds to the zone not affected by the thermal flow. The hardened depth is then characterized by the first zone, called the hard zone. The resulting hardness values are relatively similar for all the experimental tests with about 60 HRC. The achieved results for the hardened depth and the hardened width for different scanning patterns are characterized by a good concordance between numerical simulations and experimental results. Table 6 illustrates the numerical and experimental results of the hardened depth and the hardened width for the various performed tests. Figure 10 and Figure 11 show a comparison between experimental and numerical results for HD and HW respectively. These results reveal globally that the developed numerical models can estimate HD and HW accurately with maximum errors not exceeding 9.1%. The model accuracy in confirmed by coefficients of determination (R 2 ) between 95% and 98%. The maximum error is observed in the case of linear and sinusoidal scanning patterns. The trochoid pattern presents the best results with 3.7%, 2.3% and 1.00% as maximum, average and minimum errors respectively.
Moreover, ANOVA is applied to the simulation data and the obtained results are compared to those found by ANOVA applied to the experimental data. The comparison of the parameters contribution between experimentation and simulation data for all considered scanning patterns is presented in Table 7. It can be observed that laser power and scanning speed effects on the variation of simulated HD and HW, in term of percentage contributions (% C), are practically the same as for the experimental data excepting for HW in the case of the triangular scanning pattern [18].       The largest difference between experimental and simulated results is about 12% observed for the power contribution in HW variation when using triangular scanning pattern. This is due to some HW measurement errors and probably to the fact of not considering the effect of the laser beam diameter variation in the simulations. The evaluation of the laser power and scanning speed interaction shows also non-significant effects as established in the experimentation investigations [18]. The contributions of various interactions do not exceed 5%. These results are corroborated by the graphs of effects presented in Figures 12-15. An F-test is also executed to confirm the laser parameter contributions and effects significance. Therefore, the proposed numerical model is sufficiently accurate for been used in this numerical investigation and ultimately to generate additional data that can be very useful to integrate in an eventual predictive modelling approach.

Simulation Results and Discussions
Numerous simulations are performed under various heating parameters and conditions using different scanning patterns. The main objective of these exhaustive simulation efforts is to collect sufficient data to evaluate the effects of heating parameters and patterns design parameters on the hardened surface This will significantly help to identify the most appropriate scanning patterns in conjunction with heating parameters to reach a uniform temperature distribution and maximum treated surface dimensions. Additional simulations using a linear trajectory are performed to evaluate the results accomplished by different scanning patterns compared to the traditional single tract laser hardening.    During the simulations, the temperature appears, as anticipated, always at its maximum under the laser beam spot and decreases progressively as the spot moves away. Outside the focal spot, the part absorbs the heat by conduction and the temperature drops gradually. In addition, some of the heat is lost by convection in the ambient air and by radiation through the room. However, due to its good thermal absorption coefficient, the material absorbs most of the heat emitted by the laser. The part volume plays a very important role with a direct cooling of the heated areas after reaching the austenization temperature without the need for forced convection [13]. These observations are illustrated in Figure 16. It was also observed that under the same heating parameters, the temperature distribution is related to the used scanning patterns. Consequently, the dimension of the treated zone depends on heating parameters and scanning pattern design parameters. Figure 17 shows the maximum temperature achieved by different scanning patterns as a function of y-axis observed in the section located at x = (5; 30 and 55 mm). These temperature curves are obtained using the same heating parameters in order to compare the effect of the scanning patterns on the temperature distribution. This temperature distribution visualization looks like a cross section that represents the width of the treated zone and the reached maximum temperature represents practically the hardened depth. According to this figure, it appears that the sinusoidal pattern produces the largest treated width. The highest maximum temperature is achieved using the linear pattern that concentrates the heat flow in restricted area in the middle of the part. Linear pattern generates a temperature of about 1500˚C exceeding the melting temperature of 4340 steel. At y = 2 and 4 and along x-axis, the sinusoidal path presents the highest temperature magnitude (about 1100˚C) but the lowest temperature in the part center (y = 3) with (about 900˚C). The temperatures reached using triangular and trochoid patterns are equivalent with some variations of the order of 30˚C to 100˚C. The temperature distribution is practically symmetrical. Figure 18 shows the temperature distribution over the entire surface. This helps to not only detect the variation of the heat flow distribution, but also to compare the temperature distribution uniformity for different scanning patterns and to visualize the maximum temperature reached at each point of the surface. Consequently, the dimensions of the treated zone can be estimated accurately. These figures reveal that the linear pattern is characterized by the most uniform distribution longitudinally but heat treat only a limited portion of the surface. In the case of linear multi-tracks pattern, the distribution can be improved but remains irregular. The sinusoidal pattern produces a relatively uniform temperature distribution longitudinally similar to that generated by the triangular pattern. The trochoid pattern produce a temperature distribution characterized by some periodicity in temperature peaks due to the intersection points in the laser beam path. Some points of the surface are heated two times.
The sinusoidal pattern produces a good distribution with a relatively high temperature that covers the majority of the treated surface offering a maximum width compared to other patterns under the same hardening conditions (P = 800 W, SS = 20 mm/s). For the triangular pattern, the heat flux is more concentrated in the part centerline and the temperature remains relatively low at the margins (y = 2 and 4 mm). The trochoid pattern is characterized by temperature peaks generated by the double passage of the laser beam over some points on the treated surface, leading to a non-uniform and non-symmetrical distribution of the heat. Generally, these temperature peaks are undesirable since they can usually lead to melting zone and surface deformation.   The identification of the temperature distribution on the surface indicates the delimitation of the hardened zone and consequently allows to estimate the hardened depth and the hardened width for each scanning pattern. Figure 19 shows contour plots representing the temperature distributions according to the part surface points coordinates and as a function of the part cross-section points coordinates. These plots provide a good indication how the temperature is distributed not only on the surface but also on different depths. This allows furthermore to estimate the volume of transformed material using the same energy for different scanning patterns. Consequently, HD and HW can be estimated easily. Figure 19 confirm the results observed from Figure 18. Although sinusoidal pattern does not produce a transformed material volume as high as that produced using the trochoid pattern (105 vs. 108 mm 3 ), it presents the best compromise in terms of HD and HW. Indeed, sinusoidal pattern produces a relatively uniform temperature distribution with a reasonably high temperature covering the entire surface and offering the maximum width compared to other patterns under the same processing conditions. For specific applications, maximizing the surface of the hardened zone requires the selection of the appropriate scanning pattern. An adequate scanning pattern is the one that produce a maximum uniform hardened width while respecting a minimum limit of hardened depth. In such a case, if a minimum is set for HD, the HW/HD ratio can be used as an indicator to evaluate the hardened zone surface according to the laser parameters for evaluating different scanning patterns. HW/HD ratio can provide also useful indication on how the generated heat is absorbed by the surface. Table 8 presents the HW/HD ratio for various simulated combination of laser parameters and scanning patterns. It can be observed that the general average HW/HD ratio is higher than that calculated ex-cluding tests where HD is lower than 200 μm. Sinusoidal and triangular scanning patterns present the highest HW/HD ratio allowing to achieve the largest hardened surfaces. With high HD and small HW, linear pattern produces a low and relatively constant HW/HD ratio. Trochoid pattern represents an improved version of linear pattern with approximatively the same HD but with a HW more than 2.5 times higher. Then, if the laser parameters can be set to control a constant minimal HD to 200 μm, by adopting the maximum ratio, HW can reach a minimum more than 3500 μm for sinusoidal and triangular scanning pattern. If the total generated heating energy is adjusted to maintain a constant minimal HD, HW can reach higher values.

Conflicts of Interest
The authors declare no conflicts of interest regarding the publication of this paper.