Two-dimensional simulation of grain structure growth within selective laser melted AA-2024

A two-dimensional Cellular Automata (CA) – Finite Element (FE) (CA-FE) coupled model has been developed to predictthemicrostructuresformedduringthelasermeltingofapowderedAA-2024feedstockusingtheAdditive Manufacturing(AM)processSelectiveLaserMelting(SLM).ThepresentedCAmodeliscoupledwithathermalFE model, which computes the heat ﬂ ow characteristics of the SLM process. The developed model considers the powder-to-liquid-to-solid transformation, tracks the interaction between several melt pools within a melted track, and several tracks within various layers. It was found that the calculated temperature pro ﬁ les as well as the simulated microstructures bared close resemblance with SLM fabricated AA-2024 samples. The developed model was capable of predicting melt pool cooling and solidi ﬁ cation rates, the type of microstructure obtained, the size of the melt pool (with 14% error) and the heat affected zone, average grain size number (with 12% error) and the growth competition present in microstructures of components manufactured via SLM.


Introduction
Selective Laser Melting (SLM) is an Additive Manufacturing (AM) technology that in recent years has experienced notable increases in industrial uptake for the manufacture of end-use engineering components. The SLM technology can be used to process a wide range of metallic alloys (e.g. nickel, titanium, aluminium based etc.). The melting process is rapid, fusing multiple layers successively together creates complex thermal histories within the material. According to Verhaeghe et al. [25] it is crucial to understand the physical phenomena involved within the fabrication process in order to accurately control it. Materials and Design 113 (2017)   Numerical methods that simulate the SLM process have been undertaken by several researchers (Shiomi, 1999; [3,6,7,10,12,13,16,18,19,23,26]). Each of these numerical approaches attempts to develop an improved understanding of the physical phenomena that occur during the laser processing of a powder bed (thermal history, Marangoni flows, solidification front, etc.).

Materials and Design
Understanding the thermal history of an SLM component has been a main area of investigation. Shiomi et al. [22] developed an FE simulation that calculated the temperature distribution within metallic powders exposed to a pulsed laser and experimentally validated the calculated results. It was found that the maximum temperature reached by the metallic powder was affected by the peak laser power rather than the duration of the laser irradiation. Matsumoto et al. [13] proposed a method to calculate the temperature and stress distribution within a solidified layer within the SLM process using the FE method. Although the effect of a substrate plate was not considered, neither was the detailed thermal dependent properties of the material of interest discussed, it was one of the first research to compute the change from powder-to-liquid-to-solid. Overall, these models are restricted by the complexity of the problem, and the consideration of homogeneous conditions within the model. Roberts et al. [18] developed a threedimensional FEM model in order to understand the thermal history during layer-by-layer processing, while taking into account the nonlinearities produced by the temperature-dependent material properties and phase changes. Even though their results agree with experiments, it was identified that a more detailed model is needed in order to compute the solidification phenomena within SLM. Loh et al. [12] developed a single layer FEM model using a sacrificial layer (which vaporises) in order to obtain accurate results of the generated temperature profiles within the SLM process, however this approach is not considered as suitable for a multilayer process. Foroozmehr et al. [3] used the optical penetration depth of a laser beam [2] and developed a 3D single layer powder bed model that predicts the temperature profile. Detailed thermal dependent properties (excluding mushy zone properties) were considered, and experimental values were used to calibrate the optical penetration depth. However even though the results are considered accurate, the model does not consider the interaction between layers. In an attempt to predict optimal processing parameters during SLM, Song et al. [23] simulated the process on a three-dimensional FE model and experimentally validated the results, highlighting the importance of a FEM simulation of the SLM process. Numerical models of the SLM process are important in order to achieve a certain degree of control/optimisation of the process [19].
The heat transfer phenomenon within a melt pool formed by the SLM process is highly influenced by the fluid flow [15], solely modelling laser melting process without fluid flow consideration will cause inaccuracies. Khairallah and Anderson [10] demonstrated via a three-dimensional mesoscopic micrometre scale model the importance of including the stochastic nature of the powder bed. Further to this it was found that the physics of the process is driven by the surface tension of the melt pools and subsequently effects heat transfer and topology of the solidified melt pools. Other models such as that developed by Pengpeng and Dongdong [16] used Computational Fluid Dynamics (CFD) to accurately predict the melt pool geometry and temperature distributions present within the process. This work developed a three-dimensional model to simulate the temperature evolution behaviour and the effects of the melt pool dynamics during the SLM process and validated the model with experimental trials. In order to simplify and reduce the simulation processing time, a number of researchers have used the enhanced thermal conductivity approach to account for melt pool convection. Safdar et al. [19] used a proposed enhanced thermal conductivity approach for the SLM process and experimentally validated his results. This investigator states that the enhanced thermal conductivity approach is able to artificially simulate the melt pool convection during the processing of materials in SLM, accurately modelling the melt pool profile and temperature distribution without the need of CFD models.
These studies on thermal behaviour have assisted in improving the understanding of stress formation, melt pool topology and surface tension within the SLM process. The majority of SLM microstructural studies have focussed on observations of experimentally fabricated components. Work undertaken by Yin and Felicelli [26] developed a numerical model of the microstructural development during the LENS laser powder blown process with focus only on a micro region inside the mushy zone of the melt pool. The model developed by Yin and Felicelli [26] does not consider convection on the top surface of the layer, and the obtained results are only relevant for the deposition of a single layer.
Numerical simulation has been used to understand grain growth and develop optimum processing conditions in other metal processing techniques (i.e. casting, forging, etc.) to improve efficiency. Despite the benefits numerical simulation has to offer, the development of an appropriate numerical simulation to model microstructural evolution within powder bed SLM has not yet been reported in literature.

Modelling methodology
The purpose of this research is to develop a "first of its kind" microstructural evolution model of the SLM process. The model developed in this work is based on the CA-FE (Cellular Automata-Finite Element) method developed by Gandin and Rappaz [5]. The CA model is generally used to describe the formation of grains during the solidification process, while the FE method is used to calculate the heat flow present within the process. The coupling CA-FE leads to a multiscale model, in which the FE considers the higher scale problem (temperature profiles) and the CA the smaller scale problem (grain growth). The CA-FE technique is used in order to develop a new model which is able to capture the evolution of the microstructural formation during the meltingsolidification of various melt pools within several layers of the SLM process.

Governing equations
The SLM process uses a localised laser beam in order to heat and melt the powder bed, as a result heat transfer plays an important role in this process. Generally, the spatial and temporal distribution of the temperature is governed by the heat conduction equation, which can be expressed as: where T is the temperature, t is the time, x, y and z are the spatial coordinates, k xx , k yy and k zz are the thermal conductivities, ρ is the density, C p is the specific heat and ϕ ⃛ is the heat source term. The heat source is modelled using a Gaussian model, which is the most widely adopted model that uses the symmetrical distribution of laser irradiance across the beam (assuming the irradiance is symmetrical about their propagation direction). This work uses an approximation of the heat source used by Shi et al. [21], expressed as: where P is the power of the laser beam, α is the laser energy absorptance of the material and r is the spot radius. The calculations on the time-dependent temperature distribution during the SLM process were performed with the FE software ANSYS Ver. 14.0. Ansys solves the general energy balance equation in the defined closed system in accordance with the first law of thermodynamics, which is expressed as: R are the vector heat quantities for the laser heat flux, conduction, convection losses and radiation losses respectively. Appropriate values of thermo-physical quantities (e.g. thermal conductivities, specific heat, etc.) must be used in order to attain an accurate and realistic result of the thermal distribution present during the process.

Temperature dependent material properties
The thermal properties of most of the standard materials (e.g. density, thermal conductivity, enthalpy, etc.) are typically available within literature. However, most of the values found in literature correspond to the bulk material properties, for powders different values are found.
In the current work the material investigated is AA-2024, selected due to its associated use and interest in industries such as aerospace. The temperature dependent properties of the bulk material considered are those reported by Mills [14]. In the case of the powdered AA-2024 material, several approaches can be taken, however within this work the discrete powder bed is substituted by a continuum material that possess equivalent material properties. This approach was undertaken in order to compute the effective material properties of the powder bed as that proposed by Sumin Sih and Barlow [24] in which the effective emissivity and thermal conductivity of the powder bed is calculated. The density of the powder at room temperature was experimentally obtained by measuring the mass of a container of a known volume filled with powder and temperature dependent values established by considering the expected behaviour of the powder when heated.
A main part of this research is to develop an accurate simulation of the solidification process that takes place during powder bed laser melting. It is therefore of great importance to consider in detail the physical characteristics (thermophysical properties) for the mushy region of the material. It is known that the thermophysical properties of the liquid phase differ from those of the solid alloy (e.g. thermal conductivity and density) thus the value of the thermophysical property in the mushy region will be dependent on the amount of liquid and solid (i.e. solid fraction). In order to calculate the required property (A) at a certain temperature (T), the expression proposed by Mills [14] is used: where f s(T) is the solid fraction at the desired temperature, A Tsol and A Tliq are the values of the property at the solidus and liquidus temperature respectively. This expression can be used to calculate the density, heat capacity, enthalpy of fusion, thermal conductivity, diffusivity and emissivity in the mushy region of a material. In order to mimic the fluid flow present in a melt pool generated during the SLM process, an enhanced thermal conductivity approach was used [19], where the thermal conductivity λ ii ′′ is defined as: where k is the normal isotropic thermal conductivity value at the corresponding temperature, ii represent the spatial coordinate and γ ii ′′ is the anisotropic enhancement factor for the respective spatial coordinate, which is defined as: In this simulation the anisotropic enhancement factor was used to artificially simulate the melt pool convection during laser melting. Following some experimental trials changing the value for this factor and observing its effects on the melt pool size, it was found that the anisotropic enhancement factor has a relationship with the laser power employed during processing. The values used in this investigation are detailed Section 6. The material data used for the FEM is listed in Table 1.

Simulation model description
The two-dimensional geometry and mesh used is shown in Fig. 1. The model was meshed with rectangular elements (SOLID 77 was used). A 2D thermal conduction solid element with 8 nodes and a single degree of freedom with temperature at each node. The element size was defined (1 μm) after performing a mesh convergence study shown in Fig. 2. A thermal transient analysis is considered, a full Newton-Raphson solver to implicitly solve the stiffness matrix values of the dynamic analysis in Ansys is used. The model consists of four 50 μm high layers of AA-2024 in the form of powder (powder material properties defined) and a block of 0.5 mm of height and 2 mm of width of solid form AA-2024 (solid material properties defined). In order to simulate the heat dissipation along the powder bed (without interfering with the temperature distribution generated by the laser heat source in the scanned region), temperature boundary conditions (simulations of a full powder bed were run in order to determine the temperature boundary conditions, and it was defined as 380 K) were defined (applied to the bottom, left and right walls of the model) in which the laser heat source and convection as a result of air flow is not defined.
As the model reproduces the layer-by-layer building methodology employed by SLM, a recoating time (12 s) (convection boundary condition over the powder bed surface) is considered as the time taken to deposit a new layer of powder prior to subsequent laser scanning. The chamber temperature is maintained within the model at 80°C. The scanning pattern considered is shown in Fig. 3.
Once the simulation commences, the laser (with diameter of 50 μm) is applied as a heat flux ( ϕ⃛ ⃛ ) on the powder bed surface at the corresponding location of the determined laser spot for a defined time (dictated by the exposure time), a convection surface boundary (convection coefficient of air at room temperature, 20 W/m 2 K) was applied to the rest of the powder bed in order to simulate the flow of inert gas inside the chamber. The layers are deposited in the model using the element birth and death technique. As the laser moves through the powder bed, the temperature is monitored to determine which elements have melted and then the material properties of those elements are simultaneously revised and fed back into the model in order to  simulate the solidification process. After the irradiation time of a layer has finished, a recoating time is simulated (this includes thermal convection considerations). This routine is continued layer-by-layer until the last layer is complete.

Cellular Automata model
The nucleation and growth kinetics involved during the solidification process are integrated within the CA grid. The main objective is to mimic the development of grain structures by changing the state indices of the cells in the CA grid. As the nucleation and the growth kinetics are temperature dependent, the temperature values are obtained from the calculated FE nodes.
All the cells (v) of the CA grid are initialised with values of the state index corresponding to the values indicated by the initial time step of the FE nodes. Once the material is in its liquid state and starts to cool down and the local temperature of a cell (ΔT v t ) becomes lower than the critical temperature for the nucleation site (ΔT nucl ), a new grain is created. These undercooling temperatures follow a Gaussian distribution as described by Rappaz [17]. Subsequently a unique grain number is attributed to the cell (v) for each nucleation event and its state index is updated corresponding to a non-liquid state.
In the case of the growth algorithm, the shape associated with the cell (for FCC crystals in 2D, the growing shape is associated with a square [1]) can eventually become large enough to encompass the centre of a neighbouring cell, the growth only occurs if the local temperature corresponds to the defined undercooling temperature. The index state of the neighbouring cell is then switched to a value that corresponds to the growing grain structure. The objective is to preserve the local microstructure-liquid growth front. In the present algorithm, all liquid cells are successively scanned for capture by their direct neighbours, i.e. 4 neighbours in 2D. If a liquid cell is captured by several neighbouring cells during the same time step, a randomly selected neighbouring cell will capture that liquid cell and transform it into solid.

Cellular Automatacoupling with Finite Element
The local temperature or undercooling temperatures of the cells are calculated through FE and are a key parameter of the CA nucleation-growth algorithm. The CA-FE method superimposes the FE mesh to the CA lattice. Usually the grid used to calculate the FE solution is coarser than that of the CA cell size, however, the cell size of both CA and FE can be a similar size if the computing performance is not significantly affected, such is the case for the present work.
In the CA-FE method a weak or a full coupling mode can be used, as described by Gandin et al. [4]. In this work the weak coupling mode was used, in which a unique solidification path (e.g. the Gulliver-Sheil micro-segregation path) can be used on the FE calculation. The variation of enthalpy is then a simple function of temperature variation alone and the temperature field is directly solved on the macroscopic scale. The CA rules defined uses the values of temperature distributions calculated by the FE at certain time intervals.

Thermal history simulation validation
In order to determine the suitability of the developed model, experimental work was conducted in order to compare the results obtained by the model. Several experiments with constant layer thickness (50 μm) were undertaken, varying processing parameters (as shown in Table 2) in both the simulation and the experimental trials. Consequently, the melt pool size (height and diameter) in different locations of the benchmark samples was measured and compared to that of the obtained simulation in order to calibrate the model. SLM samples were then created to validate the results of the numerical simulation.
The SLM samples used for metallographic inspection were grinded and polished according to standard procedures and etched using Keller's reagent (solution consisting of 190 ml distilled water, 5 ml HNO 3 , 3 ml HCl, 2 ml HF). A Carl Zeiss inverted optical microscope was used to observe and measure the formed melt pool geometries. Melt pools were measured in different locations of the sample and an average size of melt pool was determined per specimen. Fig. 4a shows the measured melt pool size of sample number 1, the dashed lines in Fig. 4a represent an approximate location of the deposited layer, approximately 3 to 5 measurements (shown in different colours in Fig.  4a) to the clearest melt pools per micrograph were performed to the ten micrographs obtained per sample. Fig. 4b shows the predicted melt pool size obtained with the developed numerical model for sample number 1 processing conditions. Fig. 5 shows a comparison between the experimental and the simulated data for melt pool geometry. The presented data suggests that a 14%   error is present between the experimental and modelling results, thus the prediction of the melt pool dimensions of parts produced in the SLM system can be performed within this limit.

Thermal history analysis
From the 2D layer-by-layer FEM model, valuable data can be extracted. Information such as cooling and solidification rates of the melt pool and porosity can be used.
Jacobson and McKittrick [9] estimated that the cooling rates for RSP are in the range of 10 5 to 10 6 K/s. Kurz and Trivedi [11] outlined the solidification conditions of a process similar to SLM, laser surface processing; and established a relationship between the cooling rate (j _ Tj), the thermal gradient in the liquid ahead of the solid-liquid interface (G) and the interface growth rate (V), which is defined as j _ Tj ¼ G Ã V. It was determined that the solidification conditions of that of laser processing will in most the cases, lead to a columnar (directional) growth.
In Table 3 the cooling rate, the thermal gradient in the liquid, the solidification rate and the calculated mean V is shown.
The data extracted from the developed FEM model suggests that the resulting microstructure of the sample, will be in all of the cases a columnar growth. The simulated results also agree with calculated values of cooling rates determined by Harrison et al. [8], which was determined using measurements of primary dendrite arm spacing of parts produced with SLM. Using the values extracted from the developed FEM model, the primary dendrite arm spacing on parts produced by SLM could also be predicted.

Microstructural evolution, simulation and validation
A similar validation process to that used in Section 6.1 is used to determine the suitability of the developed simulation model. Using the set of parameters shown in Table 2, benchmark specimens were produced and results were inputted into the simulation, validation SLM samples were experimentally fabricated and compared with the simulation results.
The SLM samples were grinded and polished according to standard procedures and anodized with Baker's reagent (solution of 1.8% HBF 4 in water, using 20 VDC for 80 s and an agitation velocity of 10 RPM at 22°C). An inverted optical microscope with polarised light was used to observe and measure the formed grains in the melt pools of the SLM processed specimens.
The average grain size number was measured in random locations of the specimens along the transverse direction of the elongated grains in order to determine an average grain size number per specimen according to the intercept procedures described in ASTM E112 standard. Measurements along the length direction of the grains were not performed, since the numerical model will be only representative for four layers and the experimental values will be representative for the whole sample size. The orientation of the grains as well as the development of a grain within layers was observed and visually compared with the obtained simulated results. Fig. 6a shows an example of the obtained micrograph of sample number 1 as well as of how the average grain size measurements were performed. An average of 50 measurements per sample was performed in order to determine the average grain size number. The calculations performed by the developed CA-FE model were exported into a .bmp file in order to have a visual comparison with experimental data, different colours (red, green, blue and pink) will determine different grains. An example of the calculations performed using the data numerically calculated for sample number 1 is shown on Fig. 6b.
The average grain size number was then measured in each of the predicted microstructures in order to compare the results with   Fig. 7 compares the experimental average grain size number with the predicted average grain size number. The shown data suggests that a 12% error between measurements and predictions is present. The average grain size number of parts produced via SLM can be determined considering this limit.

Microstructural analysis results of the multilayer model.
The developed CA-FE model uses the calculated temperature profiles by the developed FEM model in order to predict the microstructural evolution of a simulated part produced by SLM. The developed CA-FE model uses data provided by fabricated benchmark samples in order to successfully predict the microstructure. The data inputted to the developed code is in form of a probability of nucleation based on the measured grain size. The probability of nucleation is a parameter that depends of the material to be used, in the case of the present work this probability was established as 0.025 based on experimental observations and trial and error runs of the code. The probability of nucleation was also calculated using a theory proposed by Sands [20], however the calculated value is several orders of magnitude lower (0.0003) than the above-mentioned value. Using the probability of nucleation calculated with the theory proposed by Sands [20] on the developed CA-FE model, the obtained microstructures do not agree with experimental information.
The average grain size of the predictions made by the developed CA-FE model agrees with experimental information obtained from experiments. Nevertheless, making a close comparison of the predicted microstructure with the microstructure obtained from experiments (see Fig. 8), columnar grains that grow between layers with smaller equiaxed type interspersed grains can be observed in both microstructures.
New grain boundaries intersecting primary columnar grains are occasionally formed at the limits of each melt pool within the prediction, comparable to the observed phenomenon in the experiments performed, as well as in the research undertaken by Harrison et al. [8]. New formed grains continue to grow and competitive growth has then an important role in the layer-by-layer process as the heat flux reduces.

Conclusions
A novel approach for a 2D FEM layer-by-layer model was developed. The developed model considers as inputs the most important processing parameters required in the manufacturing of parts via SLM (e.g. laser power, point distance, exposure time, etc.) and material properties (e.g. absorptance, density, thermal conductivity, etc.) in order to predict the temperature profiles generated within the process.
The temperature profiles calculated by the novel 2D FEM layer-bylayer model in this investigation serves as an input for the CA-FE coupling. The developed CA-FE model uses the temperature profiles calculated by the developed FEM approach in order to simulate the solidification phenomenon present within the SLM process. After the calculations were performed by the CA-FE model, the calculated microstructures were found to be similar to the microstructures of components manufactured via SLM.
Benchmarking samples were produced in order to experimentally validate both FEM and CA-FE models. The melt pool dimensions and the average grain size for the benchmark samples were measured and the obtained results were used as an input to calculate both the anisotropic enhancement factors (FEM) and the nucleation density (CA-FE) of the models. Simulations were performed and the results were compared with those of the validation samples. The melt pool dimensions and the average grain size number were used in order to validate the calculated temperature profiles and microstructures. It was found that melt pool dimensions predicted by the FEM model had an error of approximately 14% compared to experimental values, the average grain size number had a 12% error. After the validations were performed, it was concluded that both FEM and CA-FE models successfully predict both the temperature profiles and microstructures of components manufactured via SLM.
From the developed 2D FEM layer-by-layer model, important data was extracted. From the simulation, cooling and solidification rates were calculated, confirming that the process is considered as a RSP technique. Using these rates and the GV microstructure selection map, it was determined that the formed microstructure would compose of dendritic growth. Porosity, a common defect within parts manufactured via SLM, was predicted by the developed FEM model as a result of lack of fusion. These predictions can be used to avoid such defects by running multiple simulations with variable processing conditions in order to minimise porosity. As further work, further validations on porosity predictions is required in order to fully rely on the results obtained from the FEM model. The calculated microstructures by the CA-FE model held similarities with those components manufactured via SLM. Calculated microstructures had similar growth phenomena to those found in experiments. The formation of small equiaxed grains at the liquid interface of the formed melt pools were predicted, which then will compete with larger columnar grains and grow towards the general heat flux. This growth competition led to the appearance of either small dispersed grains or columnar grains that halted growth between layers due to these dispersed grains. These same effects were observed in experiments, concluding that the calculated microstructure by the developed CA-FE model agree with experiments.
Given that the developed models agreed with experimental data, it is expected that such models could be used with other materials and expanded to a three-dimensional space in order to predict generated temperature profiles and microstructure of components manufactured with SLM. Understanding paths for the development of customised microstructures could be a further benefitting from the use of this developed model.