Predicting Mesoscale Microstructural Evolution in Electron Beam Welding

Using the kinetic Monte Carlo simulator, Stochastic Parallel PARticle Kinetic Simulator, from Sandia National Laboratories, a user routine has been developed to simulate mesoscale predictions of a grain structure near a moving heat source. Here, we demonstrate the use of this user routine to produce voxelized, synthetic, three-dimensional microstructures for electron-beam welding by comparing them with experimentally produced microstructures. When simulation input parameters are matched to experimental process parameters, qualitative and quantitative agreement for both grain size and grain morphology are achieved. The method is capable of simulating both single- and multipass welds. The simulations provide an opportunity for not only accelerated design but also the integration of simulation and experiments in design such that simulations can receive parameter bounds from experiments and, in turn, provide predictions of a resultant microstructure.


INTRODUCTION
The Potts Monte Carlo model, a generalization of the Ising model for magnetic systems to systems with an arbitrary number of spins or grain identifiers, has been used for several decades to study grain growth in polycrystalline materials. [1][2][3][4] The simplicity of the method has led to the development of many in-house codes for specific applications but few widespread, general-use, publically available ones. Stochastic Parallel PARticle Kinetic Simulator (SPPARKS) is the result of an effort to develop a general-use mesoscale Monte Carlo suite at Sandia National Laboratories, much like the LAMMPS package 5 for atomic-scale simulations, thereby providing a framework for the broad implementation of a variety of mesoscale Monte Carlo solvers. By using SPPARKS, the Potts model has been extended to several novel applications, including hybrid approaches incorporating cellular automata and phase field models. [6][7][8][9] Since the 1990s, several studies have attempted to simulate the influence of materials processing methods with moving heat sources on the microstructure. Dress et al. 10,11 used cellular automata to produce an accurate grain structure within the fusion zone (FZ) of an autogenous weld. However, the relatively small scale and two-dimensional (2D) nature of the method limited further application or quantitative analysis. Significant advances were made by Debroy et al. with a three-dimensional (3D) weld simulation. [12][13][14][15] Their studies demonstrated good agreement with experimental grain size distributions at various distances from the centerline of the FZ. Debroy and co-workers focused their analysis on the heat-affected zone (HAZ) surrounding the FZ and excluded the FZ from their simulation domain.
Here, we introduce a Monte Carlo Potts-based method for the simulation of melting, solidification, and grain growth during welding at the mesoscale. The model simulates the widely known dependence of solidification behavior on thermal gradient (G) and solidification front velocity (V) in metals, rather than on specific material system properties. As shown in Chapter 4 of Kurz and Fisher, 16  solidification based on the combination of the thermal gradient, G, and solidification velocity, V. The model presented in this work simulates melting in the FZ and grain growth in the FZ and HAZ to match grain morphology predictions based on G and V. A 3D steady-state temperature profile with a given G is rastered with a velocity V. Since these simulations are performed at a length scale that is much larger than required to resolve the formation and growth of dendrites, the authors suggest this approach enables studying large numbers of grains and their evolution, which in turn, provides greater flexibility and a potentially broader range of applicability.

Grain Growth
The Monte Carlo Potts model evolves spins (or grain identifiers) on a discrete lattice to simulate microstructural evolution. To initialize the simulation, the starting microstructure is digitized on a 3D lattice by assigning spins to each lattice site. Since the driving force for curvature-driven grain growth is the reduction in grain boundary energy, only grain boundary energy is considered and is given by the sum of the bond energies between neighboring lattice sites of unlike spins: where N is the total number of lattice sites in the simulation, q i and q j are the spins at lattices sites i and j, and L is the number of neighbors of each lattice site (26 for a 3D cubic lattice used here). The 1/2 prefactor eliminates double counting in the summation. With this formulation, each pair of unlike neighboring sites contributes one unit of energy to the total system energy. Grain growth is simulated by selecting a lattice site and attempting to change its spin to that of a randomly selected neighbor with a different spin. The total energy of the system is then recalculated with Eq. 1, and the probability of accepting the change is determined by the Metropolis function: 18 where DE is the change in system energy, k B is Boltzmann's constant, and T s is the Potts simulation temperature (which does not correspond directly to physical system temperature 8 ). With this formulation, all configurational changes that decrease the global system energy are accepted, while those that increase system energy are accepted based on the Boltzmann distribution.
When all lattice sites have attempted one spin change, simulation time is advanced by one Monte Carlo step (MCS) and the process repeats. This formulation has been shown to simulate curvaturedriven grain growth with correct kinetics and topology. 3

Model Modification for Simulation of Welding
Several adaptations to simulate welding have been made to the standard Potts model to incorporate a localized, moving heat source 9 that melts a small fusion zone, FZ, and heats the surrounding region called the heat-affected zone, HAZ. The temperature T in the moving FZ and HAZ are input into the Potts model as a function of time, t, and location, i, T(t, i). Note that this temperature, T, is the physical temperature not to be confused with the simulation temperature, T s , of Eq. 2. The FZ is the region where T(t, i) > T m , the melting temperature, and melting within the FZ is simulated by disordering the spins, so that neighboring sites do not have the same spin. This effectively simulates a molten zone with high configurational energy given by Eq. 1 in which no sites are bonded to one other.
In the solid HAZ, grain growth is simulated as described earlier with the mobility, M, of grain boundaries now made a function of temperature, T(t, i), given by: where Q is the activation energy for grain boundary motion, and M o is a user-defined prefactor. This nonuniform mobility modifies the Monte Carlo acceptance probabilities at each lattice site given by Eq. 2 by adding M as a prefactor: where M is determined by Eq. 3.
In the HAZ, grain growth is simulated with grain boundary mobility being highest near the FZ and decreasing with temperature away from the FZ until M = 0 far from the weld where temperature remains at the ambient temperature. As the heat source moves, previously solid regions with a grain structure are melted by disordering their spins in the FZ. As the heat source continues to move, regions in the FZ cool below T m and grains grow with grain boundary mobility, M(T). In this way, we simulate melting, solidification, and grain growth during welding.
Simulation time is then advanced in a rejection kinetic Monte Carlo scheme 9 such that during each MCS, every site has attempted a change. In this model, the velocity of the heat source is defined as the number of lattice sites traveled by the FZ per MCS.

Experimental Process Parameters as Model Inputs
The newly developed method also allows the user to specify several process parameters, heat source geometry, size and shape of the FZ and HAZ, scan velocity, and number of heat source passes. For common autogenous weld conditions, melt pools are roughly ellipsoidal in shape. 19 Prior experimental investigations have shown at low scan velocities, the surface tension of the liquid interface is minimized, and the melt pool has a symmetric, near-hemispherical shape. As velocities increase, an elongated tail develops behind the melt pool and the overall shape becomes more ''comet-like.'' At high beam powers, local temperatures near the center of the melt pool can surpass the boiling point of the base metal, and a vapor cavity can develop. This leads to deep penetration of the incident beam into the material and the formation of a ''keyhole'' shape where a wide melt pool near the surface gives way to a slimmer, extended melt through the thickness of the component. 20,21 These variations in the characteristic mode of the molten zone are simulated in the model through three-dimensional description of the molten zone geometry.
The welding model described earlier was applied to a starting microstructure consisting of equiaxed grains. The 3D melt pool and surrounding temperature gradient along with cross-sectional schematics of its geometry are shown in Fig. 1. The geometry was designed to model a keyhole weld, and it was created by the union of two ellipsoids. The simulation parameters were chosen to match the keyhole conditions in welding experiments on Ni-200 and are given in Table I. One ellipsoid was broad and shallow, and its three axes were defined by the shallow melt pool length, width, and depth. The second ellipsoid was defined by the deep melt pool length and width and by the height of the simulation domain. The extent of the temperature gradient surrounding the melt pool was defined by two concentric ellipsoids (one shallow and one deep) centered on the shallow melt pool's midpoint. The deep melt pool was positioned so that its leading edge coincided with that of the shallow melt pool; however, the deep temperature gradient shells are centered on the shallow midpoint. This results in an asymmetry in the deep HAZ are as shown in Fig. 1d. Previous work has shown that little temperature variation occurs in front of the melt pool during keyhole welding. 21 The resulting microstructures were analyzed and compared with experimental microstructures obtained by autogenous welding of Ni-200.

ELECTRON-BEAM WELDING VALIDATION EXPERIMENTS
To validate the weld model, autogenous electronbeam butt-welds were made on plates of Ni-200. Beam power density and diameter were measured by using the Enhanced Modified Faraday Cup technique. 22 The beam diameter had a full width at half maximum of 0.16 mm and a 1/e 2 value of 0.42 mm. The beam was 5 mA over-focused. Additional experimental parameters are summarized in Table I. The beam conditions were chosen to enable keyhole-mode welding with penetration of the melt through the entire plate. The weld samples were sectioned, polished, and analyzed via optical microscopy. The resulting experimental microstructures were then used as guidance for the size and shape of the simulation's heat source. Single-pass welds were made, characterized, and compared to simulation results.

Qualitative Comparison
Top-down views of the experimental and simulated electron-beam welds are shown in Fig. 2. Both images show the three primary regions of the welded microstructure: base metal, HAZ, and FZ. In both experiment and simulation, grains in the HAZ isotropically coarsen while in the FZ, grains grow much larger and elongate in the scan direction.
Transverse views of the experimental and simulated weld microstructures are shown in Fig. 3. In the FZ, grains are much larger than in the base metal but not as elongated or curved as in the top- down view. The keyhole shape of the FZ is apparent as it has a maximum width at the top surface, which decreases at lower heights through the plate thickness, while HAZ width remains constant through the depth of the plate.

Quantitative Comparison
Experimental micrographs were segmented and skeletonized to isolate grain boundaries with the Fiji image analysis software suite. 23,24 Quantitative analysis of experimental and simulated digitized microstructures was performed by using Dream3D. 25 The grain equivalent circle diameter (ECD) 26 and aspect ratio distributions of grains seen in the top view of Fig. 1 within the FZ are shown in Fig. 4. The ECD was determined by measuring each grain's cross-sectional area and then calculating the diameter of a circle with the identical area. The aspect ratios were determined by fitting ellipses to the grain cross sections and calculating the ratio of their major and minor axes. The ECD and aspect ratio distributions of grains within the FZ are shown, respectively, in Fig. 4a and b. Both experimental and simulated size distributions are heavily skewed to the larger sizes with few large, outlier grains. These large grains are found at the center of the weld pass and can be readily seen in Fig. 2.
The grain aspect ratio distributions shown in Fig. 4b also show good agreement between simulation and experiment. A circular grain in the cross section would have a value of 1, while a needle-like  grain would have a value near 0. The distributions show that most FZ grains are far from equiaxed with a peak between 0.2 and 0.4. The details of 3D grain geometry can be easily determined for the simulated microstructure, but they require costly and time-consuming serial sectioning or comparable analysis techniques for experimental samples. 27 As an example, the simulated FZ grains were fitted with ellipsoids, and the ratios for the three ellipsoid axes are plotted in Fig. 5. Each data point on the plot represents the morphology of an individual grain. This plot shows the distribution of 3D grain shapes in which grains occurring in the upper right corner are spherical, those in the bottom right are infinitely thin disks, and grains in the lower left corner are needle-like. Due to the convention adopted here, c < b < a, all grains must exist in the space represented by c/b < 1. Overall, the FZ grain shapes are elongated to needle-like shape and not flatted like disks. 26

Influence of Initial Grain Size
In addition to process-related variables, the initial grain structure may influence simulation behavior. We investigate the influence of the initial grain size on the weld simulation microstructure by varying the initial average cross-sectional grain area from 1 to 66 voxels while holding all other welding variables constant. The resulting microstructures were cross-sectioned along the scan direction in the FZ and HAZ, and the grain was cross-sectional area measured. The results are shown in Fig. 6. The average FZ and HAZ grain area remained roughly constant for small initial sizes, and it increased linearly with the initial area when the initial grain size was comparable to the weld dimensions. In welding experiments, initial grain sizes are typically much smaller than those of the solidified structure. In simulated welds, initial grain size relative to the weld dimensions should be made as small as possible to achieve realistic results.
Examining the simulated microstructures shown in Fig. 7 provides an explanation for the increase in solidified grain size. At the conditions simulated, grain growth is driven by surface curvature, but grain boundary mobility varies locally within the HAZ temperature gradient to give different rates of grain growth. The final sizes of HAZ grains are controlled by the residence time in and the temperature of the HAZ. 28 HAZ grains that grow large enough to encroach on the FZ serve as nucleation sites for elongated FZ grains. This is shown by the termination of FZ grains in isotropic surfaces at the HAZ/FZ boundary. When initial grain size is comparable to the resulting solidified microstructure, the initial size is large enough to influence the initial growth of the FZ and HAZ grains. A large initial grain size reduces the number of unique grains in the HAZ and limits available nucleation sites for FZ grains. However, competition between the grains ensures that further growth follows similar trends regardless of the initial size. Thus, the initial grain size does influence final grain sizes, but it has little influence upon grain shapes.

DISCUSSION
Significant materials modeling advances have been made at the nano (atomistic) 29 and continuum length scales. 30 However, limited advancements have been made at the mesoscale. Phase field modeling provides an excellent tool to model material behavior at larger-than-atomic length scales and across interfaces. 31,32 However, the computational demands of phase field limit the size of the Fig. 6. The influence of initial grain size on the resulting simulated FZ and HAZ microstructure. All measurements are in voxels. Fig. 7. Simulated weld microstructures with varying initial grain sizes. Grain morphology in the fusion zone and heat-affected zone were relatively unaffected until the initial grain size was comparable to that of the solidified structure. simulation to 10 or 100s of grains. 33,34 The Potts Monte Carlo model can simulate mesoscale behavior for 1000s, and even 10,000s, of grains. The computational efficiency of the method creates the opportunity of modeling large experimental domains to understand, predict, and design the welding process.
Finite element crystal plasticity and other mechanical response models can incorporate grain shape and texture in continuum-scale simulations through the use of representative volume elements (RVEs). 30,35 These methods require an input microstructure determined by other methods, either from experimental measurements (EBSDs), 36 mathematically generated structures (e.g., Voronoi tessellation), 37 or simplified ideal grain shapes. 38 The Monte Carlo Potts model can provide not only an additional route to bridge microstructure to continuum scales but also can provide realistic microstructures generated from specific weld process simulation. Atomistic-scale methods can be used to determine input parameters, and phase field models can generate detailed initial microstructures, while MC outputs can be used as inputs to continuumscale models.

CONCLUSION
The Monte Carlo Potts model has been adapted to simulate microstructural evolution near a moving heat source with localized melting and has been demonstrated by applying it to electron beam welds. The 3D temperature profile around the melt pool and the solidification front velocity serve as inputs to this model. Grain growth is modeled by varying the grain boundary mobility as a function of the local time-dependent temperature in the heat-affected zone around the moving heat source to simulate the unique elongated and curved grains that result from welding. Quantitative comparison of grain size and grain shape of electron beam welds with experimental cross sections show good agreement giving confidence in the predictive capability of this model. The flexibility of the method was demonstrated by showing the effect of the initial grain size on final welding grain size and shape. Furthermore, the simulations provide full 3D microstructures that can serve as baseline surrogates for 3D measures that would be either costly or destructive to procure.

ACKNOWLEDGEMENT
Sandia National Laboratories is a multiprogram laboratory managed and operated by Sandia Corporation, a wholly owned subsidiary of Lockheed Martin Corporation, for the U.S. Department of Energy's National Nuclear Security Administration under Contract DE-AC04-94AL85000. The authors would also like to thank A. Kilgo for metallographic preparation, M. Winters for e-beam welding assistance, and S. Williams for metallography.

OPEN ACCESS
This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.