Numerical Analysis of Phase Decomposition in A-B Binary Alloys Using Cahn-Hilliard Equations

The analysis of phase decomposition was carried out using the nonlinear and linear Cahn-Hilliard equations in a hypothetical A-B alloy system with a miscibility gap. These equations were solved by the explicit finite difference method assuming a regular solution model. The supersaturated solid solution and decomposed phases were considered to have an fcc structure. Different aging temperatures and thermodynamic interaction parameters Ω A-B were used to simulate different alloy systems. The numerical simulation results showed that the growth kinetics of phase decomposition in the alloy with 30at.% A was slower than that of 50 at.% A. Additionally, the start time and modulation wavelength of phase decomposition are strongly affected by the thermodynamic interaction parameter Ω A-B value. The numerical simulation results showed that the growth kinetics of phase decomposition with the linear equation is slower than that with the nonlinear one.


Introduction
The microstructure evolution during the heat treatments of industrial alloys plays an important role in the mechanical properties of final products.Nowadays, the use of numerical methods has become a good alternative for the modelation and simulation of microstructural evolution during the heating of alloys either as a part of the heat treatment or during the service-operation of industrial components.Since these methods have permitted to simulate different heat treating or operating conditions, which are sometimes difficult to obtain in a practical way for instance, very prolonged times 1 .Additionally, the use of numerical methods can be used to understand in more detail the mechanism and growth kinetics of different phase transformations, which occur during the heating of alloys 2,3 .
Several numerical methods 4 have been used to analyze the phase transformations in alloys.One of these methods is the called phase-field method which usually is based on a solution of the nonlinear Cahn-Hilliard equation using mainly thermodynamic and atomic diffusion data 5 .The phase-field method has been used to simulate phase transformations in different alloy systems 2,3,6,7 .
One of the typical applications of the phase-field method is the microstructural simulation of the phase decomposition in alloy systems with a miscibility gap 8 .This kind of alloy systems is important for different industrial alloys such as, Cu-Ni, Al-Zn, Fe-Cr, Cu-Ni-Fe, Cu-Ni-Sn, Cu-Ni-Cr, Fe-Cr-Co, etc.Thus, it is necessary to study the effect of the different parameters of Cahn-Hilliard equations on the growth kinetics and the microstructural evolution of this type of spinodally-decomposed alloys.
Additionally, it has been pointed out 9,10 that the linear Cahn-Hilliard equation is suitable to analyze the early stage of phase decomposition during aging of these alloy systems.It is also interesting to analyze the application of the linear equation to the microstructural simulation for this kind of phase transformation and for the comparison to the nonlinear equation.It is important to mention that there is no such comparison work reported in the literature.Furthermore, most of the simulation works 11,12 reported in the literature, using the phase field method based on the nonlinear Cahn-Hilliard equation for hypothetic binary A-B alloys, use mathematical functions to model the spinodal curve in a free energy vs. composition diagram, but these ones do not represent thermodynamic parameters which can be associated with a specific alloy system.Therefore, it is necessary to use free energy functions based on thermodynamic solution models in order to approach the microstructure simulation to that taking place in real alloy systems 13 .
Thus, the goal of this work is to analyze the phase decomposition process in a binary A-B alloy system using both the linear and nonlinear Cahn-Hilliard equations, in order to understand the effect of the equation parameters on the microstructural evolution and growth kinetics for this transformation.

Cahn-Hilliard equations
The linear Cahn-Hilliard equation 9,10 is the basis for the theory of the spinodal decomposition in alloys, developed by Cahn and Hilliard.This theory has been used to analyze the phase decomposition in numerous alloys 8 .The most common mathematical expression is the following: where c is the concentration of either A or B elements as a function of a position vector and time t, M is the atomic mobility, f the local free energy, and K the gradient energy coefficient.
In contrast, the nonlinear Cahn-Hilliard equation has been the base for numerical simulation of different metallurgical phenomena 5 such as, solidification, recrystallization, phase decomposition, etc.This equation is expressed as follows: The same type of parameters of the linear equation is also involved in the nonlinear one.The linearization of equation, Equation 1, can be obtained from Equation 2, if c is assumed to be only slightly different from its average value 9 .
Both equations are partial differential equations and therefore, they can be solved using the finite difference method 4 .

Thermodynamical parameters
As stated above, one of the important parameters of the nonlinear Cahn-Hilliard equation corresponds to the local free energy f which can be defined in a simple way using the strict regular solution model 13 for a binary alloy system as follows: ( ) where R is the gas constant, T is the absolute temperature.X A and X B are the mole fractions of A and B, respectively.f A and f B are the molar free energy of pure element A and B, respectively, and Ω A-B is the interaction parameter 13 between the A and B atoms.Ω A-B >0 when A and B atoms are repulsive.On the contrary, Ω A-B <0 when A and B atoms are attractive.Figure 1 shows the plots of free energy f versus composition c for the A-B binary system at different temperatures considering the value of Ω A-B as a multiple value of gas constant R, 1000 R, 1500R and 2000R.The values of f A and f B were assumed to be equal to R (T-T A ) and R (T-T B ), respectively.T A and T B correspond to the melting point of A and B, respectively, shown in Table 1.The f vs. c curve at high temperatures is convex which indicates that the single phase is stable over all the entire composition range.At the lowest temperature, the f vs. c curve shows two minima and one maximum.These points are defined where df/dc=0 and all the minima, corresponding to each temperature, represent the miscibility gap or binodal curve shown in Figure 2. The inflection points, where the f vs. c curve changes from convex to concave downward, are defined by the second derivative df 2 /dc 2 =0.All the inflection points for each temperature form the spinodal curve, also shown in Figure 2. The spinodal curve is parabolic and it is inscribed in the binodal curve at its vertex (c=0.5).Therefore, the critical temperature 13 of the two-phase decomposition T C can be expressed by:

Diffusion parameters
The diffusion parameter of both Cahn-Hilliard equations is the atomic mobility M, which is usually defined 10 as a function of the interdiffusion coefficient D by the following expression:

Parameter Value
Interaction parameter Ω A-B [J mol -1 ] 8314, 12471, 16628   D can be defined 8 for a binary alloy system as follows: (1 ) where D A and D B are the diffusion coefficients of A and B, respectively, and they can be expressed by an Arrhenius equation as follows 8 : where D 0 is the frequency factor and it depends on the crystalline structure and Q the activation energy for atomic diffusion and it is a function of the melting point 8 .These parameters are also shown in Table 1.

Other parameters
The gradient energy coefficient K can be defined as proposed by Hilliard 9 : where H M is the heat of mixing per volume unit and r 0 is the nearest-neighbor distance.The heat of mixing H M was determined according to the following equation 10 : (1 ) The nearest-neighbor distance r o was estimated considering that the solution treated and decomposed phases had an fcc crystalline structure.The lattice parameter a is shown in Table 1.
Additionally, the elastic-strain energy due to the coherency between the decomposed phases was neglected since it was considered a similar lattice parameter of the decomposed phases in order to have a very low lattice misfit, which is related directly to the coherency elasticstain energy 10 .

Computer program
Both Can-Hilliard equations, Equations 1 and 2, were solved numerically in two dimensions using the explicit finite difference method with 101 × 101 points-square grit with a mesh size of 0.25 nm and a time-step size of 1 s.The computer program was coded in FORTRAN 95.Table 1 shows the numerical values of variables used for this simulation.The selected chemical compositions were 30 and 50 at.% of a hypothetic element A, which correspond to the asymmetrical and symmetrical positions in the miscibility gap of Figure 2. The aging temperatures and times were 400-900 K and 1000-100000000 s, respectively.

Concentration profiles
Figures 3 and 4 show the concentration profiles, the plot of the concentration of A element versus distance, for the B-50 and 30at.%Aalloys, respectively, aged at 650 K for the case of Ω A-B = 16628 J/mol (2000 R).The increase in the amplitude of the initial fluctuation with the aging time can be observed in both alloys.This fact confirms the occurrence of the phase by the spinodal decomposition mechanism in the A-B alloy system 9 .The phase separation reaction took place as follows: Supersaturated solid solution → A-rich phase + B-rich phase (10)   The initial fluctuation has to get a minimum wavelength value in order to start the spinodal decomposition in both cases.This is in agreement with the theory of spinodal decomposition which states that there is a minimum wavelength in order to increase the modulation amplitude with time 9 .The evolution of the composition modulations was faster in the case of the symmetric alloy, 50 %at.A, than that of the asymmetric alloy, 30 at.% A. According to the spinodal decomposition theory of Cahn and Hilliard 9 , the amplitude of modulation, c-c o , is given by 0 ( , ) i x c c A t e − = β β (11)   in which c o is the average composition and A(β,t) is the amplitude of the Fourier component of the wave number β (β=2π/λ) at time t expressed in terms of the initial amplitude at time t=0 as follows: Where R(β) is the amplification factor given by Inside the spinodal ∂ 2 f/∂c 2 <0 and R(β) >0 for all values of β.Thus, any modulation will grow, according to Equation 12.In order to obtain higher and faster fluctuation amplitudes, the value of the driving force for spinodal decomposition, ∂ 2 f/∂c 2 , must be also high, according to Equations 11 and 13.The variation of the second derivative with composition is shown in Figure 5 for the case of Ω A-B = 16628 J/mol at T= 650 K.This figure shows that the symmetric alloy has a higher value of the driving force, ∂ 2 f/∂c 2 , than that of the asymmetric one and thus the kinetic behavior of phase decomposition kinetics is faster in the former alloy.
Figures 6 and 7 show the concentration profiles for the A-50at.%Balloy aged at 650 K for the case of Ω A-B = 6628 J/mol using the nonlinear and linear Cahn-Hilliard equations, respectively.The increase in the amplitude of the initial fluctuation with aging time can be also observed in both cases.This increase confirms the occurrence of the phase separation by the spinodal decomposition mechanism in the A-B alloy system 9 .
It is interesting to notice that the composition modulation is more sinusoidal and periodic in the case of the linear equation.In contrast, the coarsening process of the composition fluctuations is more evident with the nonlinear equation than in the linear one.
The initial fluctuation has to get a minimum wavelength value in order to start the spinodal decomposition in both cases.The simulation with the nonlinear equation requires a smaller minimum value of wavelength than that with the linear one.Besides, the modulation amplitude reached with the nonlinear equation is higher than that corresponding to the linear one.According to the spinodal decomposition theory of Cahn and Hilliard 9 , the minimum wavelength λ min for spinodal decomposition is given by: In order to obtain a lower value of minimum wavelength and higher fluctuation amplitude, the value of driving force for spinodal decomposition, ∂ 2 f/∂c 2 , must be also high, according to Equation 14.This fact suggests that the nonlinear equation involves a higher driving force than that of the linear equation.That is, the driving force is represented by the first derivative of the free energy with respect to the composition, ∂f/∂c, in the nonlinear equation, Equation 2, and by the second derivative, ∂ 2 f/∂c 2 in the linear equation, Equation 1.The variation of the first derivative and second derivative are shown in Figures 8 and 5, respectively, for the case of Ω A-B = 16628 J/mol at T= 650 K.In general, the change in the first derivative with composition is higher between the inflexion points, c= 0.2 and 0.8, than that of the second derivative.This fact means that the driving force for the nonlinear equation is higher than that in the linear one.for Ω A-B = 16628 J/mol.case of Ω A-B =16628 J/mol at 650 K.The black and gray zones correspond to the A-rich and B-rich phases, respectively.The microstructure was determined using the two-dimension concentration results.The morphology of decomposed phases is irregular an interconnected for both alloys.This type of morphology is a characteristic at the early stages of the spinodal decomposition 9 .The interconnected microstructure is more clearly observed in the case of the alloy with 50 at.%Athan at that corresponding to the other alloy composition.This fact can also be attributed to the higher driving force for the phase decomposition in the former alloy.

Microstructural evolution
Additionally, the morphology of decomposed phases was irregular an interconnected for both equations in the case of Ω A-B =16628 J/mol at 650 K, as expected in the early stages of spinodal decomposition 9 .Nevertheless, the size of decomposed phases predicted by the nonlinear equation is much smaller than that of the linear one.

Decomposition kinetics
To determine the modulation wavelength of the composition fluctuations, the concentration profiles of the present work were analyzed using the autocorrelation analysis 12 .Figure 11 shows the growth kinetics of spinodal decomposition, plot of modulation wavelength λ versus aging time t, for the B-30 and 50 at.%A alloys for Ω A-B =16628 J/mol at 650 K.
In both cases, the modulation wavelength remains almost constant with time, which is also a characteristic of the spinodal decomposition process 9 .The start of phase decomposition occurred first for the symmetric alloy than for the asymmetric one.This can also be attributed to the higher driving force decomposition of the former alloy.
According to the spinodal decomposition theory of Cahn and Hilliard 9 , the minimum wavelength λ min for spinodal decomposition is given by Equation 14.It is evident that the modulation wavelength λ min , is shorter if the alloy composition has either a higher driving force ∂ 2 f/∂c 2 or a lower gradient energy coefficient K.
The modulation wavelengths λ min , shown in Figure 11, are very close for the symmetric and asymmetric alloy compositions.This can be associated with very close values of K and ∂ 2 f/∂c 2 for these compositions.
Figure 12 shows the plot of modulation wavelength λ versus aging time t for B-50at.%Aalloy aged at 650 K for Ω A-B values of 12471 and 16628 J/mol.The modulation wavelength λ min is shorter as Ω A-B increases.Besides, the start time is shorter for the phase decomposition with the highest value of Ω A-B parameter, which can be attributed to a higher driving force.It is important to mention that the kinetic behavior of phase decomposition increases as the aging time increases in Figure 12.This fact can be related to the coarsening stage of the decomposed phases 6 .
The effect of the aging temperature on the modulation growth kinetics is shown in Figure 13 for B-50at.%Aalloy with a Ω A-B value of 16628 J/mol.The increase in the aging temperature caused a lower driving force, see Figure 1; nevertheless, the increase in the atomic mobility, due to the higher atomic diffusion, promoted that the phase decomposition started first at a higher aging temperature.The minimum modulation wavelength is very close for the two aging temperatures since they have the same interaction parameter.
Figure 14 shows the growth kinetics of spinodal decomposition for the A-50 at.%B alloy with Ω A-B values of 12471 J/mol (1500 R) and 16628 J/mol (2000 R) using both equations at 650 K.In both cases, the modulation wavelength remains almost constant with time, which is also a characteristic of the spinodal decomposition process 9 .It is evident that the modulation wavelength λ, is shorter for the nonlinear equation than for the linear one.In the case of the nonlinear equation, a higher value of Ω A-B caused that the start of the phase decomposition growth kinetics was faster.This is attributed to the higher values of the first derivative as Ω A-B increases.A higher first derivative value, or driving force, can be adopted as the reason for the shorter minimum wavelength λ min for spinodal decomposition using the nonlinear equation.The higher values for minimum wavelength λ min in the linear equation seem to be associated with the lower driving force, lower values of second derivative, as explained in section 3.1.Therefore, a longer aging time is needed to reach the longer wavelength λ min for the start of the spinodal decomposition.

Comparison with actual alloy systems
The use of low interaction Ω A-B parameters causes that the T c temperature of the miscibility gap to be located at low temperatures, see Equation 4. This promotes a low atomic diffusion which makes the phase decomposition process be very slow.For instance, Cu-Ni alloys 12 have a miscibility gap at low temperatures, lower than 573 K and its kinetic behavior of phase decomposition is very slow due to its low atomic diffusion at these temperatures.Furthermore, it has been reported that the phase decomposition takes longer aging times than several years.In contrast, the use of high Ω A-B values makes the T c temperature be higher which promotes the atomic diffusion and thus the appearance the phase decomposition.For instance, the phase decomposition in Fe-Cr alloys occurs at higher temperatures 6 and the phase decomposition process is much more rapid than that observed in Cu-Ni alloys.Table 2 shows a comparison of simulated results, the modulation wavelength and morphology of decomposed phases, for the present work using nonlinear equation with other simulated and experimental results for Cu-Ni and Fe-Cr alloys [14][15][16][17] .The simulation results corresponding to Ω A-B values of 8314 and 12471 J/mol have the T c temperature of the miscibility gap located at temperatures between 500 and 750 K.This temperature range is close to that of the miscibility gap in Cu-Ni alloys [14][15] .Furthermore, the minimum value of the modulation wavelength for the present work is in good agreement with the experimental and calculated values reported in the literature for Cu-Ni alloys.In contrast, the Ω A-B and T c values of 16628 J/mol are close to those in Fe-Cr alloys.The minimum value of modulation wavelength for the A-B alloy shows a good agreement with those reported in the literature for Fe-Cr alloys 16,17 .Finally, the morphology prediction of the decomposed phases to that reported in the literature.

Conclusions
Both linear and nonlinear Cahn-Hilliard equations reproduced the basic characteristics of the morphology and growth kinetics of the spinodal decomposition in the early    Cu-Ni [12,15] 8366 573 2 1-2 Irregular and interconected

Irregular and interconnected
Fe-Cr [16,17]  stages of aging for hypothetical A-B alloys.The increase in the thermodynamic interaction parameter Ω A-B caused a shorter start time for the phase decomposition.The higher aging temperature promoted the faster phase decomposition due to the higher atomic diffusion.The kinetic behavior for the phase decomposition calculated by the linear equation is slower than that of the nonlinear one which can be attributed to the difference in the driving force parameter of each equation.The microstructure evolution and kinetics of phase decomposition of Cu-Ni and Fe-Cr alloys showed a good agreement with the present work results.

1200 Figure 1 .
Figure 1.Plot of free energy f, vs. composition c.

Figure 3 .
Figure 3. Concentration profile of the A-50at.% B alloy aged at 650 K for different times for Ω A-B = 16628 J/mol.

Figure 4 .
Figure 4. Concentration profile of the A-30at.%Balloy aged at 650 K for different times for Ω A-B = 16628 J/mol.

Figures 9
Figures 9 and 10 show the microstructural evolution predicted for the B-50 at.%A and 30at.%Aalloys for the

Figure 5 .
Figure 5. Plot of second derivative of free energy f vs. composition c for Ω A-B = 16628 J/mol at 650 K.

Figure 6 .
Figure 6.Concentration profile of the A-50at.% B alloy aged at 650 K for different times for Ω A-B = 16628 J/mol with nonlinear equation.

Figure 7 .
Figure 7. Concentration profile of the A-50at.%Balloy aged at 650 K for different times for Ω A-B = 16628 J/mol with linear equation.

Figure 8 .
Figure 8. Plot of first derivative of free energy f vs. composition c

Figure 12 .
Figure 12.Plots of modulation wavelength λ vs. time t for B-50at.%Aalloy at 650 K for Ω A-B values of 12471 and 16628 J/mol.

Figure 14 .
Figure 14.Plots of modulation wavelength λ vs. time t for linear and nonlinear equations for Ω A-B values of 12471 and 16628 J/mol.

Table 1 .
Values of parameters for simulation.

Table 2 .
Comparison of present work simulated results with other works.