Simulation of the solidification of pure nickel via the phase-field method

REGULAR ARTICLES

Simulation of the solidification of pure nickel via the phase-field method

Alexandre Furtado Ferreira^*; Alexandre José da Silva; José Adilson de Castro

Program in Metallurgical Engineering, Universidade Federal Fluminense – UFF/EEIMVR, Av. dos Trabalhadores, n. 420, 27255-125 Volta Redonda - RJ, Brazil

ABSTRACT

The Phase-Field method was applied to simulate the solidification of pure nickel dendrites and the results compared with those predicted by the solidification theory and with experimental data reported in the literature. The model's behavior was tested with respect to some initial and boundary conditions. For an initial condition without supercooling, the smooth interface of the solid phase nucleated at the edges of the domain grew uniformly into the liquid region, without branching. In an initially supercooled melt, the interface became unstable under 260 K supercooling, generating ramifications into the liquid region. The phase-field results for dendrite tip velocity were close to experimental results reported in the literature for supercooling above 50 K, but they failed to describe correctly the nonlinear behavior predicted by the collision-limited growth theory and confirmed by experimental data for low supercooling levels.

Keywords: solidification, phase-field, mathematical modeling, nickel

1. Introduction

Understanding dendritic solidification is highly important because the microstructural scales of the dendrite determine the proprieties of the material. Although significant advances have been achieved in understanding dendritic structures in past decades, our knowledge concerning dendritic growth phenomena is largely based on experiments and idealized theoretical models.

The first well-known paper about the Phase-field Method, proposed by Caginalp¹, dates back to 1980. Since then, the method has proved to be suitable for describing the complex growth of the unstable solid/liquid interface during solidification phenomena. This method is especially attractive because it considers all the governing equations in the whole domain, without direct tracking of the interface position during numerical calculations. Several phase-field models have recently been developed, mainly for the solidification of pure materials^1-6, and have also been extended to the solidification of binary alloys^7-15.

For this type of work involving solidification of pure materials, we can cite Kim et al.⁶. In that paper, the authors present a double-grid method that significantly reduces the CPU time. This method is based on the difference in phase-field diffusivity and thermal diffusivity in a pure material, enabling one to use a large time step. The phase variable (f) is calculated adaptively within the solid/liquid interface and the thermal field is calculated within the thermal boundary layer. The results obtained by Kim et al. point to a well-developed nickel dendrite with tertiary and secondary arms forming side branches. Kobayashi² proposes a similar work for pure materials, showing the influence of the isotropic mode on the solid/liquid interface in a system cooled from the boundaries. In one simulation, the solid/liquid plane interface is stable, with no branching, and advances uniformly into the liquid region. In another, Kobayashi studies the effect of latent heat on the shapes of the solid/liquid interface. For low latent heat in a supercooled system, the plane solid/liquid interface is slightly destabilized. At high latent heat, also in a supercooled system, a cellular structure is visible. One can also see that the branches compete with each other. The anisotropy constant (d_M) is also studied in this paper². Kobayashi shows the relations between d_M and dendrite shapes. In the first computation, where d_M = 0, patterns similar to the amorphous shapes are noted. The results for anisotropy constants different from zero (d_M = 0.010) show one typical dendritic structure with side branches.

The simulation of the solidification process is an intriguing problem. The formation of complex microstructures during the solidification of metals and alloys in the liquid phase have fascinated researchers in materials science and related areas for hundreds of years. A question a researcher in this area might come across is: How does the phase-field method generate results that resemble the morphology of a dendrite? To answer this question, this paper presents a study of some fundamental features of the method under different boundary and initial conditions. In the results section, we first present the equations of the phase-field method solved by finite volumes versus finite differences. This is a strategy adopted to validate our numerical procedure. The second section presents a simulation of solidification without supercooling. In this case, the plane interface is stable. In the third section discusses the simulation of the solidification of a pure metal with supercooling. In this condition, the solid/liquid interface becomes unstable, generating ramifications into the liquid. In the fourth section, the Phase-Field Method is used to predict solidification rates. We calculated the solidification for different supercooling levels from 10 to 160 K. The model successfully reproduced the general behavior and the predicted velocities were of the same order of magnitude found in experimental data. The phase-field model, however, was not able to reproduce the nonlinear dependence between solidification velocities and supercooling levels. The last section presents the development of the thermal field for the system. For the formation of a dendritic morphology, it is important that the thermal diffusivity become greater than its phase diffusivity so that asymmetrical secondary dendrites can be formed.

2. Phase-field Modeling for Pure Materials

The phase-field model is based on the simultaneous solution of energy and phase-field equations. Phase-field modeling assumes the growth of seeds in the liquid phase. According to this hypothesis, there are three regions to be considered: the solid nucleus, the liquid phase and the solid/liquid interface. The state of the entire domain is represented by the distribution of a single variable known as the order parameter, f, or phase-field variable. The solid material is represented by f = 1 while, in the liquid region, f = 0. The region in which f changes from 1 to 0 is defined as the solid/liquid interface. The time evolution equation of the phase-field f is described by⁶:

where M is defined as the solid/liquid interface mobility⁷, the angle q is given by the orientation of a vector perpendicular to the solid/liquid interface, e.g., Ñf. DH is the latent heat and T_m the melting temperature. The function g'(f) that multiplies w determines the distribution of the excess free-energy at the interface. h(f) is a function that satisfies the condition h'(0) = h'(1 ) = 0. As in Ref. 7, we chose

Equations 2 and 3 are widely used in combination with the Phase-field Method. Note that the term disappears when f = 0, in which case only liquid is present. Similarly, when f = +1, only solid is present. As expected, this term is different from zero only in the presence of both the solid and liquid phases. The method most widely used^6,7 to include anisotropy for two-dimensional calculations is to assume that e in Equation 1 depends on q, the orientation of the normal to the interface with respect to the x-axis:

where d_M is the anisotropy constant. The value of j controls the number of preferential directions of the material's anisotropy, equalling 0 for the isotropic cases, 4 for anisotropy of 4 directions, and so on. The constant q_o is the interface orientation with respect to the maximum anisotropy, while e and w are parameters associated with the interfacial energy (s) and interface thickness (l), as proposed by Boettinger et al.⁷:

For the interface mobility, we follow Refs. 6 and 7:

where µ is the linear interface kinetic coefficient. The Phase-field Method particularized for a pure material, subjected to a nonuniform thermal field, includes an energy transport equation:

In Equation 8, D is the thermal diffusivity, DH is the phase-change latent heat, considered positive for solidification, r is the material's density, assumed to be the same for both solid and liquid, and C_p is the specific heat. The parameters and physical properties are presented in the results section.

3. Numerical Solutions

Equations 1 and 8 were solved using the Finite Volume Method¹⁶, in a mesh sufficiently refined to describe details of the dendrites. The energy equation was solved by an implicit scheme, which ensures convergence for longer time steps. For the phase-field equation, an explicit scheme was used. Considering the case of pure nickel, e(q)² < D, convergence is then already possible for values of Dt < Dx² / 4e(q)². However, in order to observe the growth of thermal dendrites, the calculation must be done according to the time scale of the thermal diffusion. For this reason, it was necessary to use Dt < Dx² / 4D.

The anti-symmetrical side branching from primary arms around the dendrite tip is known to be possible only with the existence of a noise source in the phase-field Equation 1. Therefore, random noise was added to Equation 1, in the same way as described by Warren and Boettinger:⁸

r is a random number generated between + 1 and – 1 and a is a noise amplitude factor. An analysis of Equation 9 indicates that the maximum noise is at f = 0.5, with zero at f = 0 and f = + 1. The calculations were made on a Pentium 4, 1200MHz, with 1 GB RAM.

4. Results and Discussion

4.1. The phase-field method by finite volumes vs. phase-field by finite differences

The Phase-field Method has been developed over the two last decades. The equations pertaining to the method are usually solved by finite differences. In this work, Equations 1 and 8 were solved by the finite volume technique¹⁶. The strategy adopted to validate the numerical procedure developed here was to compare our results against those of Kobayashi² for the numerical simulation of a dendritic crystal. The parameters and properties adopted for this simulation were the same as those used by Kobayashi² and are shown in Table 1 and 2. Figure 1 shows the simulation results for a latent heat of 1.6. The boundary condition for f is zero-flux condition. An adiabatic process was assumed for the heat flux. The initial temperature of the domain is equal to - 1.0, which is below the melting temperature (T_m = 0).

In Figure 1, the similarity between the dendrite obtained in this work a) and the one obtained by Kobayashi² b) under the same conditions is evident. Small differences between the secondary arms are due to the random characteristics of the noise source.

Figure 2a depicts the simulation of a nickel dendrite calculated here with a numerical grid of 1200 x 1200 points, while Figure 2b depicts the results obtained by Kim et al.⁶. A "seed" was set at the top left-hand corner, as done before by Kim et al.⁶. Both figures show: a) the secondary and tertiary arms; b) the secondary arm increases with the distance behind the primary dendrite tip; and c) the asymmetry in the side branch found in the secondary and tertiary arms. The initial and boundaries conditions we used in this simulation are identical to those adopted by Kim et al.⁶

4.2. Simulation of the solidification of a pure metal without supercooling

Table 3 summarizes the values of the model's parameters used here to simulate a nickel dendrite, while Table 4 lists the values of the material's properties.

The initial temperature of liquid nickel was set to equal its melting temperature (T_m = 1728 K). The domain discretization was a mesh of 600 x 600 points. The four boundaries were cooled at identical rates. The temperature in the region neighboring the domain is 300 K. The cooling rate obeys the convection law for heat transfer. The regions closer to domain boundaries present lower temperatures, since the system is being cooled at the boundaries at a prescribed heat transfer coefficient of h = 45 Wm^-1 K^-1.

A thin layer of solid (f = 1) was placed on the left, right and top borders of the domain, while a rectangular layer of solid was placed on the bottom border (Figure 3a). The objective of this simulation is to verify that, in a non-supercooled melt, the interface advances uniformly into the liquid region, with no branching. This type of plane interface is stable. These results are illustrated by the time evolution of the interface shown in Figures 3b to 3d.

 

The phase-field model is also able to reproduce solidification conditions in which the solid/liquid interface is stable and no dendrites are formed. Such is the case of a pure metal that undergoes solidification due to heat extraction through the walls of a mould, with no supercooling of the bulk of the liquid. Under these circumstances, even large perturbations in the solid/liquid interface are damped and the solidification front becomes flat. Stabilization of a planar interface in a domain with no supercooling is also reported in Chalmers¹⁷ and is exemplified in Figure 3a-d, which shows the solid/liquid interface at four different moments during solidification. Figure 3a shows the initial condition, where a pattern of rectangles was set as a seed for the solidification front. In Figure 3c-d, this pattern evolves to a planar front.

4.3. Solidification of a pure metal with supercooling effects

In this simulation, a fine solid layer was added to the domain boundaries and there is 260 K supercooling. The objective of this simulation is to demonstrate that the plane interface becomes unstable in a supercooled domain, with branching advancing into the liquid region. Figure 4a-d illustrates the advance of the interface. Note that regions close to domain corners solidify more slowly than other regions, which is due to the fact that there is a greater release of latent heat in the corners, causing the local temperature to increase and thus delaying solidification in this region. Figure 4d shows the primary arms; note that some of the arms located in the central region are well developed, including some secondary arms.

 

4.4. Solidification speed

The phase-field model is usually applied as a tool to calculate the solidification process, but it can also be used to predict solidification rates or speeds. The boundary condition for f is a zero-flux condition. As for the heat flux, an adiabatic process was assumed. We explored this possibility to assess dependence of the solidification speed on supercooling, calculating the solidification of pure nickel at different supercooling levels ranging from 10 to 160 K. For each transient calculation, the displacement of a dendrite tip was measured and the speed of the tip of a dendrite was determined for each supercooling level. The results of these calculations are shown in Figure 5b, while the experimental results are depicted in Figure 5a. In both Figures 5a and 5b, the supercooling considered was the difference in temperature between the equilibrium and the initial liquid. The phase-field model predicted a linear dependence between solidification speed and supercooling levels, while experimental results in the literature¹⁷ indicate a nonlinear behavior. We surmise that the discrepancy can be explained by the fact that the linear behavior is "built into" the formulation of the phase-field model. The linear dependence on supercooling first appears in the source term in Equation 1, but it is also a consequence of the constant mobility M, which is derived based on the assumption of a linear dependence between the solidification speed and supercooling¹⁷. A consequence of this assumption is the introduction of a constant linear kinetic coefficient, µ^k, in Equation 7. The coefficient µ^k should be considered constant only for small ranges of supercooling. Disregarding this dependence could cause the linear behavior shown in Figure 5b. However, this assumption requires further investigation and will be the focus of the upcoming Ref. 18. However, the model predicted a velocity magnitude very close to that of the experimental results, especially at higher supercooling, where the experimental results scatter and do not rule out a linear behavior.

4.5. Development of the thermal field

Figure 6a shows the solidification of initially supercooled (DT = 260 K) pure nickel with adiabatic boundary conditions and zero flux for f, while Figure 6b shows the thermal field. In Figure 6a the dendrite started to grow from a solid nucleus at the top left-hand corner toward the bottom right-hand corner of the domain. One of the interesting phenomena is the asymmetry in the side branching found in the secondary and tertiary arms; side branching occurs only at one side of the arms. Such asymmetry seems to be related with the thermal-field distribution.

Figure 7 shows temperature profiles and the phase-field variable across the interface at the dendrite tip, as indicated in Figure 6a. The transient response of phase-field equations is controlled by the product M.e₀². This parameter acts in the phase-field method similarly to the thermal diffusivity D in the thermal energy equation. As Kim et al.⁶ pointed out, in the formation of a dendritic morphology it is important that the thermal diffusivity becomes greater than its similar term, M.e₀². This can be explained by analyzing Figure 7. The greater value of D forced the thermal front to be always ahead of the solidification interface. Hence, there was always a thermal gradient ahead of a small perturbation at the solidification front. This thermal gradient favored the perturbations, which penetrated this region and then grew more rapidly. Without this thermal gradient, no perturbation would hold and the dendritic pattern would not form.

5. Conclusions

The phase-field model is a powerful tool for simulating details of solidification of crystalline materials. The model can be described as a transport-like equation for the solid phase formulated in terms of a phase-variable, f, which determines whether the phase is solid or liquid. This transport equation was solved by a finite-volume scheme, explicit in time, together with a finite-volume scheme and an implicit-in-time scheme for the energy equation. The simulations of pure nickel dendrites carried out in this work were first compared with previous simulations reported in the literature to validate the present computational code. Calculations of the solidification of pure nickel were then used to investigate the model's ability to reproduce both stable and unstable solid/liquid interfaces. Predictions of the phase-field method for the solidification velocities at dendrite tips were compared with data from the literature¹⁷. Our model was able to reproduce the general behavior and the predicted velocities were of the same order of magnitude found in experimental data. The phase-field model, however, could not reproduce the nonlinear dependence between solidification velocities and supercooling levels. The likely reason for this discrepancy was "built into" the phase-field model due to the use of a constant interface mobility, M, through a kinetic coefficient µ^k, originally postulated as a constant.

Received: April 25, 2005; Revised: November 30, 2006

Appendix

From the hypothesis of equilibrium, a solution exists when T = T_m for a one-dimensional transition zone between liquid (f = 0) and solid (f = 0), where f varies in the x direction normal to the interface. In this case, Equation 1 becomes:

g(f) is here represented as follows:

and

Inserting Equation 12 into 10, one has:

The solution is:

The goal here is to show that Equation 14 is a solution of Equation 13 and that it satisfies the boundary conditions, where f = 1 (solid) when x ® - ¥ and f = 0 (liquid) when x ® + ¥. Rearranging Equation 14,

The relation between the hyperbolic functions is:

Taking Equation 15 into Equation 16,

Differentiating Equation 14 twice:

Using Equations 15 and 16 in 17 and rearranging the terms, one has:

Thus, Equation 19 is equivalent to Equation 13.

• 1. Caginalp G, Fife P. Phase-Field methods for interface boundaries. Physical Review B. 1986; 33:7792-7794.
• 2. Kobayashi R. Modeling and numerical simulations of dendritic crystal growth. Physical D. 1993; 63:410-423.
• 3. Andersson C. Phase-Field simulation of dendritic solidification. [Unpublished PhD Thesis]. Stockholm: Royal Institute of Technology; 2002.
• 4. McFadden GB, Wheeler AA, Braun RJ, Coeriell SR. Phase-Field model for anisotropic interfaces. Physical Review E 1993; 48:2016-2024.
• 5. Wheeler AA, Murray BT, Schaefer RJ. Computation of dendrites using a phase-field model. Physical D 1993; 66:243-262.
• 6. Kim SG, Kim WT, Lee JS, Ode M, Suzuki T. Large scale simulation of dendritic growth in pure undercooled melt by phase-field model. ISIJ International 1999; 39(4):335-340.
• 7. Boettinger WJ, Warren JA, Beckermann C, Karma A. Phase-field simulation of solidification. Annual Review of Materials Research 2002; 32:163-194.
• 8. Warren JA, Boettinger WJ. Prediction of dendritic growth and microsegregation patterns in a binary alloy using the phase-field method. Acta Metallurgical Materials 1995; 43(2):689-703.
• 9. Lee JS, Suzuki T. Numerical simulation of isothermal dendritic growth by phase-field model. ISIJ International 1999; 39(3):246-252.
• 10. Ode M, Suzuki T. Numerical simulation of initial microstructure evolution Fe-C alloys using a phase-field model. ISIJ International 2002; 42(4):368-374.
• 11. Ode M, Lee JS, Kim SG, Kim WT, Suzuki T. Phase-field model for solidification of ternary alloys. ISIJ International 2000; 40(9):870-876.
• 12. Ode M, Kim SG, Kim WT. Suzuki T. Numerical prediction of the secondary dendrite arm spacing using a phase-field model. ISIJ International 2001; 42(4):345-349.
• 13. Ode M, Suzuki T, Kim SG, Kim WT. Phase-field model for solidification of Fe-C alloys. Science and Technology of Advanced Materials 2000; 1:43-49.
• 14. Kim SG, Kim WT, Suzuki T. Phase-field model for binary alloys. Physical Review E 1999; 60(6):7186-7197.
• 15. Kim SG, Kim WT, Suzuki T. Interfacial compositions of solid and liquid in a phase-field model with finite interface thickness for isothermal solidification in binary alloys. Physical Review E 1998; 58(3):3316-3323.
• 16. Patankar SV. Numerical Heat Transfer and Fluid Flow McGraw-Hill-Hemisphere Publication, Minnesota, USA, 1979.
• 17. Chalmers B. Principles of Solidification, USA: John Wiley & Sons; 1964.
• 18. Furtado AF, Castro JA, Silva AJ. Simulation of solidification speed via the phase-field method. [in preparation].

Publication Dates

• Publication in this collection
  22 Mar 2007
• Date of issue
  Dec 2006

History

• Reviewed
  30 Nov 2006
• Received
  25 Apr 2005