Unusual domain growth behavior in the compressible ising model

Unusual domain growth behavior in the compressible ising model

S. J. Mitchell^*; Luiz F. C. Pereira^†; D. P. Landau^‡

The Center for Simulational Physics, and the Department of Physics, University of Georgia, Athens, GA 30602-2451, USA

ABSTRACT

Large scale Monte Carlo simulations have been used to study long-time domain growth behavior in a compressible, two-dimensional Ising model undergoing phase separation. The system is quenched below the transition temperature from a random spin state, and we investigated the late-time domain size growth law, R(t) = A + Btⁿ. For "lattice mismatched" systems, we found n = 0.224 ± 0.004 which deviates significantly from the Lifshitz-Slyozov value of n = 1/3 for late-time growth . For a compressible model with no mismatch, we find only a slight deviation from n = 1/3. These results strongly suggest that we do not yet fully understand domain growth.

Keywords: Ising Model; Compressible Model; Phase Separation; Domain Growth

I. INTRODUCTION

The phenomenon of phase separation is extremely common in diverse condensed matter systems including magnets [1,2], alloys [3–14], fluids [15–20], and polymers [21]. Phase separation may be observed whenever a system is quenched from a homogeneous disordered phase into an ordered phase where multiple domains coexist. During phase separation, diffusion and surface tension cause the domain size to increase with time. For an infinite system, the domain size is predicted to grow as a power law [22],

where R is the domain size, t is the time after the quench occurs, A and B are constants which depend upon system specific details, and n is the domain growth exponent [1]. Stochastic domain growth with a conserved order parameter in the absence of hydrodynamic modes (Model B [23]) is believed to be in a class of domain growth with n = 1/3, independent of the dimensionality of the system. This value for the growth exponent, first predicted by Lifshitz and Slyozov [24], has since been observed in numerous simulation results of lattice models, most notably Ref. [1]. For reviews, see Refs. [25–27]. Thus, it would seem as though domain growth under these conditions would not show the rich diversity of behavior that is seen near critical points in different systems.

Monte Carlo simulations are well suited for the study of domain growth under these conditions, but because large system sizes are required to access the asymptotic growth limit region and multiple long runs are needed to measure n, few high precision computational studies exist that could verify the universality of the n = 1/3 growth law. Furthermore, until recently, high precision studies of compressible systems were outside the limits of available computational resources. Studies of Ising binary alloys [28–32] have now shown that compressibility can affect static critical phenomena and can have a non-trivial effect on the interfaces between domains [17]. There is reason to believe that domain growth may also be affected by compressibility [4, 6]; however, we are unaware of any high precision computational study that has confirmed any significant deviation from the theoretical prediction of n = 1/3.

II. MODEL AND METHOD

In order to investigate the universality of the n = 1/3 growth law, we have considered the simplest possible compressible model with coherent interfaces, a two-dimensional Ising model in which L² Ising spins have continuous positions within a periodic box of size L_x × L_y. Although the Ising model is most often used to describe magnetic systems, it is equivalent to a binary alloy in which the two spin values represent the two different species. The Ising model Hamiltonian is

where å_áijñ is a sum over nearest-neighbor pairs of particles, r_ij is the distance in dimensionless units between nearest-neighbor particles, s_i is the spin value of the ith spin with possible values of ±1, f(r,a,b) is the nearest-neighbor interaction potential (given below), J_q = 50 is the bond angle stiffness in dimensionless units, and å_áijkñ is the sum over bond angles q_ijk (four per particle), where i and k are nearest-neighbors of j. The second term stabilizes an ordered square lattice structure.

In this model the nearest-neighbor interaction is given by a Lennard-Jones-like potential:

where r is the displacement between two nearest-neighbors with spin values a and b, l_ab is the preferred bond length, J₊₊ = J_-- = 30, J_+- = J_-+ = J₊₊-2, and l_+- = l_-+ = 1. With these interactions, when l_ab = 1 and J₊₊ and J_q ® ¥, the Hamiltonian reduces to that of the common ferromagnetic, rigid (lattice) Ising model, allowing comparison with previous studies [1]. We will give all spatial units in terms of the l_+- bond length, and all energies and temperatures in units where J₊₊-J_+- = 2. The well-known Onsager critical temperature for the rigid (lattice) model is thus T_c = 2.269.

The model was investigated using Monte Carlo simulations with the standard Metropolis acceptance criterion within the canonical constant pressure ensemble [33]. Three types of Monte Carlo moves were used: a spin exchange (Kawasaki dynamic), in which a nearest neighbor pair of spins were randomly chosen and their spin values exchanged, lateral displacement, in which a particle was randomly selected and its position displaced by a small, random amount, and a global volume rescale (expansion or contraction), which was needed to maintain a constant pressure. The energy change for the global rescale had an additional effective term not shown in Eq. 2, which was needed to correctly reproduce a constant pressure ensemble [28]. Monte Carlo time was measured in units of Monte Carlo steps per site (MCS), where one MCS consisted of one attempted volume change followed by L² exchange or lateral displacement attempts, where the probability to choose exchange or displacement was 50%. For the rigid (lattice) model, only spin exchange moves could be implemented, and then one MCS consisted of L² exchange attempts.

Each simulation began with a state in which the value of each spin was randomly chosen to be +1 or -1 with 50% probability, and our Monte Carlo algorithm conserved the total magnetization at this initial value. The system was then equilibrated at a temperature (T = 7.0) which was well above the critical temperature (T_c). Following equilibration the temperature was quenched to 0.6T_c, and the configurations were recorded and analyzed every 10³ MCS.

III. RESULTS

We investigated three specific variations of the model, the rigid two-dimensional Ising model (i.e. the traditional Ising square lattice), a symmetric elastic net, where l_ab = 1, referred to as the 0% mismatch model, and a 4% lattice mismatch model, where l₊₊ = 1.02 and l_-- = 0.98. The 4% mismatch value was chosen to give insight into phase separation in SiGe alloys [28–30]. To decrease the computational effort of the simulations, the compressible systems were treated as distortable square nets, a reasonable approach for the low temperatures and relatively stiff models considered here.

All simulations were performed in the spinodal decomposition regime, as the configurations in Fig. 1 indicate, and all three variations of the model clearly show phase separation behavior [34]. No noticeable qualitative differences in the overall structure can be seen in the "snapshots" of the rigid and 0% mismatch models; however, the 4% mismatch model shows clear, diagonally oriented domains. Anisotropic configurations are often observed for systems containing a lattice mismatch [10]. A system with coherent interfaces, such as the one studied here, might be expected to stop phase separating at very late times (not reached in these simulations) [14], where the domain size would be finite due to interfacial energy.

Simulations were performed for various system sizes between L = 64 and L = 512. These showed that system sizes of L = 512 (512² particles)and time scales of 10⁶ MCS are needed to measure the asymptotic growth behavior. For the compressible models, each 10⁶ MCS run required one to two weeks of run time on a 2 GHz AMD Opteron processor, with a combined computational effort of ~ 5 cpu years, and over a terabyte was required to store and then analyze the resulting time series. The rigid model (lattice) simulations accounted for less than 5% of the total computational effort and are included for comparison to earlier work [1] that used a vectorized, multispin coding algorithm. Although we would not expect any difference in n, the vectorizing algorithm examined the spins in a different order and could yield different values of A and B. Note that larger system sizes and longer runs are not feasible for the compressible models with currently available resources, since L = 1024 only reaches 10⁵ MCS after two weeks of run time. On the other hand, simulations for L = 256 show finite size effects for t ~ 10⁷ MCS so longer simulations would not help for this size system.

Quantitative information about the phase separation behavior can be extracted by analyzing the spin-spin spatial correlation function; however, because the compressible models have continuous particle positions, methods developed for regularly spaced particles were not applicable. Therefore, the correlation function was calculated by direct summation over the "lattice" in only four "high symmetry" directions,

where (p,q) denotes one of the four displacement directions, i.e. the "lattice" directions, (1,0) and (0,1), and the "lattice" diagonal directions, (1,1) and (-1,1). Furthermore, s_i and s_j are the spin values of two particles with a displacement r in the (p,q) direction, and ás_is_jñ denotes an average over all such pairs of particles. For brevity, r is used here to denote an average over the actual interparticle displacements. The correlation length, x, is then defined as the first zero crossing of C(r) and is found by fitting a second order polynomial to the three points closest to the crossing.

The correlation functions were isotropic for both the rigid (lattice) model and the 0% mismatch model, and Fig. 2 shows C(r) averaged over all four directions for 0% mismatch. Both the rigid and 0% mismatch models show well-defined first zero crossings for all directions, and the zero crossing moves to greater distances as the time increases. When the correlation function is plotted against the "scaled" correlation length, see Fig. 3, the curves almost collapse onto each other. In contrast, with 4% mismatch we find quite anisotropic results, as one might expect from Fig. 1. For the (1,0) and (0,1) directions, see Fig. 4 the first zero crossings again move to greater distances with increasing time, moreover, additional zero crossings are visible on this scale. For t > 10⁴ MCS, no first zero crossing is observed in the diagonal directions, but we can perform alternative scaling analyses using the location of the first local minimum or the perimeter density (i.e. fraction of mixed bonds). These additional analyses yielded the same growth exponents as those obtained from the first-zero crossings presented below. When the correlation function is plotted against the "scaled" correlation length, see Fig. 5, the curves do not collapse onto each other, although the zero crossings occur quite close to each other.

The average correlation lengths are shown in Fig. 6. In all cases, power-law growth is observed at late times, but the form, x(t) = A + Btⁿ is not adequate for t < 10⁴. From non-linear fits to the data shown in Fig. 6 for t > 10⁴, we find that n = 0.332 ± 0.003 for the rigid model, n = 0.318 ± 0.005 for the 0% mismatch model, and n = 0.224 ± 0.004 for the 4% mismatch model. The error bars include small systematic uncertainties arising from the fitting procedure, such as the number of points, weighting of the data, and fitting range.

In addition to performing least squares fits of the long time behavior of the domain size, shown in Fig. 6, we also performed a ratio analysis [33], which is a more sensitive technique giving a "local" estimate of the growth exponent, i.e. over a relatively short time interval. The ratio estimate for n(t) is defined by

where x¢(t) = x(t)- A, A is obtained from fits shown in Fig. 6, and t and D are multiples of 10³ MCS. Because the ratio method is very sensitive to small statistical errors in x(t), bin averaging of the final sequence was necessary. Therefore, the calculation was repeated for different bin widths and different values of D and A, but such differences did not noticeably affect the estimate of n. The results from the ratio analysis are shown in Fig. 7.

The growth exponent estimates for the lattice (the rigid model) are in excellent agreement with the theoretical prediction of n = 1/3 [1, 23–25, 33]. The small deviation from n = 1/3 for 0% mismatch is within the fluctuation in the "local estimates" for n; hence we cannot firmly conclude that it is inconsistent with n = 1/3 at long times. Longer runs and larger systems would be needed to determine if this deviation is a real effect, but such simulations are beyond our current computational capabilities. However, for the 4% mismatch model the results clearly show a smaller growth exponent and anisotropic correlations. (Deviations from n = 1/3 have been seen before in systems with mismatch [4], but it was unclear if the asymptotic growth regime had been reached, and more recent studies [9] indicated n = 1/3.)

We also used large scale Monte Carlo simulations to study the static critical behavior of these models. Data for the magnetization distribution P(m) were analyzed using histogram reweighting in H - T space to adjust the estimate for the location of the critical point until the finite size scaled distributions sP(m) vs m/s, where s is the 2nd moment of the distribution, were rather well superimposed. Data for the 4% mismatch model are shown in Fig. 8. The scaled data for L = 16 appear to show small corrections to finite size scaling, and if these are omitted the quality of the scaling is improved. For larger sizes any finite size corrections are so small that it seems highly unlikely that data for still larger lattices would bring scaling curves into agreement with that for the Ising square lattice which is included in the figure for comparison. Note that field mixing was not employed, and it is possible that an even better scaling could be obtained if it were implemented; but even without field scaling the data scale extremely well.

The "universal" curve expected for systems in the same universality class as the two dimensional Ising model is dramatically different and strongly suggests that the compressible, Ising net is not in the same static universality class as the Ising lattice model.

IV. CONCLUSION

The inclusion of compressibility and lattice mismatch in the Ising model Hamiltonian can alter the domain growth exponent, and our results can be readily generalized to systems with differing sizes or differing bond lengths. We do not know if the deviations from n = 1/3 indicate a breakdown of universality or if they indicate new classes of domain growth. Obviously our current understanding of domain growth is incomplete, and further theoretical and computational consideration is needed.

Acknowledgments

The authors thank S. H. Tsai and K. Binder for helpful comments and discussions. This work was funded by NSF grants #DMR-0341874 and #DMR-0307082.

Received on 28 August, 2007

• [1] J. G. Amar, F. E. Sullivan, and R. D. Mountain, Phys. Rev. B 37, 196 (1988).
• [2] P. Fratzl, J. L. Lebowitz, O. Penrose, and J. Amar, Phys. Rev. B 44, 4794 (1991).
• [3] A. Sadiq and K. Binder, J. Stat. Phys. 35, 517 (1984).
• [4] H. Nishimori and A. Onuki, Phys. Rev. B 42, 980 (1990).
• [5] A. Maheshwari and A. J. Ardell, Phys. Rev. Lett. 70, 2305 (1993).
• [6] C. Sagui and R. C. Desai, Phys. Rev. Lett. 74, 1119 (1995).
• [7] C. A. Laberge, P. Fratzl, and J. L. Lebowitz, Phys. Rev. Lett. 75, 4448 (1995).
• [8] V. I. Gorentsveig, P. Fratzl, and J. L. Lebowitz, Phys. Rev. B 55, 2912 (1997).
• [9] P. Nielaba, P. Fratzl, and J. L. Lebowitz, J. Stat. Phys. 95, 23 (1999).
• [10] D. Orlikowski, C. Sagui, A. M. Somoza, and C. Roland, Phys. Rev. B 62, 3160 (2000).
• [11] R. Weinkamer, P. Fratzl, H. S. Gupta, O. Penrose, and J. L. Lebowitz, 4894 (1995)., Phase Transitions 77, 433 (2004).
• [12] K. Binder and D. Stauffer, Phys. Rev. Lett. 33, 1006 (1974).
• [13] P. Fraztl, O. Penrose, and J. L. Lebowitz, J. Stat. Phys. 95, 1429 (1999).
• [14] B. J. Schulz, B. Dünweg, K. Binder, and M. M uller, Phys. Rev. Lett. 95, 096101 (2005).
• [15] T. Koga and K. Kawasaki, Physica A 196, 389 (1993).
• [16] A. T. Bernardes, T. B. Liverpool, and D. Stauffer, Phys. Rev. E 54, R2220 (1996).
• [17] S. K. Das, S. Puri, J. Horbach, and K. Binder, Phys. Rev. Lett. 96, 016107 (2006).
• [18]M. Grant and K. R. Elder, Phys. Rev. Lett. 82, 14 (1998).
• [19] V. M. Kendon, J. C. Desplat, P. Bladon, and M. E. Cates, Phys. Rev. Lett. 83, 576 (1999).
• [20] P. B. Warren, Phys. Rev. Lett. 87, 225702 (2001).
• [21] T. Hashimoto, M. Itakura, and H. Hasegawa, J. Comp. Phys. 85, 6118 (1986).
• [22] D. A. Huse, Phys. Rev. B 34, 7845 (1986).
• [23] P. C. Hohenberg and B. Halperin, Rev. Mod. Phys. 49, 435 (1977).
• [24] I. M. Lifshitz and V. V. Slyozov, J. Phys. Chem. Solids 19, 35 (1961).
• [25] A. J. Bray, Adv. Phys. 43, 357 (1994).
• [26] A. Onuki, Phase Trasition Dynamics (Cambridge University Press, Cambridge, 2002).
• [27] J. D. Gunton, M. S. Miguel, and P. S. Sahni, in Phase Transitions and Critical Phenomena, Vol. 8, edited by C. Domb and J.L. Lebowitz (Academic Press, New York, 1983).
• [28] B. Dünweg and D. P. Landau, Phys. Rev. B 48, 14182 (1993).
• [29] M. Laradji, D. P. Landau, and B. D unweg, Phys. Rev. B 51, 4894 (1995).
• [30] F. Tavazza, D. P. Landau, and J. Adler, Phys. Rev. B 70, 184103 (2004).
• [31] D. P. Landau, F. Tavazza, and J. Adler, Comp. Phys. Commun. 169, 149 (2005).
• [32] E. M. Vandeworp and K. E. Newman, Phys. Rev. B 55, 14222 (1997).
• [33] D. P. Landau and K. Binder, A Guide to Monte Carlo Simulations in Statistical Physics, 2nd ed. (Cambridge University Press, New York, 2005).
• [34] S. M. Mitchell and D. P. Landau, Phys. Rev. Lett. 97, 025701 (2006).

Publication Dates

• Publication in this collection
  27 Mar 2008
• Date of issue
  Mar 2008

History

• Received
  28 Aug 2007