Comment on ‘ ‘ Kinetics of spinodal decomposition in the Ising model with vacancy diffusion ’ ’

In a recent paper @Phys. Rev. B50, 3477~1994!#, P. Fratzl and O. Penrose present the results of the Monte Carlo simulation of the spinodal decomposition problem ~phase separation ! using the vacancy dynamics mechanism. They observe that the t growth regime is reached faster than when using the standard Kawasaki dynamics. In this Comment we provide a simple explanation for the phenomenon based on the role of interface diffusion, which they claim is irrelevant for the observed behavior.

In the study of nonequilibrium growth phenomena two processes are considered as prototype. 1Both correspond to the evolution of systems quenched from a high-temperature disordered phase to a temperature below an equilibrium phase transition.The first case is the growth of ordered domains ͓nonconserved order parameter ͑NCOP͔͒ and the second the spinodal decomposition problem ͓conserved order parameter ͑COP͔͒ also called phase separation.There is general agreement 2 in the fact that in both cases the mean size of the domains grows, during the late stages of the evolution, following a power law.Such laws ͓R(t)ϳt x ͔ are characterized by the growth exponent x that takes the value x NCOP ϭ1/2 ͑Allen-Cahn law 3 ͒ for the nonconserved case and x COP ϭ1/3 ͑Lifshitz-Slyozov law 4 ͒ for the conserved case.These values are considered to be universal, independent of details of the model and dimensionality.They have been corroborated by theoretical arguments, numerical analysis of Langevin equations, and Monte Carlo simulations of lattice systems. 2 For the case of COP and concentration cϭ0.5 the first Monte Carlo results suggested that the growth exponent was x T Ӎ0.25. 5 After a wide discussion it was shown that this behavior was just a transient regime with enhanced diffusion along the interfaces 6,7 corresponding to a growth law R(t)ϭR 0 ϩt 1/3 .For very large systems a correct extrapolation revealed the existence of the xϭ1/3 law at late times, 8 dominated by the bulk diffusion.
Most Monte Carlo simulations in the literature have been carried out using the standard Kawasaki dynamics, proposing neighboring atom-atom exchanges homogeneously on the lattice, and the Metropolis acceptation probability, which ensures the evolution towards equilibrium.For a lattice with N sites, the time unit, called a Monte Carlo step ͑MCS͒, is defined as N proposals of atom-atom exchanges, independently of its acceptance ratio.Such a definition is based on the assumption that a system out of equilibrium, but in con-tact with a heat bath, is excited ͑jump proposals͒ with a rate constant in time and homogeneous in space.
Recently a more realistic vacancy dynamics has been proposed for the study of such growth processes: 9-14 a very small amount of vacancies is introduced in the system, almost not affecting its equilibrium properties.Only neighboring atom-vacancy exchanges are proposed and accepted or not according to the standard Metropolis rules.The time unit is also defined as N proposals, which are inhomogeneously distributed since the vacancy walk through the lattice is very correlated with the current order existing in the system. 11evertheless, to our knowledge, such an inhomogeneous time scale has not been analyzed by general theoretical arguments.
1][12] The physical reason for this acceleration is that the vacancy stays most of the time on the interfaces, avoiding a waste of time in the already ordered bulk.
In a recent paper 13 P. Fratzl and O. Penrose point out the usefulness of the vacancy dynamics also for the COP case in order to reach the asymptotic xϭ1/3 regime much faster than when using the standard Kawasaki dynamics, suppressing the transient xӍ0.25 regime.
We have performed Monte Carlo simulations with the vacancy dynamics of a system with the same conditions as in Ref. 13 ͑square lattice with linear size Lϭ128, cϭ0.5, and cϭ0.1, 1 vacancy and reduced temperature Tϭ0.5). 15We have measured the structure factor S(k,t) at different time steps, averaging over the (1,0) and (0,1) directions in reciprocal space.Also, averages over 43 and 20 independent runs for cϭ0.5 and cϭ0.1, respectively, have been taken.We have estimated the mean domain size by measuring the first moment of the structure factor: In Fig. 1͑a͒ we show the evolution of k 0 (t) and compare with the xϭ1/3 behavior for cϭ0.5 and cϭ0.1.Fig. 2 shows the scaled structure factors k 0 2 S(k/k 0 ) for different time steps, revealing the existence of dynamical scaling, at least, from 10 3 to 10 5 MCS for the two studied concentrations, in agreement with the results in Ref. 13.The straight line corresponds to the behavior k Ϫ3 predicted by Porod's law. 2 In order to measure the growth exponent, we have computed the logarithmic derivative of k 0 (t) defined as Results for ␣ϭ10 are shown in Fig. 1͑b͒ and compared with the value xϭ1/3 ͑dashed line͒.The error bars have been estimated from the difference between the values of k 0 obtained from the structure factor in the (1,0) and (0,1) directions.We have tested that such error bars include the deviations arising from the change in the derivation step in Eq. ͑2͒ between ␣ϭ1.11 and ␣ϭ100.
For cϭ0.5, a value of the dynamical exponent xӍ1/3 is reached after the first 10 3 MCS, in agreement with the results of Fratzl and Penrose. 13Furthermore our analysis also reveals, first, that the exponent continues increasing towards higher values (xӍ0.4)and, second, that after a maximum it decreases again towards values around xӍ1/3.
For cϭ0.1 we find that the dynamical exponent monotonously increases from very low values.After 3ϫ10 5 MCS it still remains below xϭ1/3.In our opinion this result is in agreement with those obtained by Fratzl and Penrose ͑according to their Fig.3͒ using both the atom-atom exchange and the vacancy mechanisms.On the basis of Fig. 2 of their paper, the authors claim that the vacancy mechanism accelerates the process.We think that such analysis of R vs t 1/3 only demonstrates the change of the prefactor of the growth law but not the change of the dynamical exponent.The study of the dynamical exponent is better performed on a log R vs log t as in their Fig. 3.By comparing their Fig.3a and 3b we cannot see any significant change in the slope enough to justify the change in the dynamical exponent.We believe that a change in the prefactor may well exist ͑as suggested by their Fig.2͒, but this could depend on fine details of the simulations, such as time correlation between configurations, diffusion constants, energy barriers, etc.Therefore, contrarily to their conclusions, we think that for cϭ0.1 the vacancy mechanism does not change the dynamical exponent.
The authors of Ref. 13 claim that the transient regime xӍ0.25 is avoided since the vacancy mechanism introduces a curvature dependence in the surface tension and that interface diffusion does not play any role.Such curvature dependence of the surface tension is introduced from a phenomenological point of view, and the irrelevance of the interface diffusion is justified only from the paradox that, according to Huse 6 and Vengrenovitch, 16 such diffusion enhances the xӍ0.25 regime.Here we shall provide a simple explanation, based on the existence of interface diffusion and heterogeneous excitations, for the change in the growth exponent from a value xӍ0.25 towards xӍ1/3 in the initial stages of the COP case with cϭ0.5.The argument also applies to the change from xӍ1/2 towards xӍ1 in the NCOP case. 11et us consider a quenched system for which, under homogeneous excitations, the mean domain size evolves following a power law R(t)ϳt z , and the amount of interface decreases like A(t)ϳt Ϫz .If a small amount of vacancies are introduced in the system, for energetic considerations they will concentrate on the domain interfaces destroying the most unfavorable bonds ͓as shown for both NCOP ͑Ref.11͒ and COP systems ͑Ref.13͔͒.This fact immediately implies an accumulation of the excitations on the domain boundaries compared to the excitations in the bulk.If the physical process underlying the growth takes place on the interfaces ͑as happens in the NCOP case and in the intermediate stages of the COP case for cϭ0.5) the growth will obviously be accelerated by the vacancy mechanism.We can estimate the change in the growth exponent in the limiting case in which the following two hypotheses hold: only the excitations on the interface are responsible for the evolution of the system, and the vacancies always stay on the interface.
Consider the usual Kawasaki exchange dynamics.During a Monte Carlo step (N trials͒ only A(t) are useful while NϪA(t) represent excitations in the bulk that are useless for the growth.Within the limiting case of the vacancy dynamics all the N trials would be useful ͑concentrated on the inter-face͒.This means that a configuration that has been reached after t MCS using standard Kawasaki dynamics, is reached at time tЈϽt with vacancy dynamics according to where t 0 is an initial transient time and the last aproximation holds for large enough t.With this new time scale the growth law will be RЈ͑tЈ͒ϵR͑t ͒ϭR͑ tЈ 1/͑1Ϫz ͒ ͒ϳtЈ z/͑1Ϫz ͒ .͑4͒ Thus, the exponent for the vacancy mechanism (zЈ) is related to the homogeneous one (z) according to zЈϭz/(1Ϫz).Such enhancement only applies when the dominant process for the growth is of interfacial nature, like interface diffusion in the intermediate stages of the COP case for cϭ0.5 (zӍ0.25⇒zЈӍ1/3) and interfacial ordering in the NCOP case (zӍ1/2⇒zЈӍ1).The above scenario is also in agreement with the following three facts.
͑1͒ The vacancy mechanism does not increase the exponent during the late stages of the COP case for cϭ0.5 since then the bulk diffusion becomes responsible for the growth.This might be an explanation for the observed decrease of the exponent towards xӍ1/3 in the late stages of the evolution in Fig. 1 ͑filled circles͒.Nevertheless, we cannot completely discard that such a decrease is a finite size effect.Only simulations to larger systems will clarify this point.
͑2͒ For cϭ0.1, the growth enhancement due to the vacancy mechanism almost disappears.In this case the pattern consists of isolated domains of the minority phase enbedded in a sea of the majority phase.One then expects that the interfacial diffusion is not relevant for the growth, and consequently the vacancy mechanism cannot accelerate the evolution.
͑3͒ We expect that the growth enhancement due to the vacancy mechanism is reduced when the temperature is increased.This is because thermal fluctuations compel the vacancy to scan the lattice more homogeneously. 11In our opinion results shown in Fig. 2 of the Rapid Communication by Fratzl and Penrose, are consistent with this explanation.
In conclusion, we have revised and completed the study by Fratzl and Penrose of the domain growth in the spinodal decomposition problem using the vacancy mechanism.In agreement with them we have found ͑i͒ that there is dynamical scaling and ͑ii͒ that for cϭ0.5 the system evolves with a growth exponent which reaches a value xӍ1/3 very fast.Contrarily with their results, extending the simulations up to 10 6 MCS, we obtain ͑iii͒ that for cϭ0.5 the exponent continues growing to higher values and ͑iv͒ that for cϭ0.1 the exponent takes monotonously increasing values below xϭ1/3.We have provided a very simple explanation for the observed behaviors based on the existence or absence of diffusion along the interface, which Fratzl and Penrose claim is irrelevant for the process.We acknowledge A. M. Lacasta and J. M. Sancho for fruitful discussions, the Comisio ´n Interministerial de Ciencia y Tecnologı ´a ͑CICyT͒ for financial support ͑Project number MAT92-884͒, and the Centre de Supercomputacio ´de Catalunya ͑CESCA͒ and Fundacio ´Catalana per a la Recerca ͑FCR͒ for computing facilities.C.F. also acknowledges financial support from the Comissionat per a Universitats i Recerca ͑Generalitat de Catalunya͒.

FIG. 1 .
FIG. 1. ͑a͒ Evolution of the reciprocal mean domain size as a function of time in log-log scale for cϭ0.5 (᭹) and cϭ0.1 (᭺).Data for cϭ0.5 correspond to averages over 43 different runs up to 10 5 MCS and 12 different evolutions from 10 5 to 10 6 MCS; for cϭ0.1 averages are taken over 20 different evolutions.Error bars are smaller than the symbols.The dashed line shows the xϭ1/3 behavior.͑b͒ Logarithmic derivative of the curves shown in ͑a͒ in log-linear scale displaying the effective exponent as a function of time.The dashed line indicates the value xϭ1/3.

FIG. 2 .
FIG. 2. Scaled structure factor in log-log and linear ͑insets͒ scales for the two studied concentrations.The different symbols correspond to different times as indicated in the legend.Data correspond to averages over 43 independent runs for cϭ0.5 and 20 independent runs for cϭ0.1.The line corresponds to the behavior k Ϫ3 predicted by Porod's law.