Residual Stress Modeling and Deformation Measurement in Laser Metal Deposition Process

Direct metal deposition (DMD) has become very popular within the space of rapid manufacturing and repair. Its capability of producing fully dense metal parts with complex internal geometries, which could not be easily achieved by traditional manufacturing approaches, has been well demonstrated. However, the DMD process usually comes with high thermal gradients and high heating and cooling rates, leading to residual stresses and the associated deformation, which can have negative effect on product integrity. This paper studies the features of thermal stress and deformation involved in the DMD process by constructing a 3-D, sequentially coupled, thermomechanical, finite element model to predict both the thermal and mechanical behaviors of the DMD process of Stainless Steel 304 (SS 304). A set of experiments were then conducted to validate deformation using a laser displacement sensor. Comparisons between the simulated and experimental results show good agreement. This model can be used to predict the mechanical behavior of products fabricated by the DMD process and to help with the optimization of design and manufacturing parameters.


Laser-aided direct metal deposition
Laser-aided direct metal deposition (DMD) is an advanced additive manufacturing (AM) technology which can produce fully dense, functional metal parts directly from CAD model. In its operation, the laser beam is focused onto a metallic substrate to create a melt pool, and a powder stream is continuously transported into the melt pool by the powder delivery system. The substrate is attached to a computer numerical control (CNC) multi-axis system, and by moving the substrate according to the desired route pattern, a 2-D layer can be deposited. By building successive layers on top of one another (layer by layer), a 3-D object can be formed. The DMD process has demonstrated its ability in the area of rapid manufacture, repair, and modification of metallic components. Practically, this process is most suitable for components with complex internal geometries that cannot be fabricated by traditional manufacturing methods such as casting.
Furthermore, this process is very cost-effective compared with traditional subtractive manufacturing techniques because it can produce near net shape parts with little or no machining.

Residual stress and distortion
Residual stresses exist in a part while no external loads were applied. When a part is heated evenly from its previous residual stress-free state, it expands evenly, and no thermal stress is generated. However, when a part is heated unevenly, thermal stress is produced [1].
During DMD processing, highly localized heating and cooling lead to nonuniform thermal expansion and contraction, which further results in a complicated distribution of residual stresses in the heat-affected zone (HAZ) and distortion across the entire structure. These residual stresses could be detrimental-they may cause fractures, promote fatigue, and induce unpredictable buckling during the service of deposited parts; the deformation often is harmful to the dimensional accuracies of structures. Therefore, it becomes critical to predict the two behaviors of materials after the DMD process and to optimize the manufacturing parameters to reduce the residual stresses and deformation.

Literature review
The thermal behavior of the DMD process has been investigated numerically by many scholars. In [2], a 2-D finite element model was created to simulate the temperature field during the laser cladding process. The results indicated that quasisteady thermal field could not be reached in a short time. Other scholars have chosen to experimentally investigate thermal behavior. In [3], Griffith et al. employed radiation pyrometers and thermocouples to monitor the thermal signature during laser engineered net shaping (LENS) processing. The results showed that the integrated temperature reheat had a significant effect on the microstructural evolution during the fabrication of hollow H13 tool steel parts. Utilizing a twowavelength imaging pyrometer, Wang and Felicelli [4] measured the temperature distribution in the melt pool and the area surrounding it during the LENS deposition process. It was found that the maximum temperature in the molten pool is approximately 1600. Only thermal behaviors were investigated in these papers, while no residual stresses were modeled and analyzed.
Some researchers have focused on the modeling and simulation of traditional welding processes, which share many similarities with DMD processes. Using a double-ellipsoid heat source, Gery et al. [5] generated the transient temperature distributions of the welded plates. The results demonstrated that the welding speed, energy input, and heat source distributions had important effects on the shape and boundaries of heat-affected zone. Deng [6] investigated the effects of solid-state phase transformation on the residual stress and distortion caused by welding in low carbon and medium steels. The model discovered that the solid-state phase transformation did not have a noteworthy impact on the final residual stress and the welding deformation in low carbon steel. However, the final residual stresses and the welding deformation appear to be significantly affected by the martensitic transformation in medium carbon steels. Feli et al. [7] analyzed the temperature history and the residual stress field in multi-pass, butt-welded, stainless steel pipes. It was found that in the weld zone and its vicinity, a tensile axial residual stress is produced on the inside surface and compressive axial stress at the outside surface.
Other researchers have attempted to obtain the distribution of residual stress caused by the DMD process through experiments. For example, Moat et al. [8] measured strain in three directions using a neutron diffraction beam line to calculate the stress in DMD manufactured Waspaloy blocks. They found that large tensile residual stresses exist in the longitudinal direction near the top of the structure. Zheng et al. [9] measured residual stress in PZT thin films fabricated by a pulsed laser using X-ray diffraction. Although experiments can provide relatively accurate results, their flexibility and high cost limit their ability to serve as a general method by which to solve residual stress problems.
In recent years, analyses of the residual stress involved in laser deposition processes using the FE model have been well documented in the literary. In [10], 2-D FE model was built to observe the impact of process parameters on the melt pool size, growth-direction residual stress, and material properties in laser-based deposition processes. They observed that after the deposition was completed and the wall was cooled to room temperature, large tensile stresses exist in the vertical direction at vertical free edges, which is a contrast to the observations in this study. Wang et al. [11] utilized commercial welding software SYSWELD to characterize the residual stress in LENS-deposited AISI 410 stainless steel thin-wall plates. Tensile longitudinal stresses were found near the mid-height, and compressive stresses were found near the top and bottom of the walls. In [12], Kamara et al. investigated the residual stress characteristics of a laser-deposited, multiple-layer wall of Waspaloy on an Inconel 718 substrate. The results indicated that along the length of the wall, residual stresses were almost zero at the bottom and top of the wall. Along the height of the wall, tensile stress with large magnitudes existed at the bottom of the wall while close to the top surface, near stress-free condition, seem to prevail. This matches well with the results presented in this thesis.

Simulation and experiment approach
Based on the finite element (FE) analysis package ABAQUS, a 3-D, sequentially coupled, thermomechanical model was developed to simulate the transient temperature field, residual stress, and final deformation involved in the DMD process of Stainless Steel 304 (SS 304). The numerical modeling consisted of two major steps shown in Figure 1. A transient thermal analysis was firstly carried out to produce the temperature history of the entire workpiece. Then, in the second step, using the temperature field file generated in the previous step as load, a mechanical analysis was carried out to calculate the residual stress and deformation of workpiece.
A laser displacement sensor was used in the experiment to record the vertical deflection of the workpiece resulted from thermal stresses during the deposition process. The accuracy of the numerical model was checked and validated by comparing the experimental results with the simulation results. This validated model can be applied to a multilayer DMD process of stainless steel under different process parameters and can be used for other materials.

Thermal analysis 2.1 Governing equations
During DMD processing, the stress and deformation field in the workpiece largely depend on the temperature history; however, the effect of the stress and deformation field on the temperature field is insignificant. Thus, a heat transfer analysis not coupled with mechanical effect is considered.
The transient temperature domain T x, y, z, t ð Þwas attained by solving the heat conduction equation, Eq. (1), in the substrate, along with the initial and boundary conditions discussed in Section 2.2 [13]: where T is the temperature, ρ is the density, C is the specific heat, k is the heat conductivity, and Q is the internal heat generation per unit volume. All material properties were considered temperature-dependent.

Initial and boundary conditions
To solve Eq. (1), the initial conditions were set as where T 0 is the ambient temperature. In this study, the room temperature of298 K was used. The boundary conditions, including thermal convection and radiation, are governed by Newton's law of cooling and the Stefan-Boltzmann law, respectively. The heat source parameter, Q in Eq. (1), was considered in the boundary conditions as a surface heat source (moving laser beam). The boundary conditions then could be expressed as [13] where k, T, T 0 , and Q bear their previous definitions; n is the normal vector of the surface; h c is the heat convection coefficient; ε is the emissivity which is 0.9; σ is the Stefan-Boltzmann constant which is 5:6704 Â 10 À8 W=m 2 K 4 ; Γ denotes the surfaces of the workpiece; and Λ denotes the surface area covered by the laser beam.

Assumptions and adjustments
Accurate modeling of the thermal process yields highly nonlinear coupled equations, which is time-consuming and expensive to solve. To speed up the solution process and reduce the computational time without sacrificing accuracy, the following assumptions and adjustments were considered.

Energy distribution of the deposition process
In this study, a circular-shaped laser beam with a constant and uniform power density was used. Thus, to match the experiment setup, the heat source parameter Q in Eq. (1) was considered a constant and uniformly distributed surface heat flux defined as where α is the absorption coefficient, P is the power of the continuous laser, and r is the radius of the laser beam. α was set as 0.4 according to the previous experiments conducted, and r ¼ 1:25 mm.

Movement of laser beam
The motion of the laser beam was taken into account by updating the position of the beam's center R with time t as follows where x, y, and z are the spatial coordinates and the laser beam centers, u, v, and w, are the continuous velocities the laser beam travels along x-, y-, and z-direction.
In ABAQUS, a user subroutine "DFLUX" [14] was written to simulate the motion of the laser beam.

Powder projection
When modeling, the continuous powder injection process is broken into many small discrete time steps. Using the model change method provided by ABAQUS [14], in each time step, a set of finite elements was added onto the substrate to form deposits along the center line of the substrate. The width of the deposits was assumed to be the same as the diameter of the laser beam, and the thickness of the deposits was calculated from the laser or table travel speed and the powder feed rate. An efficiency of 0:3 was assumed for the power feeding process to account for the powder that did not reach the melt pool.

Modeling the latent heat of fusion
To account for the effect of the latent heat of fusion during the melting and solidification process, the specific heat capacity is modified to generate an equivalent specific heat capacity c * p as [15] where C * p T ð Þ is the modified specific heat, C p T ð Þ is the original temperaturedependent-specific heat, L is the latent heat of fusion, T m is the melting temperature, and T 0 is the ambient temperature. The values of the latent heat of the fusion, solidus temperature, and liquidus temperature of SS 304 [16] appear in Table 1.

Marangoni effect
As discussed in [17], the temperature distribution is significantly impacted by the effect of Marangoni flow, which is caused by the thermocapillary phenomenon. To obtain an accurate thermal field solution, based on the method proposed by [18], an artificial thermal conductivity was used to account for the Marangoni effect: where k m T ð Þ is the modified thermal conductivity, T liq is the liquidus temperature, and T and k T ð Þ maintain their previous definitions.

Combined boundary conditions
The boundary conditions shown in Eq. (4) can be rewritten as where h r is the radiation coefficient expressed as Eq. (9) shows that when the temperature is low, convection is dominant in heat loss and when temperature is high, radiation becomes dominant. As shown in Eq. (10), radiation coefficient is the third-order function of temperature T, which is highly nonlinear. This would greatly increase the computational expense and time. Based on experimental data, an empirical formula combining convective and radiative heat transfer were given by [19] as where h is the combined heat transfer coefficient which is a lower-order function of temperature T compared with h r . The associated loss in accuracy using this relationship is estimated to be less than 5% [20]. In ABAQUS, a user subroutine "FILM" is written to simulate heat loss.

Dimension and parameter
As shown in Figure 2, a finite element model for a one-pass, three-layer DMD process was built. The dimension of the substrate under consideration is 50:8 Â 12:7 Â 3:175 mm (2 Â 0:5 Â 0:125 inch). Two cases were simulated with different process parameters including laser power, laser travel speed, and powder feed rate. These parameters were chosen according to the criterion that the final geometry of deposits and the total energy absorbed by the specimen be the same in each case. These process parameters are detailed in Table 2.

Material properties
Temperature-dependent thermal physical properties of SS 304, including the density, specific heat, thermal conductivity, and latent heat, were used as inputs. The values of these properties appear in Figure 3.

Element selection method
Based on the computational accuracy and cost, the type and size of the finite elements used to simulate the substrate and deposits were determined. In transient heat transfer analysis, second-order elements generally produce more accurate results; however, there is a minimum time increment. A simple guideline can be written as [14] Δt> 6c ρk Δl 2 (12) where ρ, c, and k have been defined; Δl is the element dimension; and Δt is the time increment. If the time increment Δt is smaller than this value, nonphysical   oscillations may appear in the solution. According to [14], the first-order elements can eliminate such oscillations but may lead to inaccurate solutions [13]. To take both the computational efficiency and accuracy into account, the first-order heat transfer elements, C3D8, with h-version mesh refinement, to refine the mesh by subdividing existing elements into more elements of the same order, were used for the whole domain. Fine meshes were used in zones close to the deposits, and the mesh size progressively increased with the distance from the deposits. In regions more away from the heat-affected zone, coarser meshes were implemented. A total of 14,496 elements and 17,509 nodes was created.

Increment control
In order to obtain reliable results from the mechanical analysis, the maximum nodal temperature change in each increment was set as 5 K, and the time increments were selected automatically by ABAQUS to ensure that this value was not exceeded at any node during any increment of the analysis [14].

Governing equations
The total strain ε ij can be represented generally as where ε M ij is the strain contributed by the mechanical forces and ε T ij is the strain from thermal loads. Eq. (12) can be decomposed further into five components as [6] where ε E ij is the elastic strain, ε P ij is the plastic strain, ε T ij is the thermal strain, ε ΔV ij is the strain due to the volumetric change in the phase transformation, and ε Trp ij is the strain caused by transformation plasticity. Solid-state phase transformation does not exist in stainless steel [21], so ε ΔV ij and ε Trp ij vanish. The total strain vector is then represented as The elastic stress-strain relationship is governed by isotropic Hooke's law as where D ijkl is the elastic stiffness tensor calculated from Young's modulus E and Poisson's ratio ν as [11] where δ ij is the Kronecker delta function defined as For isotropic elastic solids, Eq. (15) can be simplified as Thermal strain ε T ij can be calculated from the thermal expansion constitutive equation where α is the thermal expansion coefficient and ΔT is the temperature difference between two different material points. Rate-independent plasticity with the von Mises yield criterion and linear kinematic hardening rule [21] was utilized to model the plastic strain.
Unlike the elastic and thermal strain, no unique relationship exists between the total plastic strain and stress; when a material is subjected to a certain stress state, there exist many possible strain states. So strain increments, instead of the total accumulated strain, were considered when examining the strain-stress relationships. The total strain then was obtained by integrating the strain increments over time t. The plastic strain-stress relationship for isotropic material is governed by the Prandtl-Reuss equation [22]: where dε P ij is the plastic strain increment, λ is the plastic multiplier, and s ij is the deviatoric stress tensor defined by By substituting Eq. (18), Eq. (19), Eq. (20), and Eq. (21) into Eq. (14) and taking the derivative with respect to time, the total strain rate can be described by [

Initial and boundary conditions
The temperature history of all the nodes generated in the thermal analysis was imported as a predefined field into the mechanical analysis. The only boundary condition applied to the domain was that the substrate was fixed on one side to prevent rigid body motion. In ABAQUS, the node displacements on the left side of the substrate were set as 0.

Material properties
Temperature-dependent mechanical properties including the thermal expansion coefficient [23], Young's modulus, Poisson's ratio [21], and yield stress [16] were used to model the thermomechanical behavior of SS 304. The values of these properties appear in Figure 4.

Element selection
The order of element and integration method used in the mechanical analysis differed from those used in the thermal analysis, while the element dimension and meshing scheme remained unchanged. To ensure the computational accuracy of the residual stress and deformation, second-order elements were utilized in the heataffected zone, while first-order elements were used in other regions to reduce the computation time. Preventing shear and volumetric locking [14] requires the selection of reduced-integration elements. Therefore, elements "C3D20R" and "C3D8R" in ABAQUS were combined in use to represent the domain.
As shown in Figure 5, the 3-D 20-node element used in the mechanical analysis had 12 more nodes than the 3-D 8-node element used in the thermal analysis. Therefore, when mapping the temperature data from the thermal analysis to the mechanical analysis, interpolation had to be conducted to obtain the temperature of the 12 extra mid-side nodes (Nodes 9-20 in Figure 5(b)). Figure 6 shows the temperature field of the melt pool and surrounding areas from top view at different times in Case 1 (laser power 607 W, laser travel speed 250 mm/min, powder feed rate 6.3 g/min). Laser beam cyclically moves along +z and -z-direction. At t ¼ 0:9 s, t ¼ 2:7 s, and t ¼ 4:5 s, laser beam is located in the center of the substrate. Figure 7 shows the temperature field and isotherms of the substrate and deposits from the side view at t ¼ 4:5 s in Case 1. The peak temperature during the process was around 2350 K, while the lowest temperature was close to room temperature. The big temperature differences and small geometrical dimensions caused very large temperature gradients.

Temperature gradient
The temperature gradient involved in the DMD process was quantitatively analyzed in details. The temperature of nodes along the x 0 and y 0 (shown in Figure 8) axes in simulation Case 1 at t ¼ 4:5 s is shown in Figure 9. The x 0 -direction nodes were selected along the top surface of the substrate (bottom surface of the deposits), while the y 0 -direction nodes were selected along the height of the deposits. The temperature of the substrate's top surface reached a maximum of 1069 K just below the center of the laser beam and decreased gradually along the x 0direction. In the y 0 -direction, the temperature of the deposits reached a maximum of 2220 K on the top surface of the deposits and decreased rapidly to1069 K . The slopes of the temperature curves represent the thermal gradients along the x 0 -and y 0 -directions. Along x 0 , the temperature gradient reached a maximum of483 K=mm; along y 0 , the maximum temperature gradient occurred near the top surface of the deposits, reaching 1416 K=mm and then decreasing along the negative y 0 -direction. These steep thermal gradients induced large compressive strains within the deposits and substrates [24].

Residual stress
The nature and magnitude of residual stresses existing in final deposits would affect the integrity of the entire structure. In general conditions, compressive residual stresses are advantageous since they increase the load resistance and     prevent crack growth, while tensile residual stresses are detrimental that they reduce the load resistance and accelerate crack growth.
The residual stress (in Pa) distribution within the final deposits is shown in Figures 10 and 11, where Figure 10 shows the whole substrate and deposit, while Figure 11 shows a y-y cross section view with half of the deposits hidden to show the internal residual stress. Normal stresses σ 11 , σ 22 , and σ 33 along three spatial directions are shown in Figures 10 and 11(a)-(c), respectively, and the von Mises stress is shown in Figures 10 and 11(d). As the figures indicate, residual stresses in the lower part of the deposits were mostly tensile stresses due to the cooldown phase of the molten layers [24]. After the deposition was finished, the remelted lower part of the deposits began to shrink; this shrinkage was restricted by the underlying material, thus inducing tensile stresses. Compressive residual stresses existed at the top free surface of the deposits, caused by the steep temperature gradient. The expansion of the hotter top layer was inhibited by the underlying material, thus introducing compressive stress at the top surface.
Various experimental methods for measuring residual stress have been developed, such as destructive methods, including incremental hole drilling, layer removal, and crack compliance, and nondestructive methods including X-ray diffraction and neutron diffraction [8,9,[24][25][26][27]. These methods could be used to measure the residual stress directly with relatively good accuracy; however, they usually are not cost-effective or easy to set up. Therefore, instead of measuring the residual stress directly, a flexible indirect method has been developed for residual stress validation. A one-one relationship exists between the deflection of the substrate and residual stress; therefore, by validating the deflection of the substrate, the residual stress results can be validated indirectly.

Deformation
During the DMD process, the substrate will continuously expand and shrink, finally maintaining a deformed shape ( Figure 12). In this study, deflection along y was the main deformation under consideration and is shown in Figure 13.

Experiment setup
As shown in Figure 14, in the experiment, the substrate was clamped at the left end to prevent rigid body motion. Keyence's LK-G5000 series laser displacement sensor shown in Figure 15 was placed just below the right end of the substrate to record the displacement of the free end along the y-direction with a frequency of 25 Hz during the process. The experimental results appear in Figure 16. The entire DMD process was controlled by the "Laser Aided Material Deposition System".      Figure 16 illustrates the comparisons of the substrate deflection between the experimental and simulation results for both cases. It is obvious that the simulated deflection matched well with the experimental results. During each deposition layer, the substrate firstly bent down due to thermal expansion on the top surface and then bent up due to thermal shrinkage from the cooling process. After cooling down, the substrate still maintained its distorted shape.

Experimental and simulation results
The deflection values from simulation were 28.5 and 24.6% higher than the deflection measured from experiments for Cases 1 and 2, respectively. There are several potential factors that could cause these differences. Firstly, the experimental setup could not perfectly match the setup in the simulation. For example, in the simulation, the laser beam traveled exactly along the centerline of the substrate. While in experiments (Figure 12), this cannot be perfectly achieved. These offsets may greatly affect the deflection value since it is very sensitive to the positions of heated zone (where expansion and shrinkage mainly happen) and measuring point. Moreover, the laser displacement sensor did not track the displacement of one particular point on the workpiece. Instead, it sensed the signal reflected by an obstacle, so the positions it tracked were always changing as the substrate continuing to deform. Last but not the least, the simplifications and assumptions considered in both thermal and mechanical analyses could also contribute to the differences between the simulation and experiment. For example, although the substrate material is originally isotropic elastic solid, it may become orthotropic after the DMD process.
It is also worth noting that for Cases 1 and 2, even the total amount of energy applied to the substrate is the same and Case 2 has significantly higher distortion than Case 1, which is caused by higher laser power.

Conclusion
To investigate the features of thermal and mechanical behavior of deposited materials involved in the DMD process, a sequentially coupled, thermomechanical finite element model was developed for multilayer DMD process of Stainless Steel 304. The results revealed the characteristics of temperature distribution, residual stress, and deformation within the formed deposits and substrates. A set of experiments were conducted to validate the mechanical effects using a laser displacement sensor. This FEA model can be used to predict the mechanical behavior of products fabricated by the DMD process or similar processes with localized heat sources such as laser sintering, laser cladding, and welding.