- Email: [email protected]
Orientation selection of equiaxed dendritic growth by three-dimensional cellular automaton model
Physica B 407 (2012) 2471–2475
Contents lists available at SciVerse ScienceDirect
Physica B journal homepage: www.elsevier.com/locate/physb
Orientation selection of equiaxed dendritic growth by three-dimensional cellular automaton model Lei Wei, Xin Lin n, Meng Wang, Weidong Huang State Key Laboratory of Solidification Processing, Northwestern Polytechnical University, Xi’an 710072, PR China
a r t i c l e i n f o
abstract
Article history: Received 30 December 2011 Received in revised form 13 March 2012 Accepted 17 March 2012 Available online 22 March 2012
A three-dimensional (3-D) adaptive mesh refinement (AMR) cellular automata (CA) model is developed to simulate the equiaxed dendritic growth of pure substance. In order to reduce the mesh induced anisotropy by CA capture rules, a limited neighbor solid fraction (LNSF) method is presented. It is shown that the LNSF method reduced the mesh induced anisotropy based on the simulated morphologies for isotropic interface free energy. An expansion description using two interface free energy anisotropy parameters (e1, e2) is used in the present 3-D CA model. It is illustrated by present 3-D CA model that the positive e1 favors the dendritic growth with the /100S preferred directions, and negative e2 favors dendritic growth with the /110S preferred directions, which has a good agreement with the prediction of the spherical plot of the inverse of the interfacial stiffness. The dendritic growths with the orientation selection between /100S and /110S are also discussed using the different e1 with e2 ¼ 0.02. It is found that the simulated morphologies by present CA model are as expected from the minimum stiffness criterion. & 2012 Elsevier B.V. All rights reserved.
Keywords: Dendritic growth Solidification Cellular automaton Interfacial energy anisotropy
1. Introduction The solidification microstructures, which generally determine the final properties of the relative product to a large extent [1,2], have been focused by the physicists and material scientists for several decades. Dendrite [3–8] is the most common microstructures observed in the casting products of commercial alloys. According to the solvability theory [3–6], a dendrite tip is determined by a stability parameter, which is a function of the anisotropy parameter of the interface free energy. Recently, molecular dynamics (MD) simulations [9] have established that the parameterization of interface free energy for fcc metals requires two anisotropy parameters (e1, e2). Phase-field simulations [10] have shown that dendritic preferred growth directions vary continuously from /100S to /110S, for anisotropy parameter e1 reduced from 0.15 to 0.0 while anisotropy parameter e2 is fixed to 0.02. However, for a large region sandwiched between the /100S and /110S regions (0.05o e1 o0.12, e2 ¼ 0.02), the phase-field results are not as expected from the minimum stiffness criterion. Cellular automaton (CA) model [11–24] is another computational model for describing dendritic growth, and it has been used to simulate many solidification phenomena, such as columnar-to-
n
Corresponding author. Tel.: þ86 29 88494001; fax: þ86 29 88492374. E-mail address: [email protected] (X. Lin).
0921-4526/$ - see front matter & 2012 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.physb.2012.03.048
equiaxed transition (CET) [15] and the dendritic growth under forced flow [16]. Pan and Zhu [17] presented a 3-D cellular automaton model for describing solutal dendritic growth, in which a Cahn–Hoffman vector is used to calculate the capillarity undercooling. In his model, only dendritic growth with /100S preferred direction could be simulated, because only one anisotropy parameter e4 is used. However, he did not quantitatively examine the effects of e4 on the dendritic growth behaviors. Yin and Felicelli [18] developed a cellular automaton model for dendritic growth with six-fold symmetry. He found that the use of hexagonal elements could reduce the mesh induced anisotropy in dendritic growth with six-fold symmetry. The ‘‘rules of capturing interface cells’’ used in his CA model set six preferred directions firstly, and then take the consideration of anisotropy parameter of interface free energy. If the interface free energy is isotropic, the capture rules in his CA model may not be suitable. Physically, the dendritic growth directions are caused by the anisotropy of interface free energy, so there is no need to set the preferred growth directions in the CA model. Due to the mesh induced anisotropy, few of the CA models [11–19] have quantitatively investigated the effects of interface free energy anisotropy in dendritic growth. As a partner and coauthor of the CA model [19,20], it is realized that the mesh induced anisotropy is the most important problem in CA model. A few years later, the capture rules in CA model were crucially investigated, and a random zigzag capture was presented to reduce mesh induced anisotropy
2472
L. Wei et al. / Physica B 407 (2012) 2471–2475
[21,22]. Based on the previously developed CA model, a CA model with interface reconstruction method [23] was developed, in which all the simulated dendrites with four-fold and six-fold symmetry were caused by physical anisotropy of interface free energy, and the mesh induced anisotropy in 3-D CA model was discussed for future work. Based on the work of 3-D CA model [23], it is found that it was difficult for the 3-D CA model to use 2-D zigzag capture rule [21,22], therefore a limited neighbor solid fraction (LNSF) method was developed for the 3-D CA model to reduce mesh induced anisotropy [24]. Compared to arXiv preprint [24], in the present paper we added some more discussion on the LNSF method as well as a few more figures to illustrate the simulated results by present 3-D AMR CA model.
where r is the distancefrom the interface to the coordinate origin, ^ g is the interface free energy as a function of interface normal n. According to Taylor’s work [26], the capillarity undercooling can be related to the weighted mean curvature (wmc), which is equivalent to the negative of the surface divergence of the x-vector. " # X ^ 1 @2 gðnÞ ^ þ wmc ¼ gðnÞ ¼ divS x ð5Þ 2 Ri @y i ¼ 1,2 i
Thus, the Gibbs–Thomson equation based on the weighted mean curvature can be expressed as TI ¼ TM
TM wmc L
ð6Þ
Combining Eq. (2) with Eq. (4), the surface divergence of the
x-vector could be derived as 2. Model description 2.1. Thermal field During 3-D equiaxed dendritic growth, the thermal field is governed by the following equation ! @T @2 T @2 T @2 T L @f s ¼a ð1Þ þ þ 2 þ @t C p @t @x2 @y2 @z where T is temperature, t is time, a is the heat diffusivity, L is latent heat of solidification, Cp is specific heat and fs is solid fraction. To solve Eq. (1), a central finite difference method was used. 2.2. Anisotropic interface free energy and interface stiffness The atomistic calculations [9] have established that the parameterizations of anisotropic interface free energy require two anisotropy parameters, as expressed by Eq. (2) ! ! 3 3 X X ^ gðnÞ 3 17 ¼ 1 þ e1 n4i n4i þ66n21 n22 n23 þ e2 3 þ ::: ð2Þ 5 7 g0 i¼1 i¼1 where g0 is the averaged interface free energy, n^ ¼ ðn1 ,n2 ,n3 Þis the interface normal, e1, e2 are the two anisotropy parameters. Eq. (2) was given in this form by Hoyt et al. [9]. It is just the beginning of an expansion in orthogonal fully symmetrized cubic harmonics, as for instance stressed by Haxhimali et al. [10] in their discussion of phase-field simulations. According to extremum stiffness criterion, the minima of the interface stiffness, S¼ g þd2g/dy2, which appears in the Gibbs– Thomson condition, indicates the dendritic growth directions [10]. The interface stiffness expressed in spherical angular coordinates [10] is S ¼ 2g þ
@2 g 2
@y
þ
1 2
@2 g
sin y @f2
þcot y
@g @y
ð3Þ
where y and j are the spherical angles of the interface normal. To solve Eq. (3), it is needed to substitute n^ ¼(cosysinj, sinysinj, cosj) into Eq. (2). Then, the plot of 1/S could be drawn for different anisotropy parameters. 2.3. Gibbs–Thomson equation and the weighted mean curvature Hoffman and Cahn [25] introduced a x-vector to represent the anisotropic interface free energy of a sharp interface. The definition of Hoffman–Cahn x-vector is as follow:
x grad ðrgÞ
ð4Þ
3 17 @nx @ny @nz þ þ e2 divS x ¼ 1 e1 5 7 @x @y @z @nz 2 @nx 2 @ny þ ny þn2z þ 12ðe1 þ 3e2 Þ nx @x @x @x @nx @ny @nz þ þ 3ðe1 þ 3e2 ÞQ @x @y @z @Q @Q @Q þ ny þ nz 3ðe1 þ 3e2 Þ nx @x @y @z 1 @R 1 @R 1 @R þ þ þ 132e2 nx @x ny @y nz @z ! 1 @nx 1 @ny 1 @nz 132e2 R 2 þ 2 þ 2 nx @x ny @y nz @z @R @R @R þ ny þnz 330e2 nx @x @y @z @nx @ny @nz 330e2 R þ þ @x @y @z
ð7Þ
where nx ¼qxfs/9rfs9,ny ¼ qyfs/9rf9 and nz ¼ @z f s =9rf s 9,Q ¼ n4x þ n4n þ n4z ,R ¼ n2x n2y n2z . To solve Eq. (7), the partial derivatives of solid fraction, such as: @x f s ,@2x f s ,@2x,y f s . . .must be accurately calculated. An actuate method to solve the partial derivatives of solid fraction for 2-D CA model is discussed in our previous work [21,22]. In the present article the method is extended into 3-D application. 2.4. Cellular automaton rule For 3-D simulation, each interface cell has 6 neighbor cells in /100S directions, 12 neighbor cells in /110S directions and 8 neighbor cells in /111S directions. So, the neighborhood rule of the cell in 3-D is much more complicated than that in 2-D. It is very difficult to present a 3-D capture rule with lower mesh induced anisotropy. So in the present work the simplest capture rule—the Von Neumann type (6 cells in /100S directions) is used, and a limited neighbor solid fraction (LNSF) method was used to control the mesh induced anisotropy. The capture rule in 3-D CA model with LNSF method is described as follows: firstly, a limited neighbor solid fraction fsLNSF was set; secondly, the average solid fraction of the neighbor cells for each liquid cell fsave was calculated. If the fsave is larger than fsLNSF, and this liquid cell has a solidified neighbor cell in the Von Neumann type, then it can be captured as an interface cell. It can be seen that the capture rule in the present 3-D CA model is a combination of the Von Neumann capture rule and the LNSF method. The Von Neumann capture rule makes sure that the CA model presents a sharp interface, and the LNSF method makes sure that all the liquid cells are captured under the similar conditions, thus the mesh induced anisotropy could be reduced.
L. Wei et al. / Physica B 407 (2012) 2471–2475
2473
Fig. 1. Simulated solidifying morphologies in the undercooled pure SCN melt, DT¼ 6 K, Dx ¼0.5 mm, the value of fsLNSF for LNSF method was set to: (a) 0.0, (b) 0.15, (c) 0.225, (d) 0.25.
Fig. 2. Simulated morphologies during solidification from the undercooled pure SCN melt, DT¼ 4.5 K, Dx ¼0.5 mm: (a) e1 ¼ 0.0, e2 ¼0.0 (b) e1 ¼0.10, e2 ¼0.0, (c) e1 ¼0.0, e2 ¼ 0.02; (d) the AMR mesh of (b) in 2D slice.
2474
L. Wei et al. / Physica B 407 (2012) 2471–2475
2.5. Computational details The equiaxed dendritic growth from undercooled pure succinonitrile (SCN) melt was simulated using 3-D CA model, which was programmed using adaptive mesh refinement (AMR) technique. The thermal physical parameters of SCN used in the simulation are the same as reference [22]. The simulation domain is 0.128 0.128 0.128 mm3, and the mesh size is smaller than 0.5 mm. At the beginning, the central cell of computational domain is initialized as a solid cell and the states for the rest cells are set to liquid. The time step is calculated as follow:
Dt ¼ 0:15
Dx2
a
ð16Þ
where Dt is the time step, Dx is mesh size. The initial temperature and boundary temperature of the domain are set to T0, thus the undercooling of the melt can be expressed as DT¼Tm T0.
3. Simulation results and discussion 3.1. The testing of mesh induced anisotropy In order to test the mesh induced anisotropy, the morphologies under the isotropic interface free energy (e1 ¼0.0, e2 ¼0.0) were first simulated and the limited neighbor solid fraction fsLNSF was set from 0.0 to 0.25. Fig. 1 shows the simulated morphologies grown from the pure SCN melt with the undercooling of 6 K. All of Fig. 1(a)–(d) were simulated by the Von Neumann capture rule, and the differences between them were the limited neighbor solid fraction fsLNSF. We found that when fsLNSF ¼0.1, the simulated result was the same as with fsLNSF ¼0.0. And when fsLNSF ¼0.3, it was too large for capture, which means no liquid cells can be captured. So we did the simulation with fsLNSF ¼0.15 and fsLNSF ¼ 0.25 for comparison. It can be seen that for fsLNSF ¼0.0, as seen in Fig. 1(a), all the 24 primary branches grew in /100S directions, which is caused by the Von Neumann capture rule. And forfsLNSF ¼0.15, the simulated morphology was still grown in /100S direction with branches, as seen in Fig. 1(b). For fsLNSF ¼0.225, the simulated morphology presents 24 primary branches with the same size evenly pointing at all the directions, as seen in Fig. 1(c). Fig. 1(d) shows the simulated morphology with /111S preferred growth directions for fsLNSF ¼0.25. It is found that when fsLNSF gets smaller, the mesh induced anisotropy prefers /100S directions. Especially, if fsLNSF is equal to 0.0, the mesh induced anisotropy is the same as that caused by the Von Neumann capture rule. However, when fsLNSF gets larger, the mesh induced anisotropy prefers /111S directions. The fsLNSF ¼0.225 presents a balance between the two mesh induced anisotropies. It can be concluded that the LNSF method with fsLNSF ¼0.225 could reduce the mesh induced anisotropy. So in the following simulations, such as Figs. 2 and 3, we used the LNSF method with fsLNSF ¼0.225.
Fig. 3. Comparison of morphologies predicted by CA model and 1/S plot for different e1 with the e2 of 0.02 for the growth from undercooled pure SCN melt, DT¼ 2.0 K, Dx¼ 0.25 mm.
(both of them are simulated with isotropic interface free energy), it can be seen that as the undercooling changes, the simulated morphology changes as well. From Fig. 2(b) and (c), It can be seen that for interface anisotropy parameters e1 ¼0.10 and e2 ¼0.0, a dendrite grew with the /100S preferred direction, and for e1 ¼0.0 and e2 ¼ 0.02, a dendrite grew with the /110S preferred direction. Fig. 2 (d) show the dendrite in 2-D slice of Fig. 2(b) illustrated by AMR mesh. It can be seen that the solid–liquid interface is covered by the most refined meshes, which can describe the interface more accurately.
3.3. Dendrite orientation selection The simulated morphologies by present 3-D CA model are further compared with the prediction of the spherical plot of the inverse of the interfacial stiffness, as expressed in Eq. (3), for anisotropy parameter e1 changes from 0.0 to 0.16 with the fixed e2 of 0.02. From Fig. 3, it can be seen that the simulated morphologies by CA model are continuously changed with the e1 as expected from the minimum stiffness criterion.
4. Summary 3.2. Dendrite morphology by 3-D CA model Fig. 2 shows the simulated morphologies grown from the pure SCN melt under the undercooling of 4.5 K. The interface anisotropy parameters are set to e1 ¼0.0, e2 ¼0.0 for Fig. 2(a), e1 ¼0.10, e2 ¼0.0 for Fig. 2(b), and e1 ¼0.0, e2 ¼ 0.02 for Fig. 2(c), respectively. When e1 ¼0.0, e2 ¼ 0.0, as seen in Fig. 2(a), the simulated morphology is a kind of cube with instability in directions of /100S, /110S and /111S. Compared Fig. 2 (a) with Fig. 1(c)
A 3-D cellular automaton model for describing the solidification of a pure substance was developed, in which a limited neighbor solid fraction (LNSF) method was presented to reduce the mesh induced anisotropy. In order to simulate the orientation selection, an expansion description using two interface free energy anisotropy parameters (e1, e2) is used in the present 3-D CA model. When e1 ¼0.10, e2 ¼0.0, the dendrite grew with the /100S preferred directions;
L. Wei et al. / Physica B 407 (2012) 2471–2475
and when e1 ¼0.0, e2 ¼ 0.02, the dendrite grew with the /110S preferred directions. The dendrite orientation selection between /100S and /110S dendrites were compared with the spherical plot of the inverse of the interfacial stiffness, using the different e1 with e2 ¼ 0.02. It is found that the simulated morphologies by present CA model are as expected from the minimum stiffness criterion.
Acknowledgment This project was supported by the National Natural Science Foundation of China (Nos. 50971102 and 50901061) and the National Basic Research Program of China (No. 2011CB610402). The work was also supported by the Program of Introducing Talents of Discipline to Universities (No. 08040) and the fund of the State Key Laboratory of Solidification Processing in NWPU (Nos. 39-QZ-2009 and 02-TZ-2008). References [1] W.J. Boettinger, S.R. Coriell, A.L. Greer, A. Karma, W. Kurz, M. Rappaz, R. Trivedi, Acta Mater. 48 (2000) 43.
2475
[2] M. Asta, C. Beckermann, A. Karma, W. Kurz, R. Napolitano, M. Plapp, G. Purdy, M. Rappaz, R. Trivedi, Acta Mater. 57 (2009) 941. [3] S.C. Huang, M.E. Glicksman, Acta Metall. 29 (5) (1981) 701. [4] J.S. Langer, Rev. Mod. Phys. 52 (1) (1980) 1. [5] D.A. Kessler, J. Koplik, H. Levine, Phys. Rev. A 31 (1985) 1712. [6] J.S. Langer, Phys. Rev. A 33 (1986) 435. [7] W.D. Huang, Q.M. Wei, Y.H. Zhou, J. Cryst. Growth 100 (1–2) (1990) 26. [8] X. Lin, W.D. Huang, J. Feng, T. Li, Y.H. Zhou, Acta Mater. 47 (11) (1999) 3271. [9] J.J. Hoyt, M. Asta, A. Karma, Mater. Sci. Eng. R 41 (2003) 121. [10] T. Haxhimali, A. Karma, F. Gonzales, M. Rappaz, Nat. Mater. 5 (2006) 660. [11] L. Nastac, Acta Mater. 47 (17) (1999) 4253–4267. [12] L. Beltran-Sanchez, D.M. Stefanescu, Metall. Mater. Trans. A 34 (2) (2003) 367. [13] M.F. Zhu, D.M. Stefanescu, Acta Mater. 55 (2007) 1741. [14] W. Wang, P.D. Lee, M. McLean, Acta Mater. 51 (10) (2003) 2971. [15] H.B. Dong, P.D. Lee, Acta Mater. 53 (3) (2005) 659. [16] M.F. Zhu, S.Y. Lee, C.P. Hong, Phys. Rev. E 69 (6) (2004) 061610. [17] S.Y. Pan, M.F. Zhu, Acta Mater. 58 (2010) 340. [18] H.B. Yin, S.D. Felicelli, Modell. Simul. Mater. Sci. Eng. 17 (2009) 075011. [19] B.W. Shan, X. Lin, L. Wei, W.D. Huang, Adv. Mater. Res. 26–28 (2007) 957. [20] L. Wei, X. Lin, B. Shan, W. Huang, Proceedings of the 5th Decennial International Conference on Solidification Processing Sheffield, U.K (2007) 53. [21] X. Lin, L. Wei, M. Wang, W.D. Huang, Mater. Sci. Forum PRICM7 (2010) 1528. [22] L. Wei, X. Lin, M. Wang, W.D. Huang, Appl. Phys. A 103 (2011) 123. [23] L. Wei, X. Lin, M. Wang, W.D. Huang, Comput. Mater. Sci. 54 (2012) 66. [24] L. Wei, X. Lin, M. Wang, W.D. Huang, arXiv preprint (2011) arXiv:1111.6159. [25] D.W. Hoffman, J.W. Cahn, Surf. Sci. 31 (1972) 368. [26] J.E. Taylor, Acta Metall. Mater. 40 (1992) 1475.
Contents lists available at SciVerse ScienceDirect
Physica B journal homepage: www.elsevier.com/locate/physb
Orientation selection of equiaxed dendritic growth by three-dimensional cellular automaton model Lei Wei, Xin Lin n, Meng Wang, Weidong Huang State Key Laboratory of Solidification Processing, Northwestern Polytechnical University, Xi’an 710072, PR China
a r t i c l e i n f o
abstract
Article history: Received 30 December 2011 Received in revised form 13 March 2012 Accepted 17 March 2012 Available online 22 March 2012
A three-dimensional (3-D) adaptive mesh refinement (AMR) cellular automata (CA) model is developed to simulate the equiaxed dendritic growth of pure substance. In order to reduce the mesh induced anisotropy by CA capture rules, a limited neighbor solid fraction (LNSF) method is presented. It is shown that the LNSF method reduced the mesh induced anisotropy based on the simulated morphologies for isotropic interface free energy. An expansion description using two interface free energy anisotropy parameters (e1, e2) is used in the present 3-D CA model. It is illustrated by present 3-D CA model that the positive e1 favors the dendritic growth with the /100S preferred directions, and negative e2 favors dendritic growth with the /110S preferred directions, which has a good agreement with the prediction of the spherical plot of the inverse of the interfacial stiffness. The dendritic growths with the orientation selection between /100S and /110S are also discussed using the different e1 with e2 ¼ 0.02. It is found that the simulated morphologies by present CA model are as expected from the minimum stiffness criterion. & 2012 Elsevier B.V. All rights reserved.
Keywords: Dendritic growth Solidification Cellular automaton Interfacial energy anisotropy
1. Introduction The solidification microstructures, which generally determine the final properties of the relative product to a large extent [1,2], have been focused by the physicists and material scientists for several decades. Dendrite [3–8] is the most common microstructures observed in the casting products of commercial alloys. According to the solvability theory [3–6], a dendrite tip is determined by a stability parameter, which is a function of the anisotropy parameter of the interface free energy. Recently, molecular dynamics (MD) simulations [9] have established that the parameterization of interface free energy for fcc metals requires two anisotropy parameters (e1, e2). Phase-field simulations [10] have shown that dendritic preferred growth directions vary continuously from /100S to /110S, for anisotropy parameter e1 reduced from 0.15 to 0.0 while anisotropy parameter e2 is fixed to 0.02. However, for a large region sandwiched between the /100S and /110S regions (0.05o e1 o0.12, e2 ¼ 0.02), the phase-field results are not as expected from the minimum stiffness criterion. Cellular automaton (CA) model [11–24] is another computational model for describing dendritic growth, and it has been used to simulate many solidification phenomena, such as columnar-to-
n
Corresponding author. Tel.: þ86 29 88494001; fax: þ86 29 88492374. E-mail address: [email protected] (X. Lin).
0921-4526/$ - see front matter & 2012 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.physb.2012.03.048
equiaxed transition (CET) [15] and the dendritic growth under forced flow [16]. Pan and Zhu [17] presented a 3-D cellular automaton model for describing solutal dendritic growth, in which a Cahn–Hoffman vector is used to calculate the capillarity undercooling. In his model, only dendritic growth with /100S preferred direction could be simulated, because only one anisotropy parameter e4 is used. However, he did not quantitatively examine the effects of e4 on the dendritic growth behaviors. Yin and Felicelli [18] developed a cellular automaton model for dendritic growth with six-fold symmetry. He found that the use of hexagonal elements could reduce the mesh induced anisotropy in dendritic growth with six-fold symmetry. The ‘‘rules of capturing interface cells’’ used in his CA model set six preferred directions firstly, and then take the consideration of anisotropy parameter of interface free energy. If the interface free energy is isotropic, the capture rules in his CA model may not be suitable. Physically, the dendritic growth directions are caused by the anisotropy of interface free energy, so there is no need to set the preferred growth directions in the CA model. Due to the mesh induced anisotropy, few of the CA models [11–19] have quantitatively investigated the effects of interface free energy anisotropy in dendritic growth. As a partner and coauthor of the CA model [19,20], it is realized that the mesh induced anisotropy is the most important problem in CA model. A few years later, the capture rules in CA model were crucially investigated, and a random zigzag capture was presented to reduce mesh induced anisotropy
2472
L. Wei et al. / Physica B 407 (2012) 2471–2475
[21,22]. Based on the previously developed CA model, a CA model with interface reconstruction method [23] was developed, in which all the simulated dendrites with four-fold and six-fold symmetry were caused by physical anisotropy of interface free energy, and the mesh induced anisotropy in 3-D CA model was discussed for future work. Based on the work of 3-D CA model [23], it is found that it was difficult for the 3-D CA model to use 2-D zigzag capture rule [21,22], therefore a limited neighbor solid fraction (LNSF) method was developed for the 3-D CA model to reduce mesh induced anisotropy [24]. Compared to arXiv preprint [24], in the present paper we added some more discussion on the LNSF method as well as a few more figures to illustrate the simulated results by present 3-D AMR CA model.
where r is the distancefrom the interface to the coordinate origin, ^ g is the interface free energy as a function of interface normal n. According to Taylor’s work [26], the capillarity undercooling can be related to the weighted mean curvature (wmc), which is equivalent to the negative of the surface divergence of the x-vector. " # X ^ 1 @2 gðnÞ ^ þ wmc ¼ gðnÞ ¼ divS x ð5Þ 2 Ri @y i ¼ 1,2 i
Thus, the Gibbs–Thomson equation based on the weighted mean curvature can be expressed as TI ¼ TM
TM wmc L
ð6Þ
Combining Eq. (2) with Eq. (4), the surface divergence of the
x-vector could be derived as 2. Model description 2.1. Thermal field During 3-D equiaxed dendritic growth, the thermal field is governed by the following equation ! @T @2 T @2 T @2 T L @f s ¼a ð1Þ þ þ 2 þ @t C p @t @x2 @y2 @z where T is temperature, t is time, a is the heat diffusivity, L is latent heat of solidification, Cp is specific heat and fs is solid fraction. To solve Eq. (1), a central finite difference method was used. 2.2. Anisotropic interface free energy and interface stiffness The atomistic calculations [9] have established that the parameterizations of anisotropic interface free energy require two anisotropy parameters, as expressed by Eq. (2) ! ! 3 3 X X ^ gðnÞ 3 17 ¼ 1 þ e1 n4i n4i þ66n21 n22 n23 þ e2 3 þ ::: ð2Þ 5 7 g0 i¼1 i¼1 where g0 is the averaged interface free energy, n^ ¼ ðn1 ,n2 ,n3 Þis the interface normal, e1, e2 are the two anisotropy parameters. Eq. (2) was given in this form by Hoyt et al. [9]. It is just the beginning of an expansion in orthogonal fully symmetrized cubic harmonics, as for instance stressed by Haxhimali et al. [10] in their discussion of phase-field simulations. According to extremum stiffness criterion, the minima of the interface stiffness, S¼ g þd2g/dy2, which appears in the Gibbs– Thomson condition, indicates the dendritic growth directions [10]. The interface stiffness expressed in spherical angular coordinates [10] is S ¼ 2g þ
@2 g 2
@y
þ
1 2
@2 g
sin y @f2
þcot y
@g @y
ð3Þ
where y and j are the spherical angles of the interface normal. To solve Eq. (3), it is needed to substitute n^ ¼(cosysinj, sinysinj, cosj) into Eq. (2). Then, the plot of 1/S could be drawn for different anisotropy parameters. 2.3. Gibbs–Thomson equation and the weighted mean curvature Hoffman and Cahn [25] introduced a x-vector to represent the anisotropic interface free energy of a sharp interface. The definition of Hoffman–Cahn x-vector is as follow:
x grad ðrgÞ
ð4Þ
3 17 @nx @ny @nz þ þ e2 divS x ¼ 1 e1 5 7 @x @y @z @nz 2 @nx 2 @ny þ ny þn2z þ 12ðe1 þ 3e2 Þ nx @x @x @x @nx @ny @nz þ þ 3ðe1 þ 3e2 ÞQ @x @y @z @Q @Q @Q þ ny þ nz 3ðe1 þ 3e2 Þ nx @x @y @z 1 @R 1 @R 1 @R þ þ þ 132e2 nx @x ny @y nz @z ! 1 @nx 1 @ny 1 @nz 132e2 R 2 þ 2 þ 2 nx @x ny @y nz @z @R @R @R þ ny þnz 330e2 nx @x @y @z @nx @ny @nz 330e2 R þ þ @x @y @z
ð7Þ
where nx ¼qxfs/9rfs9,ny ¼ qyfs/9rf9 and nz ¼ @z f s =9rf s 9,Q ¼ n4x þ n4n þ n4z ,R ¼ n2x n2y n2z . To solve Eq. (7), the partial derivatives of solid fraction, such as: @x f s ,@2x f s ,@2x,y f s . . .must be accurately calculated. An actuate method to solve the partial derivatives of solid fraction for 2-D CA model is discussed in our previous work [21,22]. In the present article the method is extended into 3-D application. 2.4. Cellular automaton rule For 3-D simulation, each interface cell has 6 neighbor cells in /100S directions, 12 neighbor cells in /110S directions and 8 neighbor cells in /111S directions. So, the neighborhood rule of the cell in 3-D is much more complicated than that in 2-D. It is very difficult to present a 3-D capture rule with lower mesh induced anisotropy. So in the present work the simplest capture rule—the Von Neumann type (6 cells in /100S directions) is used, and a limited neighbor solid fraction (LNSF) method was used to control the mesh induced anisotropy. The capture rule in 3-D CA model with LNSF method is described as follows: firstly, a limited neighbor solid fraction fsLNSF was set; secondly, the average solid fraction of the neighbor cells for each liquid cell fsave was calculated. If the fsave is larger than fsLNSF, and this liquid cell has a solidified neighbor cell in the Von Neumann type, then it can be captured as an interface cell. It can be seen that the capture rule in the present 3-D CA model is a combination of the Von Neumann capture rule and the LNSF method. The Von Neumann capture rule makes sure that the CA model presents a sharp interface, and the LNSF method makes sure that all the liquid cells are captured under the similar conditions, thus the mesh induced anisotropy could be reduced.
L. Wei et al. / Physica B 407 (2012) 2471–2475
2473
Fig. 1. Simulated solidifying morphologies in the undercooled pure SCN melt, DT¼ 6 K, Dx ¼0.5 mm, the value of fsLNSF for LNSF method was set to: (a) 0.0, (b) 0.15, (c) 0.225, (d) 0.25.
Fig. 2. Simulated morphologies during solidification from the undercooled pure SCN melt, DT¼ 4.5 K, Dx ¼0.5 mm: (a) e1 ¼ 0.0, e2 ¼0.0 (b) e1 ¼0.10, e2 ¼0.0, (c) e1 ¼0.0, e2 ¼ 0.02; (d) the AMR mesh of (b) in 2D slice.
2474
L. Wei et al. / Physica B 407 (2012) 2471–2475
2.5. Computational details The equiaxed dendritic growth from undercooled pure succinonitrile (SCN) melt was simulated using 3-D CA model, which was programmed using adaptive mesh refinement (AMR) technique. The thermal physical parameters of SCN used in the simulation are the same as reference [22]. The simulation domain is 0.128 0.128 0.128 mm3, and the mesh size is smaller than 0.5 mm. At the beginning, the central cell of computational domain is initialized as a solid cell and the states for the rest cells are set to liquid. The time step is calculated as follow:
Dt ¼ 0:15
Dx2
a
ð16Þ
where Dt is the time step, Dx is mesh size. The initial temperature and boundary temperature of the domain are set to T0, thus the undercooling of the melt can be expressed as DT¼Tm T0.
3. Simulation results and discussion 3.1. The testing of mesh induced anisotropy In order to test the mesh induced anisotropy, the morphologies under the isotropic interface free energy (e1 ¼0.0, e2 ¼0.0) were first simulated and the limited neighbor solid fraction fsLNSF was set from 0.0 to 0.25. Fig. 1 shows the simulated morphologies grown from the pure SCN melt with the undercooling of 6 K. All of Fig. 1(a)–(d) were simulated by the Von Neumann capture rule, and the differences between them were the limited neighbor solid fraction fsLNSF. We found that when fsLNSF ¼0.1, the simulated result was the same as with fsLNSF ¼0.0. And when fsLNSF ¼0.3, it was too large for capture, which means no liquid cells can be captured. So we did the simulation with fsLNSF ¼0.15 and fsLNSF ¼ 0.25 for comparison. It can be seen that for fsLNSF ¼0.0, as seen in Fig. 1(a), all the 24 primary branches grew in /100S directions, which is caused by the Von Neumann capture rule. And forfsLNSF ¼0.15, the simulated morphology was still grown in /100S direction with branches, as seen in Fig. 1(b). For fsLNSF ¼0.225, the simulated morphology presents 24 primary branches with the same size evenly pointing at all the directions, as seen in Fig. 1(c). Fig. 1(d) shows the simulated morphology with /111S preferred growth directions for fsLNSF ¼0.25. It is found that when fsLNSF gets smaller, the mesh induced anisotropy prefers /100S directions. Especially, if fsLNSF is equal to 0.0, the mesh induced anisotropy is the same as that caused by the Von Neumann capture rule. However, when fsLNSF gets larger, the mesh induced anisotropy prefers /111S directions. The fsLNSF ¼0.225 presents a balance between the two mesh induced anisotropies. It can be concluded that the LNSF method with fsLNSF ¼0.225 could reduce the mesh induced anisotropy. So in the following simulations, such as Figs. 2 and 3, we used the LNSF method with fsLNSF ¼0.225.
Fig. 3. Comparison of morphologies predicted by CA model and 1/S plot for different e1 with the e2 of 0.02 for the growth from undercooled pure SCN melt, DT¼ 2.0 K, Dx¼ 0.25 mm.
(both of them are simulated with isotropic interface free energy), it can be seen that as the undercooling changes, the simulated morphology changes as well. From Fig. 2(b) and (c), It can be seen that for interface anisotropy parameters e1 ¼0.10 and e2 ¼0.0, a dendrite grew with the /100S preferred direction, and for e1 ¼0.0 and e2 ¼ 0.02, a dendrite grew with the /110S preferred direction. Fig. 2 (d) show the dendrite in 2-D slice of Fig. 2(b) illustrated by AMR mesh. It can be seen that the solid–liquid interface is covered by the most refined meshes, which can describe the interface more accurately.
3.3. Dendrite orientation selection The simulated morphologies by present 3-D CA model are further compared with the prediction of the spherical plot of the inverse of the interfacial stiffness, as expressed in Eq. (3), for anisotropy parameter e1 changes from 0.0 to 0.16 with the fixed e2 of 0.02. From Fig. 3, it can be seen that the simulated morphologies by CA model are continuously changed with the e1 as expected from the minimum stiffness criterion.
4. Summary 3.2. Dendrite morphology by 3-D CA model Fig. 2 shows the simulated morphologies grown from the pure SCN melt under the undercooling of 4.5 K. The interface anisotropy parameters are set to e1 ¼0.0, e2 ¼0.0 for Fig. 2(a), e1 ¼0.10, e2 ¼0.0 for Fig. 2(b), and e1 ¼0.0, e2 ¼ 0.02 for Fig. 2(c), respectively. When e1 ¼0.0, e2 ¼ 0.0, as seen in Fig. 2(a), the simulated morphology is a kind of cube with instability in directions of /100S, /110S and /111S. Compared Fig. 2 (a) with Fig. 1(c)
A 3-D cellular automaton model for describing the solidification of a pure substance was developed, in which a limited neighbor solid fraction (LNSF) method was presented to reduce the mesh induced anisotropy. In order to simulate the orientation selection, an expansion description using two interface free energy anisotropy parameters (e1, e2) is used in the present 3-D CA model. When e1 ¼0.10, e2 ¼0.0, the dendrite grew with the /100S preferred directions;
L. Wei et al. / Physica B 407 (2012) 2471–2475
and when e1 ¼0.0, e2 ¼ 0.02, the dendrite grew with the /110S preferred directions. The dendrite orientation selection between /100S and /110S dendrites were compared with the spherical plot of the inverse of the interfacial stiffness, using the different e1 with e2 ¼ 0.02. It is found that the simulated morphologies by present CA model are as expected from the minimum stiffness criterion.
Acknowledgment This project was supported by the National Natural Science Foundation of China (Nos. 50971102 and 50901061) and the National Basic Research Program of China (No. 2011CB610402). The work was also supported by the Program of Introducing Talents of Discipline to Universities (No. 08040) and the fund of the State Key Laboratory of Solidification Processing in NWPU (Nos. 39-QZ-2009 and 02-TZ-2008). References [1] W.J. Boettinger, S.R. Coriell, A.L. Greer, A. Karma, W. Kurz, M. Rappaz, R. Trivedi, Acta Mater. 48 (2000) 43.
2475
[2] M. Asta, C. Beckermann, A. Karma, W. Kurz, R. Napolitano, M. Plapp, G. Purdy, M. Rappaz, R. Trivedi, Acta Mater. 57 (2009) 941. [3] S.C. Huang, M.E. Glicksman, Acta Metall. 29 (5) (1981) 701. [4] J.S. Langer, Rev. Mod. Phys. 52 (1) (1980) 1. [5] D.A. Kessler, J. Koplik, H. Levine, Phys. Rev. A 31 (1985) 1712. [6] J.S. Langer, Phys. Rev. A 33 (1986) 435. [7] W.D. Huang, Q.M. Wei, Y.H. Zhou, J. Cryst. Growth 100 (1–2) (1990) 26. [8] X. Lin, W.D. Huang, J. Feng, T. Li, Y.H. Zhou, Acta Mater. 47 (11) (1999) 3271. [9] J.J. Hoyt, M. Asta, A. Karma, Mater. Sci. Eng. R 41 (2003) 121. [10] T. Haxhimali, A. Karma, F. Gonzales, M. Rappaz, Nat. Mater. 5 (2006) 660. [11] L. Nastac, Acta Mater. 47 (17) (1999) 4253–4267. [12] L. Beltran-Sanchez, D.M. Stefanescu, Metall. Mater. Trans. A 34 (2) (2003) 367. [13] M.F. Zhu, D.M. Stefanescu, Acta Mater. 55 (2007) 1741. [14] W. Wang, P.D. Lee, M. McLean, Acta Mater. 51 (10) (2003) 2971. [15] H.B. Dong, P.D. Lee, Acta Mater. 53 (3) (2005) 659. [16] M.F. Zhu, S.Y. Lee, C.P. Hong, Phys. Rev. E 69 (6) (2004) 061610. [17] S.Y. Pan, M.F. Zhu, Acta Mater. 58 (2010) 340. [18] H.B. Yin, S.D. Felicelli, Modell. Simul. Mater. Sci. Eng. 17 (2009) 075011. [19] B.W. Shan, X. Lin, L. Wei, W.D. Huang, Adv. Mater. Res. 26–28 (2007) 957. [20] L. Wei, X. Lin, B. Shan, W. Huang, Proceedings of the 5th Decennial International Conference on Solidification Processing Sheffield, U.K (2007) 53. [21] X. Lin, L. Wei, M. Wang, W.D. Huang, Mater. Sci. Forum PRICM7 (2010) 1528. [22] L. Wei, X. Lin, M. Wang, W.D. Huang, Appl. Phys. A 103 (2011) 123. [23] L. Wei, X. Lin, M. Wang, W.D. Huang, Comput. Mater. Sci. 54 (2012) 66. [24] L. Wei, X. Lin, M. Wang, W.D. Huang, arXiv preprint (2011) arXiv:1111.6159. [25] D.W. Hoffman, J.W. Cahn, Surf. Sci. 31 (1972) 368. [26] J.E. Taylor, Acta Metall. Mater. 40 (1992) 1475.