Micro-macroscopic Coupling in the Cellular Automaton Model of Solidification

A cellular automaton (CA) model to predict the formation of grain macrostructure during solidification has been implemented and the coupling between the microscopic and the macroscopic submodels has been investigated. The microscopic submodel simulates the nucleation and growth of grains, whereas the macroscopic solves the heat conduction equation. The directional solidification of an Al-7 wt. (%) Si alloy was simulated, enabling the calculation of the temperature and solid fraction profiles. The calculated temperature was used to obtain the solid fraction profile by an application of Scheil equation. This solid fraction disagrees with that calculated in the micro-macro coupling of the model, although this coupling is completely based on Scheil equation. Careful examination of the discrepancies shows that it is a result of the undercoolings for nucleation and growth of grains and also of the interpolations of enthalpy change and temperature from the finite volume mesh to the CA cell mesh.


Introduction
The cellular automaton (CA) technique is an algorithm constructed to simulate the spatial/temporal evolution of a system by applying simple transformation rules to the sites of a lattice 1 , which is a set of fixed points.The lattice is used to map some important parameter of the system, such as atomic particles, grains of a metallic alloy, etc.Each point of the lattice must have a state, selected among a set of finite states that represent important information about the phenomena under study.In the case of modeling the grain structure of a given material, the points of the lattice might have a number that indicates the orientation of the grains.Therefore, lattice points within the same grain would have the same number and, an approximate image of all the grain structure could be constructed from the states and positions of all lattice points.
In the CA algorithm, provided the state values of all lattice points are known at the beginning of the simulation, the time evolution of these values is calculated after each time step using transformation rules based on the state of a lattice and its neighbors and, sometimes, on parameters such as temperature and concentration, depending on the type of model.The CA technique has been used in different areas of materials science.Some important applications are the mathematical modeling of the microstructure evolution during recrystallization 2 , sintering 3 , and solidification 4 .Further applications are given by Raabe 1 .One important result of these mathematical models based on the CA is an image of the microstructure during its evolution and, obviously, the final microstructure.
Gandin and Rappaz 4 proposed one of the first models of solidification based on the CA technique.This model evolved into the popular CAFE (cellular automaton -finite element) model 5 , which was applied to predict grain macrostructures of alloys in the as-cast state 6 .Although modified CA models have been proposed to simulate detailed microstructures 7 , the CAFE model is still frequently used, because it is more computationally efficient when only the grain structure (without the details of dendrite arms) is required.
In most solidification models, including CAFE, the calculation of the solid fraction field depends on the coupling between microscopic and macroscopic submodels.The time evolution of the solid fraction field during solidification is of paramount importance to: the latent heat release, the solute rejection by the solid to the liquid, and the fluid flow in the mushy zone.Nevertheless, in the cellular automaton model the calculation of this field has never been examined, because the main interest has always been to obtain an image of the grain macrostructure.The objective of this work is to examine the micromacroscopic coupling and the resulting solid fraction field calculated with the cellular automaton model.

Cellular Automaton Model
A stochastic model based on CAFE 5 was developed in the present work by implementing and coupling a microscopic and a macroscopic submodel.A brief description of the model is presented here, emphasizing the important features of the coupling between the submodels (micro-macro coupling), summarized in Figure 1.
The macro-model (macroscopic submodel) consists of the numerical solution by the explicit finite volume method 8 of the following one-dimensional energy conservation equation where H is the volumetric enthalpy; T is the temperature; t is time; y is the spatial coordinate in the direction of heat flow; and κ is the thermal conductivity, defined as an average of the conductivities of the liquid (κ l ) and solid (κ s ) phases, weighted by their corresponding local volume fractions.Equation 1 was discretized by a one-dimensional mesh of equal rectangular finite volumes with a node at the centre of each volume and aligned in the y direction.From the discretized equation, a volumetric enthalpy change DH t + δt was calculated for each finite volume after a time step δt.
The micro-model (microscopic submodel) was developed mainly to simulate the nucleation and growth of dendritic grains in two dimensions during solidification.A two-dimensional microscopic mesh of CA square cells with a site at each cell centre was adopted.Each CA cell had two possible states, namely, liquid and solid.
At the beginning of simulation, all cells were liquid and particles representing substrates for heterogeneous nucleation of grains were randomly located in some of these cells.Nucleation of a new grain occurred at a particle when the liquid cell reached a predetermined local undercooling DT N = T L -T N , where T L and T N are the liquidus and nucleation temperatures, respectively.In the present work, the same DT N value was adopted for all particles, which is equivalent to assume instantaneous nucleation of grains 9 .At nucleation, the cell state is changed from liquid to solid and a square with a random orientation and centered at the cell site is created to represent a new grain envelope.
The growth of a newly created grain is simulated in the micromodel by increasing the diagonals of the square at a velocity V = A DT m , where DT is the undercooling at the CA site in relation to T L , and A and m are constants dependent on the alloy (Figure 2a).After sufficient growth, some part of the square eventually reached the position of one of the four nearest neighbor sites of the twodimensional mesh (Figure 2b).Then, the state of the corresponding neighbor cell was changed to solid and a new square was created with an orientation and size defined as in the CAFE model 5 to represent another part of the same growing envelope (Figure 2c).This whole process is shown schematically in Figure 2.
The strong coupling between the micro-model and macro-model is established in the solid fraction calculation.Immediately after nucleation at a liquid CA cell of the micro-model, its solid fraction e s is calculated by Scheil equation 10 , given below: where T f is the melting point of the pure metal and k is the solute partition coefficient.The temperature T is obtained by interpolation, at the position of the cell site, of the T values at the four nearest nodes of the finite volume mesh (macro-model).During grain growth, after each time step (δt) the change in the solid fraction (δe s t + δt ) is calculated at every CA cell by an equation derived from a combination of Scheil equation and the change in volumetric enthalpy at the cell site position, δH t + δt , as follows 4 ( 2) ( )( 1) 1 where ρ is the alloy density; c p is the specific heat; and L f is the latent heat of fusion.The value of δH t + δt at the cell site position is also obtained by a linear interpolation of DH t + δt , calculated at the nodes of the finite volume mesh (macro-model).During solidification, the interpolated cell temperature eventually reaches the eutectic temperature, T E .At this moment, the temperature is kept constant at T E and the change in solid fraction at the cell is calculated as follows 4 At a finite volume of the macroscopic model, the change in solid fraction, De s t + δt , is obtained as in CAFE 5 , i.e., by averaging the values of δe s t + δt at all cells within this finite volume, weighted by the same interpolation coefficients used to interpolate T and δH t + δt .Finally, the temperature at a node of the finite volume mesh (macro-model) is calculated by the following equation, completing the micro-macro coupling cycle (Figure 1).
Results obtained with the computer code implemented in the present work showed excellent agreement with those presented by Rappaz and Gandin 4,11 for the CAFE model.

Analysis of the Micro-macro Coupling
The unidirectional solidification of an Al-7 wt.(%) Si alloy was simulated in a rectangular two-dimensional domain of size 0.15 m in the heat flow direction (y) and 0.05 m perpendicular to it (x).An adiabatic boundary was assumed at the top of the domain (y = 0.15 m), while through the bottom (y = 0), the heat flux was given by q = h(T -T ∞ ), where h = 250 W.m -2 .K -1 is the heat transfer coefficient and T ∞ = 298 K is a reference temperature.A uniform temperature of 991 K was assumed as the initial condition.The number density of substrate particles for the nucleation of grains was defined as n T,bulk = 5 × 10 5 m -3 for the bulk liquid and n T,layer = 3 × 10 6 m -2 for the liquid layer adjacent to the lower boundary (x = 0).Since the mesh of CA cells is two-dimensional, these values were converted by 11 ( ) The alloy properties were: ρ = 2452 kg.m -3 ; c p = 1126 J.kg -1 .K -1 ; κ l = 60.5 W.m -1 .K -1 ; κ s = 137.5 W.m -1 .K -1 ; L f = 387400 J.kg -1 ; k = 0.13; T L = 891 K; T E = 850 K; T f = 933 K; A = 3 × 10 -6 m.s -1 .K -m ; m = 2.7; and DT N = 3K.Convergence of the numerical solution was achieved with a finite volume mesh (macro-model) of 30 volumes and a CA mesh (micro-model) of 300 (y) × 100 (x) cells.
The calculated grain structure and the profiles of e s and T at the finite volume nodes (indicated by dots) are shown in Figure 3.The solid fraction given by Scheil model, e s Scheil , was obtained with Equation 2, in which T is the temperature at the finite volume nodes.It is possible to see that e s Scheil does not agree with the e s profile calculated with the CA model equations (Equations 2, 3 and 4), which were based on the Scheil model.
To confirm the discrepancy between e s Scheil and e s , the implemented CA code was adjusted to simulate the unidirectional solidification problem studied by Guillemot et al. 12 , in which the growth of only one grain, positioned at the CA cell adjacent to the lower boundary, is simulated,.Therefore, after nucleation the grain envelope grows throughout the domain from the lower to the upper boundary.The calculated T and e s as a function of time are compared with Guillemot et al. 12 results in Figure 4.The solid fraction e s Scheil obtained with Scheil equation (Equation 2) using the calculated temperature at 80 mm from the lower boundary is significantly different from e s obtained with both the implemented model and Guillemot et al. 12 model, confirming the type of discrepancy observed in Figure 3.

Effect of Undercooling and Interpolation on the Micro-macro Coupling
Three types of discrepancy are observed in Figure 3 and Figure 4b: (1) e s Scheil > e s = 0 at the undercooled liquid ahead of the solid-liquid interface; (2) e s Scheil < e s in the mushy zone, down to T E ; and (3) e s < 1 below T E , although solidification should end at T E .The first discrepancy exists because the Scheil equation (Equation 2), used to obtain e s Scheil , neglects the undercooling for nucleation and growth of grains 10 .Therefore, it predicts e s Scheil > 0 immediately below T L .In the stochastic model, however, the nucleation undercooling DT N (= 3K) is partly responsible for the undercooled liquid ahead of the solid-liquid interface.The other part is the undercooling necessary for grains to grow up at a velocity V, approximately equal to the T L isotherm velocity.
To confirm these effects, the nucleation undercooling was removed (DT N = 0) and the growth undercooling was eliminated by inserting a substrate particle within every CA cell of the mesh.In these conditions, the profile of e s calculated with the CA model shows very good agreement with e s Scheil near T L (Figure 5), but some discrepancy still persists within the mushy zone and below T E .
This discrepancy is related to the interpolations of DH t + δt and δe s t + δt to calculate δH t + δt and De s t + δt , respectively (Figure 1).The Scheil model is implicitly given by Equation 3; consequently,δe s t + δt calculated at a CA cell is exactly related to δH t + δt by the Scheil model.Nevertheless, the node temperature is obtained from Equation 5, using DH t + δt and De s t + δt , which are not related by Equation 3. Actually, the relation between DH t + δt and De s t + δt is complex, because De s t + δt is obtained from the weighted average of δe s t + δt values at all    12 .Also shown is the solid fraction e s Scheil calculated from Equation 2, in which T is the temperature curve in (a) at x = 80 mm.To eliminate the interpolation effect described above, the number of finite volumes in the y direction was increased and made equal to the number of CA cells in the same direction, resulting in δe s t + δt = De s t + δt and δH t + δt = DH t + δt .Now, De s t + δt , DH t + δt , and the temperature at the finite volume node are all linked by the Scheil model.The nucleation and growth undercoolings were also removed as in the previous simulation.Figure 6 shows that the resulting e s profile closely agrees with that from the Scheil model, e s Scheil , confirming the interpolation effect.
Examining the equations used to calculate the solid fractions e s Scheil (Equation 2) and e s (Equation 3) might give futher insight into the discrepancies discussed in the previous paragraphs.Equation 3 is a differential form of Scheil equation, in which a relationship between δT and δH (based on thermodynamics) was substituted.Equation 2, which is the traditional Scheil equation, can actually be obtained after numerical integration of Equation 3, using as initial condition that solidification begins at the liquidus temperature of the alloy, i.e., e s = 0 for T = T L .When this Equation 3 is used in the CA microscopic submodel, however, this initial condition is not used, because when T = T L at a CA cell, nucleation or grain growtn might not have activated the site of this cell.Therefore, solidification does not necessarily begin locally at the liquidus temperature, a behavior that differs from the hypothesis adopted in the tradicional Scheil equation.Furhtermore, another reason for the discrepancy is that e s Scheil is calculated at the finite volume node, directly from its temperature, whereas e s is obtained from an average of solid fraction changes calculated in all cells located within the finite volume, as explained in section 2.

Summary and Conclusions
A cellular automaton (CA) model to predict the grain macrostructure formation during solidification has been implemented and the micro-macro coupling of this model investigated.In the coupling cycle, after each time step of the numerical method, a change in solid fraction is calculated at the CA cells by a Scheil model.The solid fraction (e s ) in the finite volume is calculated from an averaged change in all cells located within this volume.A second type of solid fraction (e s Scheil ) was also calculated in the finite volume using its temperature; this fraction was compared with that obtained from the micro-macro coupling (e s ).Three types of important discrepancies were observed: (1) e s Scheil > e s = 0 at the undercooled liquid ahead of the solid-liquid interface; (2) e s Scheil < e s in the mushy zone, down to the eutectic temperature (T E ); and (3) e s < 1 below T E , although solidification should end at T E .These discrepancies were examined and their cause were seen to be the undercoolings for nucleation and growth of grains and the interpolations of enthalpy change and temperature from the finite volume mesh to the CA cell mesh, required in the micro-macroscopic coupling.

Figure 1 .
Figure 1.Diagram showing the exchange of information in the coupling cycle between the micro and macro submodels (micro-macro coupling).

Figure 2 .
Figure 2. Description of the a) growth of a CA square; b) activation of a neighbor CA site; and c) creation of a new square associated with the newly activated CA site.

Figure 3 .
Figure 3. Grain macrostructure, solid fraction, and temperature profile calculated with the stochastic model during solidification (t = 670 seconds).The solid fraction (e s Scheil) is calculated using Equation2, in which T is the temperature at the finite volume nodes of the macroscopic model.

Figure 4 .
Figure 4. a) Temperature and b) solid fraction as a function of time at different distances from the lower boundary calculated with the implemented cellular automaton model and presented by Guillemot et al.12 .Also shown is the solid fraction e s Scheil calculated from Equation2, in which T is the temperature curve in (a) at x = 80 mm.

Figure 5 .
Figure 5. Grain macrostructure, solid fraction, and temperature profile calculated with the stochastic model during solidification (t = 670 seconds) when the nucleation and growth undercoolings are negligible.The solid fraction (e s Scheil) is calculated using Equation2, in which T is the temperature at the finite volume nodes of the macroscopic model.

Figure 6 .
Figure 6.Grain macrostructure, solid fraction, and temperature profile calculated with the stochastic model during solidification (t = 670 seconds) when the nucleation and growth undercoolings are negligible and using a finite volume mesh as refined as the CA cell mesh.The solid fraction (e s Scheil ) is calculated using Equation 2, in which T is the temperature at the finite volume nodes of the macroscopic model.