- Email: [email protected]
Quantitative phase field modeling of diffusion-controlled precipitate growth and dissolution in Ti–Al–V
Scripta Materialia 50 (2004) 471–476 www.actamat-journals.com
Quantitative phase field modeling of diffusion-controlled precipitate growth and dissolution in Ti–Al–V Qing Chen, Ning Ma, Kaisheng Wu, Yunzhi Wang
*
Department of Materials Science and Engineering, The Ohio State University, 2041 College Road, Columbus, OH 43221, USA Received 26 June 2003; received in revised form 23 October 2003; accepted 29 October 2003
Abstract A method for quantitative phase field modeling of diffusion-controlled phase transformation in multicomponent systems is demonstrated. With inputs from CALPHAD thermodynamic and DICTRA kinetic databases, the growth and dissolution of a precipitates in Ti–Al–V is simulated on experimentally relevant length and time scales. The results agree well with DICTRA simulations. 2003 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved. Keywords: Phase field; CALPHAD; Phase transformations; Kinetics; Ti–Al–V
1. Introduction Diffusion-controlled phase transformations in multicomponent systems have traditionally been simulated by using sharp interface approaches and the simulations are usually limited to 1D because of the necessity of front tracking and determination of local-equilibrium tie-lines at the interfaces [1]. In contrast, the phase field method (PFM) [2–4], also known as the diffuse interface approach, treats a multiphase microstructure as a whole without boundary tracking by using continuous field variables. It can automatically account for local equilibrium and the Gibbs–Thompson effect. Moreover, the coherency elastic strain effect can be included in a straightforward way by using Khachaturyan’s elasticity theory [5]. Because of these conveniences the PFM has become a method of choice to study microstructure evolutions due to but not limited to phase transformations in 2D and 3D. Recently, there has been an increasing interest in using this method to treat real binary and multicomponent alloys by linking bulk free energies in the phase field models to critically assessed thermodynamic databases [6–10]. However, there have been few attempts to fully incorporate both thermodynamic and mobility databases into the PFM and
*
Corresponding author. E-mail address: [email protected] (Y. Wang).
verify the kinetic results against that of the established sharp interface models, which is apparently a necessary step towards quantitative simulation of microstructural evolution. Another important issue that has not been addressed properly so far is how to break the inherent length scale limit in a quantitative phase field modeling where material specific free energy and interfacial energy data are used [11,12]. In this article, we intend to develop a multicomponent phase field model that makes direct use of assessed thermodynamic and mobility databases and simulate phase transformations on real length and time scales. Examples of applications are given for the b $ a transformation in the Ti–Al–V system. Simple geometries will be used to facilitate quantitative comparison of PFM predictions with that from DICTRA [1], a commercial software based on the sharp interface approach. With this crucial validation work done, complex geometries and real microstructure evolution will be considered in a forthcoming publication.
2. Phase field model To describe phase transformations in an n-component system using the PFM, we need n 1 concentration fields and a set of order parameter fields. The order parameters characterize symmetry changes accompanying the phase
1359-6462/$ - see front matter 2003 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved. doi:10.1016/j.scriptamat.2003.10.032
472
Q. Chen et al. / Scripta Materialia 50 (2004) 471–476
transformations and their choice can be either physical or phenomenological. For an order–disorder transformation, the long-range order (lro) parameters are the default order parameters [13] and the local free energy as a function of concentration and lro parameters can be obtained directly by the CALPHAD technique [9,14]. For a reconstructive phase transformation such as bcc (b) to hcp (a) in Ti, a Landau free energy expansion with respect to appropriate physical order parameters can be constructed according to the symmetry changes during the phase transformation [15]. The same approach can be applied to alloys, but then it becomes evident that the parameters in the Landau free energy must be made temperature and composition dependent. In order to have a Landau free energy consistent with the experimental or assessed equilibrium free energy data in a multicomponent system, we thus have to face a formidable task to fit the expression in a multidimensional space at different temperatures. An alternative approach is to define a phenomenological order parameter that assumes certain different values for phases of different symmetries. The local free energy is then constructed in such a way that the equilibrium free energy of individual phases can be directly inserted into the expression, and the equilibrium phase relationship in the temperature–composition projection can be sustained in the temperature–composition–order parameter space. A convenient choice of such an expression is due to Wang et al. [16] and has been used widely in solidification modeling [17]. Adopting this choice, we write the local molar Gibbs free energy Gm as a function of temperature T , composition Xi (i ¼ 1; 2; . . . ; n 1), and order parameter g: Gm ðT ; Xi ; gÞ ¼ ½1 pðgÞGam ðT ; Xi Þ þ pðgÞGbm ðT ; Xi Þ þ qðgÞ ð1Þ where pðgÞ ¼ g3 ð10 15g þ 6g2 Þ and qðgÞ ¼ xg2 ð1 2 gÞ . The parameter x is the height of the imposed double-well hump, which, along with the gradient energy coefficients ji and jg shown below in Eq. (2), can be determined from interfacial energy, c, and interface thickness, k. Gam and Gbm are the molar Gibbs free energies of the a and b phases, respectively. For a chemically and structurally non-uniform system under the assumption of constant molar volume Vm , the total Gibbs free energy G can be expressed by # Z " n1 X 1 ji jg 2 2 G¼ jrXi j þ jrgj dV Gm ðT ; Xi; gÞ þ Vm V 2 2 i¼1 ð2Þ where ji and jg are the gradient-energy coefficients for concentration and order parameter inhomogeneities, respectively [18,19]. The temporal evolutions of the field
variables are governed by the time-dependent Ginzburg–Landau equations [20] and the generalized Cahn– Hilliard diffusion equations [21] on the basis of the phenomenological Fick–Onsager equations [22]: og dG ¼ Mg ot dg
ð3Þ
n1 X 1 oXk dG ¼ r Mkj ðT ; Xi ; gÞr Vm2 ot dX j j¼1
ð4Þ
where Mg is the mobility of the order parameter and can be directly related to the interface mobility in the sharp interface approach. The parameters Mkj are the so-called chemical mobilities in the volume-fixed frame of reference. In a single phase p (p ¼ a; b), according gren [23], the chemical mobilities to Andersson and A p Mkj are related to atomic mobilities Mlp (l ¼ 1; . . . ; n) by Mkjp ¼
n 1 X ðdjl Xj Þðdlk Xk ÞXl Mlp Vm l¼1
ð5Þ
where djl and dlk are the Kronecker delta and the composition dependence of Mlp can be modeled in a CALPHAD type fashion [1]. In a structurally and compositionally non-uniform system, the same relation should hold between Mkj and Ml , the latter is the atomic mobilities assumed to be dependent on the order parameter g by g
Ml ¼ Mla þ Mlb ðMla Þ ðMlb Þ
ð1gÞ
ð6Þ
The choice of Eq. (6) ensures that the atomic mobilities in the interface region will have a positive deviation from the simple linear interpolation. Inserting Eq. (2) into Eqs. (3) and (4), we can obtain the following dimensionless governing equations ~ og e 2 g oG m e g j~g r ¼M os og n1 X oXk e e e kj r ¼r M os j¼1
!
~m oG e 2 Xj j~j r oXj
ð7Þ ! ð8Þ
e¼ by introducing the following reduced quantities: r ~ e ½o=oðx=lÞ; o=oðy=lÞ; Gm ¼ Gm =DGm ; M ki ¼ Vm Mki =M; j~g ¼ jg =ðDGm l2 Þ; j~i ¼ ji =ðDGm l2 Þ; s ¼ ðMDGm =l2 Þt; e g ¼ Mg l2 =ðMVm Þ, where l is the mesh size, DGm and M M are normalization factors for molar Gibbs free energy and atomic mobility, respectively. This dimensionless version of the governing equations is particularly convenient for numerical calculations and is very useful in rescaling the space and time for diffusion-controlled phase transformations.
Q. Chen et al. / Scripta Materialia 50 (2004) 471–476
3. Length and time scales The length scale of a quantitative phase field modeling in a uniform mesh scheme is inherently limited to a few or a few tenth of micrometers because the physical thickness of interfaces is in the order of nanometers or ngstroms. In order to treat microstructures in even A tens or hundreds of microns, we have to make the interface more diffuse by adjusting certain model parameters and at the same time keeping fixed the driving forces of the process. As discussed in [11], for precipitate growth or dissolution processes, if the Gibbs–Thompson effect can be ignored, we can simply increase the gradient energy coefficients by n2 times, i.e. j0g ¼ n2 jg and j0i ¼ n2 ji , and get a larger interface thickness k0 ¼ nk. Using the same number of mesh points to discretize the interface, we can now use a larger mesh size l0 ¼ nl and thus treat a larger system. According to Allen and Cahn [19], the velocity of a moving g profile is proportional to jg Mg . To keep the same velocity for a more diffuse interface, we should then use a scaled mobility of the order parameter, i.e., Mg0 ¼ Mg =n2 . Using the same procedure as that given in the preceding section to nondimensionalize the governing equations, we get all the reduced quantities the same e 0 ¼ ½o=oðx=nlÞ; o=oðy=nlÞ and s0 ¼ as before except r 2 2 ðMDGm =n l Þt. This result indicates that in order to treat a larger system size, we just need to rescale the length and time obtained for the system size dictated by actual interface thickness. It should be emphasized that by doing so the interfacial energy has been increased by n times. Obviously, this scheme is useful to treat planar interface problems where the Gibbs–Thompson effect is not a factor in the growth or dissolution kinetics. For spherical particles, however, this scheme is limited to the range of particle sizes over which the Gibbs–Thompson effect is negligible. Scaling up from particle sizes below this range will retain the Gibbs–Thompson effect acting on the growth or dissolution of the particles because both interfacial energy and particle size have been increased by the same amount, i.e., n times. In this case, we should first use the scheme for concurrent growth and coarsening processes [12] to increase the length scale to the extent where this effect is negligible by increasing the gradient energy coefficients and simultaneously decreasing the double well hump so that the interface is getting more diffuse without altering the interfacial energy. Based on this newly achieved length scale, the simple scaling scheme can then be used readily.
4. Application to b $ a transformation in Ti–Al–V The Gibbs free energies as a function of temperature and composition were adopted for the b and a phases in
473
the Ti–Al–V system from a Ti-base thermodynamic database developed by CompuTherm [24] using the CALPHAD technique. Based on this thermodynamic information and with the help of DICTRA, a set of selfconsistent parameters describing the atomic mobility of Ti, Al, and V in the two phases were obtained by assessing the experimental diffusivity in the ternary and its constituent binary systems [25]. After directly incorporating these thermodynamic and kinetic data into Eqs. (7) and (8), model parameters related to interface properties, i.e. x, ji and jg , were fitted to the assumed interfacial energy c ¼ 0:5 J/m2 (typical value for incoherent phase interfaces) and interface thickness k ¼ 5 109 m (similar in order of magnitude to grain boundary thickness) by performing a one-dimensional phase field dynamical relaxation of a system with sharp interface and equilibrium composition. In this fitting procedure, a grid size of 109 m has been chosen, which means that 5 grid points fall into the interfacial region. For simplification, the gradient energy coefficients for composition, ji , were set to zero as this will not affect the kinetic results. Under these conditions, we obtained x ¼ 30 kJ/mol and jg ¼ 6 1014 J m2 /mol. Taking l ¼ 109 m and Vm ¼ 105 m3 / mol as well as the normalizing quantities DGm ¼ 50 kJ/ mol and M ¼ 1018 mol m2 /s J, we have j~g ¼ 1:2. Finally e g ¼ 6 to warrantee a diffusion-controlled we chose M process. 4.1. Thickening of a plate We consider an a precipitate growing with a planar interface into supersaturated b at 1173 K. The initial thickness of the a plate was chosen to be 0.2 lm, and its composition was set to the equilibrium value: 11.3 at.%Al, 1.575 at.%V. The initial composition of b was 10.19 at.%Al and 3.6 at.%V. The total system size was chosen as 10 lm. The phase field simulation was performed with 500 grid points and a grid size of 109 m in one dimension. In order to describe the desired system size, the length scale of the phase field modeling must be increased by 20 times, which means the actual time must be scaled up by 400 (see Section 3). Sharp interface simulation on real length and time scales has also been carried out by using DICTRA. The two results are compared in Fig. 1. It is clear that the phase field simulation results are in good agreement with that of DICTRA simulation. As expected and clearly shown in Fig. 1a, the thickening of the plate follows the parabolic law initially. The gradual slow down of the thickening kinetics in the later stage is due to the soft impingement (see Fig. 1b and c). It is worth mentioning that the local equilibrium comes out automatically in the phase field method without explicit calculation provided that Mg and Ml within the interface are high enough. Because the
474
Q. Chen et al. / Scripta Materialia 50 (2004) 471–476
Fig. 1. PFM and DICTRA results for the growth of a precipitate plate: (a) growth kinetics; (b) composition profile of Al; and (c) composition profile of V.
diffusivities of Al and V are different in both phases, the local equilibrium tie line determined by the actual flux balance is different from the equilibrium one. As soft impingement develops, the local equilibrium tie line moves toward the equilibrium one. This is more obvious if the two composition profiles are superimposed in the isothermal phase diagram of the system at 1173 K. As can be seen from Fig. 2, first the system finds a tie line (the one for 0.1 and 0.4 ks) below the equilibrium one (dotted line marked by the triangle symbol), and then a series of tie lines approaching and slightly passing the equilibrium one. If annealing time is long enough and homogeneity is reached within each phase, the equilibrium tie line will be assumed eventually. The system described above corresponds to the growth of grain boundary a precipitates except that in reality sideplates will form shortly either by interface instability [26,27] or sympathetic nucleation mechanism [28] and inhibit the further growth of the layer of grain boundary a. However, the predicted initial thickening kinetics should be applicable if no grain boundary diffusion enhancement, i.e. the ‘‘collector plate’’ mecha-
Fig. 2. Diffusion paths corresponding to Fig. 1b and c during the growth of a precipitate plate.
nism for Al and ‘‘rejector plate’’ mechanism for V [29], is involved. 4.2. Dissolution of globular a We examine now the dissolution of a globular a that was originally at equilibrium with the b matrix at 1173 K and is raised instantly to 1223 K. For the phase field simulation, a mesh of 500 · 500 and a quarter of a circle situated at the lower left corner were used with the zero-flux Neumann boundary conditions along both dimensions. In this case, if we need to account for an actual system of 100 · 100 lm2 , we cannot, as pointed out in Section 3, apply the mesh size l ¼ 109 m and then scale it up 200 times because the significant Gibbs– Thompson effect present on the nanometer scale will be retained due to the artificial increase of interfacial energy during the simple scaling and hence alter the dissolution kinetics on the length scale of interest. Instead, we should first diffuse the interface as much as possible and at the same time keep the interfacial energy unchanged. This can be done by increasing the gradient energy coefficient and decreasing the double-well hump simultaneously. For the Ti–Al–V system, the mesh size can be increased by 50 times (l0 ¼ 5 108 m) while maintaining the same interfacial energy by setting x0 ¼ 100 J/mol and j~0g ¼ e 0 ¼ 300. After performing 0:024. Accordingly, we have M g simulation on this base system, a simple straight scaling (see Section 3) on length (by 4 times) and time (by 16 times) was carried out to match the actual size 100 · 100 lm2 . The results are compared with DICTRA simulation in Fig. 3. Apparently, an excellent agreement between the two has been obtained. It is interesting to note that the dissolution process in this case is not simply a reverse of growth process, and the dissolution kinetics does not follow the parabolic law, see Fig. 3a, which confirms Aaron and Kotler’s analysis [30]. Due to the finite system size, soft impingement occurs after 5 ks and gradually slows down the dissolution process. The Al composition spikes in the a phase region obtained from the PFM is less sharp
Q. Chen et al. / Scripta Materialia 50 (2004) 471–476
475
Fig. 3. PFM and DICTRA results for the dissolution of globular a: (a) dissolution kinetics; (b) composition profile of Al; and (c) composition profile of V.
transformations, in which both thermodynamic and kinetic data from existing databases can be inserted directly into the phase field model and the length and time of the dynamic system can be readily scaled. Applications to phase transformation between a and b in Ti–Al– V in simple geometries have been proved very successful by comparing results with that of DICTRA simulations. The power of the phase field method relies on its capability to handle arbitrarily complex geometries and this shall be demonstrated in a forthcoming publication. Some assumptions have been made in determining model parameters in the method in order to have a volume diffusion controlled process. Different assumptions may have to be evoked if solute drag, solute trapping or interface-controlled process is considered. Fig. 4. Diffusion path at 1 ks during the dissolution of globular a.
than that from DICTRA. This is understandable from the point of view of gradient thermodynamics. The local equilibrium at the interface is obtained automatically in the phase field method. Again, the system found a different tie line from the one across the original composition at 1223 K. Fig. 4 depicts the diffusion path at 1 ks in the Ti–Al–V Gibbs triangle where two sets of tie lines at 1173 and 1223 K are drawn. The original compositions of a and b are located at the ends of the tie line across the overall composition (denoted by the triangle symbol) at 1173 K, respectively. The composition in the b phase changes gradually by following a curved path towards the interface where the local equilibrium prevails at a 1223 K tie line far above the one determined by the initial composition of the a phase. The horn-shaped path forms in the a phase because the spike in the Al profile is much sharper than that for V, which has a relatively larger mobility at 1223 K.
5. Summary We have devised a scheme for quantitative phase field modeling of multicomponent diffusion-controlled phase
Acknowledgements This work is supported by AFRL under MAI grant P035558 (QC, NM and YW) and NSF under grant DMR 0139705 (KW) and Focused Research Group grant DMR-0080766 (YW). We thank Wei Wang and Sandy Ye for their valuable help with the MAI project.
References gren J. J Phase [1] Borgenstam A, Engstr€ om A, H€ oglund L, A Equilibria 2000;21:269. [2] Wang Y, Chen LQ. Simulation of microstructural evolution using the phase field method. In: Methods in materials research, a current protocols. John Wiley & Sons, Inc.; 2000. [3] Chen LQ. Ann Rev Mater Res 2002;32:113. [4] Boettinger WJ, Warren JA, Beckermann C, Karma A. Ann Rev Mater Res 2002;32:163. [5] Khachaturyan AG. Theory of structural transformations in solids. New York: John Wiley & Sons; 1983. [6] Grafe U, Botteger B, Tiaden J, Fries SG. Scripta Mater 2000;42:1179. [7] Cha P-R, Yeon D-H, Yoon J-K. Acta Mater 2001;49:3295. [8] Yeon D-H, Cha P-R, Yoon J-K. Scripta Mater 2000;45:661. [9] Zhu JZ, Liu ZK, Vaithyanathan V, Chen LQ. Scripta Mater 2002;46:401.
476
Q. Chen et al. / Scripta Materialia 50 (2004) 471–476
gren J. Acta Mater [10] Loginova I, Odqvist J, Amberg G, A 2003;51:1327. [11] Shen C, Chen Q, Wang Y, Wen YH, Simmons JP. Scripta Mater, submitted for publication. [12] Shen C, Chen Q, Wang Y, Wen YH, Simmons JP. Scripta Mater, submitted for publication. [13] Wang Y, Banerjee D, Su CC, Khachaturyan AG. Acta Mater 1998;46:2983. [14] Ansara I, Dupin N, Lukas HL, Sundman B. J Alloy Compd 1997;247:20. [15] Toledano P, Dimitriev V. Reconstructive phase transitions: in crystals and quasicrystals. Singapore: World Scientific; 1996. [16] Wang S-L, Sekerka RF, Wheeler AA, Murray BT, Coriell SR, Braun RJ, et al. Physica D 1993;69:189. [17] Warren JA, Boettinger WJ. Acta Metall Mater 1995;43:689. [18] Cahn JW, Hilliard JE. J Chem Phys 1958;28:258. [19] Allen SM, Cahn JW. Acta Metall 1979;27:1085.
[20] Gunton JD, Miguel MS, Sahni PS. The dynamics of first-order phase transitions. In: Domb C, Lebowitz JL, editors. Phase transitions and critical phenomena, vol. 8. New York: Academic Press; 1983. [21] Cahn JW. Acta Metall 1961;9:795. [22] Kirkaldy JS, Young DJ. Diffusion in the condensed state. London: Institute of Metals; 1987. gren J. J Appl Phys 1992;72:1350. [23] Andersson J-O, A [24] Zhang F. CompuTherm, LLC, USA, 2002, private communication. [25] Chen Q, Wang Y. in preparation (data available upon request). [26] Mullins WW, Sekerka RF. J Appl Phys 1963;34:323. [27] Townsend RD, Kirkaldy JS. Trans ASM 1968;61:605. [28] Aaronson HI, Spanos G, Masamura RA, Vardiman RG, Moon DW, Menon ESK, et al. Mater Sci Eng, B 1995;32:107. [29] Menon ESK, Aaronson HI. Metall Trans A 1986;17A:1703. [30] Aaron HB, Kotler GR. Metall Trans 1971;2:393.
Quantitative phase field modeling of diffusion-controlled precipitate growth and dissolution in Ti–Al–V Qing Chen, Ning Ma, Kaisheng Wu, Yunzhi Wang
*
Department of Materials Science and Engineering, The Ohio State University, 2041 College Road, Columbus, OH 43221, USA Received 26 June 2003; received in revised form 23 October 2003; accepted 29 October 2003
Abstract A method for quantitative phase field modeling of diffusion-controlled phase transformation in multicomponent systems is demonstrated. With inputs from CALPHAD thermodynamic and DICTRA kinetic databases, the growth and dissolution of a precipitates in Ti–Al–V is simulated on experimentally relevant length and time scales. The results agree well with DICTRA simulations. 2003 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved. Keywords: Phase field; CALPHAD; Phase transformations; Kinetics; Ti–Al–V
1. Introduction Diffusion-controlled phase transformations in multicomponent systems have traditionally been simulated by using sharp interface approaches and the simulations are usually limited to 1D because of the necessity of front tracking and determination of local-equilibrium tie-lines at the interfaces [1]. In contrast, the phase field method (PFM) [2–4], also known as the diffuse interface approach, treats a multiphase microstructure as a whole without boundary tracking by using continuous field variables. It can automatically account for local equilibrium and the Gibbs–Thompson effect. Moreover, the coherency elastic strain effect can be included in a straightforward way by using Khachaturyan’s elasticity theory [5]. Because of these conveniences the PFM has become a method of choice to study microstructure evolutions due to but not limited to phase transformations in 2D and 3D. Recently, there has been an increasing interest in using this method to treat real binary and multicomponent alloys by linking bulk free energies in the phase field models to critically assessed thermodynamic databases [6–10]. However, there have been few attempts to fully incorporate both thermodynamic and mobility databases into the PFM and
*
Corresponding author. E-mail address: [email protected] (Y. Wang).
verify the kinetic results against that of the established sharp interface models, which is apparently a necessary step towards quantitative simulation of microstructural evolution. Another important issue that has not been addressed properly so far is how to break the inherent length scale limit in a quantitative phase field modeling where material specific free energy and interfacial energy data are used [11,12]. In this article, we intend to develop a multicomponent phase field model that makes direct use of assessed thermodynamic and mobility databases and simulate phase transformations on real length and time scales. Examples of applications are given for the b $ a transformation in the Ti–Al–V system. Simple geometries will be used to facilitate quantitative comparison of PFM predictions with that from DICTRA [1], a commercial software based on the sharp interface approach. With this crucial validation work done, complex geometries and real microstructure evolution will be considered in a forthcoming publication.
2. Phase field model To describe phase transformations in an n-component system using the PFM, we need n 1 concentration fields and a set of order parameter fields. The order parameters characterize symmetry changes accompanying the phase
1359-6462/$ - see front matter 2003 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved. doi:10.1016/j.scriptamat.2003.10.032
472
Q. Chen et al. / Scripta Materialia 50 (2004) 471–476
transformations and their choice can be either physical or phenomenological. For an order–disorder transformation, the long-range order (lro) parameters are the default order parameters [13] and the local free energy as a function of concentration and lro parameters can be obtained directly by the CALPHAD technique [9,14]. For a reconstructive phase transformation such as bcc (b) to hcp (a) in Ti, a Landau free energy expansion with respect to appropriate physical order parameters can be constructed according to the symmetry changes during the phase transformation [15]. The same approach can be applied to alloys, but then it becomes evident that the parameters in the Landau free energy must be made temperature and composition dependent. In order to have a Landau free energy consistent with the experimental or assessed equilibrium free energy data in a multicomponent system, we thus have to face a formidable task to fit the expression in a multidimensional space at different temperatures. An alternative approach is to define a phenomenological order parameter that assumes certain different values for phases of different symmetries. The local free energy is then constructed in such a way that the equilibrium free energy of individual phases can be directly inserted into the expression, and the equilibrium phase relationship in the temperature–composition projection can be sustained in the temperature–composition–order parameter space. A convenient choice of such an expression is due to Wang et al. [16] and has been used widely in solidification modeling [17]. Adopting this choice, we write the local molar Gibbs free energy Gm as a function of temperature T , composition Xi (i ¼ 1; 2; . . . ; n 1), and order parameter g: Gm ðT ; Xi ; gÞ ¼ ½1 pðgÞGam ðT ; Xi Þ þ pðgÞGbm ðT ; Xi Þ þ qðgÞ ð1Þ where pðgÞ ¼ g3 ð10 15g þ 6g2 Þ and qðgÞ ¼ xg2 ð1 2 gÞ . The parameter x is the height of the imposed double-well hump, which, along with the gradient energy coefficients ji and jg shown below in Eq. (2), can be determined from interfacial energy, c, and interface thickness, k. Gam and Gbm are the molar Gibbs free energies of the a and b phases, respectively. For a chemically and structurally non-uniform system under the assumption of constant molar volume Vm , the total Gibbs free energy G can be expressed by # Z " n1 X 1 ji jg 2 2 G¼ jrXi j þ jrgj dV Gm ðT ; Xi; gÞ þ Vm V 2 2 i¼1 ð2Þ where ji and jg are the gradient-energy coefficients for concentration and order parameter inhomogeneities, respectively [18,19]. The temporal evolutions of the field
variables are governed by the time-dependent Ginzburg–Landau equations [20] and the generalized Cahn– Hilliard diffusion equations [21] on the basis of the phenomenological Fick–Onsager equations [22]: og dG ¼ Mg ot dg
ð3Þ
n1 X 1 oXk dG ¼ r Mkj ðT ; Xi ; gÞr Vm2 ot dX j j¼1
ð4Þ
where Mg is the mobility of the order parameter and can be directly related to the interface mobility in the sharp interface approach. The parameters Mkj are the so-called chemical mobilities in the volume-fixed frame of reference. In a single phase p (p ¼ a; b), according gren [23], the chemical mobilities to Andersson and A p Mkj are related to atomic mobilities Mlp (l ¼ 1; . . . ; n) by Mkjp ¼
n 1 X ðdjl Xj Þðdlk Xk ÞXl Mlp Vm l¼1
ð5Þ
where djl and dlk are the Kronecker delta and the composition dependence of Mlp can be modeled in a CALPHAD type fashion [1]. In a structurally and compositionally non-uniform system, the same relation should hold between Mkj and Ml , the latter is the atomic mobilities assumed to be dependent on the order parameter g by g
Ml ¼ Mla þ Mlb ðMla Þ ðMlb Þ
ð1gÞ
ð6Þ
The choice of Eq. (6) ensures that the atomic mobilities in the interface region will have a positive deviation from the simple linear interpolation. Inserting Eq. (2) into Eqs. (3) and (4), we can obtain the following dimensionless governing equations ~ og e 2 g oG m e g j~g r ¼M os og n1 X oXk e e e kj r ¼r M os j¼1
!
~m oG e 2 Xj j~j r oXj
ð7Þ ! ð8Þ
e¼ by introducing the following reduced quantities: r ~ e ½o=oðx=lÞ; o=oðy=lÞ; Gm ¼ Gm =DGm ; M ki ¼ Vm Mki =M; j~g ¼ jg =ðDGm l2 Þ; j~i ¼ ji =ðDGm l2 Þ; s ¼ ðMDGm =l2 Þt; e g ¼ Mg l2 =ðMVm Þ, where l is the mesh size, DGm and M M are normalization factors for molar Gibbs free energy and atomic mobility, respectively. This dimensionless version of the governing equations is particularly convenient for numerical calculations and is very useful in rescaling the space and time for diffusion-controlled phase transformations.
Q. Chen et al. / Scripta Materialia 50 (2004) 471–476
3. Length and time scales The length scale of a quantitative phase field modeling in a uniform mesh scheme is inherently limited to a few or a few tenth of micrometers because the physical thickness of interfaces is in the order of nanometers or ngstroms. In order to treat microstructures in even A tens or hundreds of microns, we have to make the interface more diffuse by adjusting certain model parameters and at the same time keeping fixed the driving forces of the process. As discussed in [11], for precipitate growth or dissolution processes, if the Gibbs–Thompson effect can be ignored, we can simply increase the gradient energy coefficients by n2 times, i.e. j0g ¼ n2 jg and j0i ¼ n2 ji , and get a larger interface thickness k0 ¼ nk. Using the same number of mesh points to discretize the interface, we can now use a larger mesh size l0 ¼ nl and thus treat a larger system. According to Allen and Cahn [19], the velocity of a moving g profile is proportional to jg Mg . To keep the same velocity for a more diffuse interface, we should then use a scaled mobility of the order parameter, i.e., Mg0 ¼ Mg =n2 . Using the same procedure as that given in the preceding section to nondimensionalize the governing equations, we get all the reduced quantities the same e 0 ¼ ½o=oðx=nlÞ; o=oðy=nlÞ and s0 ¼ as before except r 2 2 ðMDGm =n l Þt. This result indicates that in order to treat a larger system size, we just need to rescale the length and time obtained for the system size dictated by actual interface thickness. It should be emphasized that by doing so the interfacial energy has been increased by n times. Obviously, this scheme is useful to treat planar interface problems where the Gibbs–Thompson effect is not a factor in the growth or dissolution kinetics. For spherical particles, however, this scheme is limited to the range of particle sizes over which the Gibbs–Thompson effect is negligible. Scaling up from particle sizes below this range will retain the Gibbs–Thompson effect acting on the growth or dissolution of the particles because both interfacial energy and particle size have been increased by the same amount, i.e., n times. In this case, we should first use the scheme for concurrent growth and coarsening processes [12] to increase the length scale to the extent where this effect is negligible by increasing the gradient energy coefficients and simultaneously decreasing the double well hump so that the interface is getting more diffuse without altering the interfacial energy. Based on this newly achieved length scale, the simple scaling scheme can then be used readily.
4. Application to b $ a transformation in Ti–Al–V The Gibbs free energies as a function of temperature and composition were adopted for the b and a phases in
473
the Ti–Al–V system from a Ti-base thermodynamic database developed by CompuTherm [24] using the CALPHAD technique. Based on this thermodynamic information and with the help of DICTRA, a set of selfconsistent parameters describing the atomic mobility of Ti, Al, and V in the two phases were obtained by assessing the experimental diffusivity in the ternary and its constituent binary systems [25]. After directly incorporating these thermodynamic and kinetic data into Eqs. (7) and (8), model parameters related to interface properties, i.e. x, ji and jg , were fitted to the assumed interfacial energy c ¼ 0:5 J/m2 (typical value for incoherent phase interfaces) and interface thickness k ¼ 5 109 m (similar in order of magnitude to grain boundary thickness) by performing a one-dimensional phase field dynamical relaxation of a system with sharp interface and equilibrium composition. In this fitting procedure, a grid size of 109 m has been chosen, which means that 5 grid points fall into the interfacial region. For simplification, the gradient energy coefficients for composition, ji , were set to zero as this will not affect the kinetic results. Under these conditions, we obtained x ¼ 30 kJ/mol and jg ¼ 6 1014 J m2 /mol. Taking l ¼ 109 m and Vm ¼ 105 m3 / mol as well as the normalizing quantities DGm ¼ 50 kJ/ mol and M ¼ 1018 mol m2 /s J, we have j~g ¼ 1:2. Finally e g ¼ 6 to warrantee a diffusion-controlled we chose M process. 4.1. Thickening of a plate We consider an a precipitate growing with a planar interface into supersaturated b at 1173 K. The initial thickness of the a plate was chosen to be 0.2 lm, and its composition was set to the equilibrium value: 11.3 at.%Al, 1.575 at.%V. The initial composition of b was 10.19 at.%Al and 3.6 at.%V. The total system size was chosen as 10 lm. The phase field simulation was performed with 500 grid points and a grid size of 109 m in one dimension. In order to describe the desired system size, the length scale of the phase field modeling must be increased by 20 times, which means the actual time must be scaled up by 400 (see Section 3). Sharp interface simulation on real length and time scales has also been carried out by using DICTRA. The two results are compared in Fig. 1. It is clear that the phase field simulation results are in good agreement with that of DICTRA simulation. As expected and clearly shown in Fig. 1a, the thickening of the plate follows the parabolic law initially. The gradual slow down of the thickening kinetics in the later stage is due to the soft impingement (see Fig. 1b and c). It is worth mentioning that the local equilibrium comes out automatically in the phase field method without explicit calculation provided that Mg and Ml within the interface are high enough. Because the
474
Q. Chen et al. / Scripta Materialia 50 (2004) 471–476
Fig. 1. PFM and DICTRA results for the growth of a precipitate plate: (a) growth kinetics; (b) composition profile of Al; and (c) composition profile of V.
diffusivities of Al and V are different in both phases, the local equilibrium tie line determined by the actual flux balance is different from the equilibrium one. As soft impingement develops, the local equilibrium tie line moves toward the equilibrium one. This is more obvious if the two composition profiles are superimposed in the isothermal phase diagram of the system at 1173 K. As can be seen from Fig. 2, first the system finds a tie line (the one for 0.1 and 0.4 ks) below the equilibrium one (dotted line marked by the triangle symbol), and then a series of tie lines approaching and slightly passing the equilibrium one. If annealing time is long enough and homogeneity is reached within each phase, the equilibrium tie line will be assumed eventually. The system described above corresponds to the growth of grain boundary a precipitates except that in reality sideplates will form shortly either by interface instability [26,27] or sympathetic nucleation mechanism [28] and inhibit the further growth of the layer of grain boundary a. However, the predicted initial thickening kinetics should be applicable if no grain boundary diffusion enhancement, i.e. the ‘‘collector plate’’ mecha-
Fig. 2. Diffusion paths corresponding to Fig. 1b and c during the growth of a precipitate plate.
nism for Al and ‘‘rejector plate’’ mechanism for V [29], is involved. 4.2. Dissolution of globular a We examine now the dissolution of a globular a that was originally at equilibrium with the b matrix at 1173 K and is raised instantly to 1223 K. For the phase field simulation, a mesh of 500 · 500 and a quarter of a circle situated at the lower left corner were used with the zero-flux Neumann boundary conditions along both dimensions. In this case, if we need to account for an actual system of 100 · 100 lm2 , we cannot, as pointed out in Section 3, apply the mesh size l ¼ 109 m and then scale it up 200 times because the significant Gibbs– Thompson effect present on the nanometer scale will be retained due to the artificial increase of interfacial energy during the simple scaling and hence alter the dissolution kinetics on the length scale of interest. Instead, we should first diffuse the interface as much as possible and at the same time keep the interfacial energy unchanged. This can be done by increasing the gradient energy coefficient and decreasing the double-well hump simultaneously. For the Ti–Al–V system, the mesh size can be increased by 50 times (l0 ¼ 5 108 m) while maintaining the same interfacial energy by setting x0 ¼ 100 J/mol and j~0g ¼ e 0 ¼ 300. After performing 0:024. Accordingly, we have M g simulation on this base system, a simple straight scaling (see Section 3) on length (by 4 times) and time (by 16 times) was carried out to match the actual size 100 · 100 lm2 . The results are compared with DICTRA simulation in Fig. 3. Apparently, an excellent agreement between the two has been obtained. It is interesting to note that the dissolution process in this case is not simply a reverse of growth process, and the dissolution kinetics does not follow the parabolic law, see Fig. 3a, which confirms Aaron and Kotler’s analysis [30]. Due to the finite system size, soft impingement occurs after 5 ks and gradually slows down the dissolution process. The Al composition spikes in the a phase region obtained from the PFM is less sharp
Q. Chen et al. / Scripta Materialia 50 (2004) 471–476
475
Fig. 3. PFM and DICTRA results for the dissolution of globular a: (a) dissolution kinetics; (b) composition profile of Al; and (c) composition profile of V.
transformations, in which both thermodynamic and kinetic data from existing databases can be inserted directly into the phase field model and the length and time of the dynamic system can be readily scaled. Applications to phase transformation between a and b in Ti–Al– V in simple geometries have been proved very successful by comparing results with that of DICTRA simulations. The power of the phase field method relies on its capability to handle arbitrarily complex geometries and this shall be demonstrated in a forthcoming publication. Some assumptions have been made in determining model parameters in the method in order to have a volume diffusion controlled process. Different assumptions may have to be evoked if solute drag, solute trapping or interface-controlled process is considered. Fig. 4. Diffusion path at 1 ks during the dissolution of globular a.
than that from DICTRA. This is understandable from the point of view of gradient thermodynamics. The local equilibrium at the interface is obtained automatically in the phase field method. Again, the system found a different tie line from the one across the original composition at 1223 K. Fig. 4 depicts the diffusion path at 1 ks in the Ti–Al–V Gibbs triangle where two sets of tie lines at 1173 and 1223 K are drawn. The original compositions of a and b are located at the ends of the tie line across the overall composition (denoted by the triangle symbol) at 1173 K, respectively. The composition in the b phase changes gradually by following a curved path towards the interface where the local equilibrium prevails at a 1223 K tie line far above the one determined by the initial composition of the a phase. The horn-shaped path forms in the a phase because the spike in the Al profile is much sharper than that for V, which has a relatively larger mobility at 1223 K.
5. Summary We have devised a scheme for quantitative phase field modeling of multicomponent diffusion-controlled phase
Acknowledgements This work is supported by AFRL under MAI grant P035558 (QC, NM and YW) and NSF under grant DMR 0139705 (KW) and Focused Research Group grant DMR-0080766 (YW). We thank Wei Wang and Sandy Ye for their valuable help with the MAI project.
References gren J. J Phase [1] Borgenstam A, Engstr€ om A, H€ oglund L, A Equilibria 2000;21:269. [2] Wang Y, Chen LQ. Simulation of microstructural evolution using the phase field method. In: Methods in materials research, a current protocols. John Wiley & Sons, Inc.; 2000. [3] Chen LQ. Ann Rev Mater Res 2002;32:113. [4] Boettinger WJ, Warren JA, Beckermann C, Karma A. Ann Rev Mater Res 2002;32:163. [5] Khachaturyan AG. Theory of structural transformations in solids. New York: John Wiley & Sons; 1983. [6] Grafe U, Botteger B, Tiaden J, Fries SG. Scripta Mater 2000;42:1179. [7] Cha P-R, Yeon D-H, Yoon J-K. Acta Mater 2001;49:3295. [8] Yeon D-H, Cha P-R, Yoon J-K. Scripta Mater 2000;45:661. [9] Zhu JZ, Liu ZK, Vaithyanathan V, Chen LQ. Scripta Mater 2002;46:401.
476
Q. Chen et al. / Scripta Materialia 50 (2004) 471–476
gren J. Acta Mater [10] Loginova I, Odqvist J, Amberg G, A 2003;51:1327. [11] Shen C, Chen Q, Wang Y, Wen YH, Simmons JP. Scripta Mater, submitted for publication. [12] Shen C, Chen Q, Wang Y, Wen YH, Simmons JP. Scripta Mater, submitted for publication. [13] Wang Y, Banerjee D, Su CC, Khachaturyan AG. Acta Mater 1998;46:2983. [14] Ansara I, Dupin N, Lukas HL, Sundman B. J Alloy Compd 1997;247:20. [15] Toledano P, Dimitriev V. Reconstructive phase transitions: in crystals and quasicrystals. Singapore: World Scientific; 1996. [16] Wang S-L, Sekerka RF, Wheeler AA, Murray BT, Coriell SR, Braun RJ, et al. Physica D 1993;69:189. [17] Warren JA, Boettinger WJ. Acta Metall Mater 1995;43:689. [18] Cahn JW, Hilliard JE. J Chem Phys 1958;28:258. [19] Allen SM, Cahn JW. Acta Metall 1979;27:1085.
[20] Gunton JD, Miguel MS, Sahni PS. The dynamics of first-order phase transitions. In: Domb C, Lebowitz JL, editors. Phase transitions and critical phenomena, vol. 8. New York: Academic Press; 1983. [21] Cahn JW. Acta Metall 1961;9:795. [22] Kirkaldy JS, Young DJ. Diffusion in the condensed state. London: Institute of Metals; 1987. gren J. J Appl Phys 1992;72:1350. [23] Andersson J-O, A [24] Zhang F. CompuTherm, LLC, USA, 2002, private communication. [25] Chen Q, Wang Y. in preparation (data available upon request). [26] Mullins WW, Sekerka RF. J Appl Phys 1963;34:323. [27] Townsend RD, Kirkaldy JS. Trans ASM 1968;61:605. [28] Aaronson HI, Spanos G, Masamura RA, Vardiman RG, Moon DW, Menon ESK, et al. Mater Sci Eng, B 1995;32:107. [29] Menon ESK, Aaronson HI. Metall Trans A 1986;17A:1703. [30] Aaron HB, Kotler GR. Metall Trans 1971;2:393.