- Email: [email protected]
A binary phase field crystal study for phase segregation of liquid phase heteroepitaxial growth
Computational Materials Science 123 (2016) 65–69
Contents lists available at ScienceDirect
Computational Materials Science journal homepage: www.elsevier.com/locate/commatsci
A binary phase field crystal study for phase segregation of liquid phase heteroepitaxial growth Yingying Peng ⇑, Yanli Lu ⇑, Zheng Chen, Genggeng Yu State Key Lab of Solidification Processing, Northwestern Polytechnical University, 710072 Xi’an, PR China
a r t i c l e
i n f o
Article history: Received 26 March 2016 Received in revised form 16 June 2016 Accepted 18 June 2016
Keywords: Binary PFC Atomic sizes Mobility differences Liquid phase heteroepitaxial
a b s t r a c t The binary phase field crystal (PFC) model is employed to investigate influences of atomic sizes and mobility differences on the liquid phase heteroepitaxial. It was found that large size atoms are driven toward regions of tensile stress which correspond to peaks in a compressively strained film but to valleys in a film with tensile strain. Small size atoms are on the contrary in contrast for large size atoms. Due to the existence of vertical separation and lateral separation resulting from atomic size differences, the epitaxial layer exhibits double-island phenomenon: light and dark islands. In the presence of different mobilities, atoms with greater mobility accumulate at the film surface. And in the process of epitaxial growth, there are misfit dislocations nucleating in valleys where the strain is highest. Ó 2016 Elsevier B.V. All rights reserved.
1. Introduction Due to unique optical and electronic properties, with potential applications in optoelectronics and semiconductor devices, self-assembled of nanostructure in semiconductor heteroepitaxy has been intensively studied [1–4]. These structures include junctions, quantum dot islands and superlattices in which planar interfaces are highly desire [5–7]. However, the difference of lattice constant between the deposited material and substrate results in mismatch strain. Morphological instability induced by the strain leads to nonplanar surfaces and defects in the film. This stressdriven morphological instability [8,9] has been studied extensively, and it is well established that the morphological instability can change the growth mode of a thin film from a layer-by-layer mechanism to a three-dimensional island mechanism. However, previous works mainly focus on the study of pure materials [10–12], with the development of alloy thin films, more and more investigations about alloy heteroepitaxy have been reported. A latest report on alloy heteroepitaxy from Elder [13] et al has shown influences of substrate thickness on the growth of the single island and multiple islands. Compared with pure materials, alloy thin films have more complicated compositions. As for epitaxial growth of alloy thin films, alloy components may be prone to phase separation and according to Vegard’s law, the lattice constant of the film is generally a function of alloy compositions. Therefore, the growth of alloy thin films needs further study. Recently, alloy phase ⇑ Corresponding authors. E-mail addresses: [email protected] (Y. Peng), [email protected] (Y. Lu). http://dx.doi.org/10.1016/j.commatsci.2016.06.017 0927-0256/Ó 2016 Elsevier B.V. All rights reserved.
separation [14–16] has captured many scholars’ attentions. The main issues motivated us are that how atomic sizes and mobility differences affect processes of phase separation and further influence the stability of epitaxial thin films. During the actual production process of thin films, it is very difficult to control the thickness of thin films because of the fast growth rate. In order to get nanoscale films, we find another way known as phase field crystal model to predict the growth of thin films. Liquid phase epitaxy as a mature technology has been used in the production of semiconductor optoelectronic devices, as well as magnetic garnets, superconductors, ferroelectrics and other optical materials. Therefore, the purpose of this paper is to illustrate how the two-dimensional binary PFC model [17,18] addresses compositional effects in alloy liquid phase heteroepitaxy, focusing on influences of atom size differences and mobility differences on phase separation. 2. Simulation method 2.1. Model The binary phase field crystal model proposed by Elder [19] recently has not only incorporated elasticity and plasticity on atomic length and diffusive timescales and multiple crystal orientations which are caused by the periodic oscillations of the crystal lattice, but also involves the solidification, phase segregation, solute diffusion and other important features. For a binary alloy made up of A and B atoms, to define the total number density q ¼ qA þ qB . To simplify calculations it is convenient to fist
66
Y. Peng et al. / Computational Materials Science 123 (2016) 65–69
, introduce the following dimensionless fields: nA ¼ ðqA qA Þ=q . Also, it indicates that the atomic number density nB ¼ ðqB qB Þ=q field n and the concentration field u (To show the atomic density difference, it can be approximated as the concentration field), which gives
dynamics of the fields n and u decouple, the respective equations of motions have the form:
@n dF ¼ M r2 @t dn
ð5Þ
Þ=q n ¼ nA þ nB ¼ ðq q u ¼ ðnB nA Þ þ ðq B q A Þ ¼ ðqB qA Þ=q
@u dF ¼ M r2 du @t
ð6Þ
ð1Þ
The numerical discretization of Eqs. (5) and (6) is implemented using the mixed explicit-implicit Fourier transformations and operator splitting by Tegze et al. [20].
where F is the free energy functional given by
Z F¼
n t v w u K ! dx ^ n n3 þ n4 þ cu þ u2 þ u4 þ j r u2 j 4 2 3 2 4 2 ð2Þ
where ^ Bl0 þ Bl2 u2 þ Bs0 ð2r2 þ r4 Þ þ 4auBs0 ðr2 þ r4 Þ. Briefly outlining the physical meaning of the model parameters is as follows: the Bl0 , Bl2 , Bs0 , t, v , w, K, c, a and u are constants, whose values depend on the material parameters. t, c and v are model parameters associated with the Taylor coefficient. The length scale of phase separation is decided by the interplay of u, K and w. Bl2 determines liquid volume modulus with the changing of atomic concentration field. Bl0 is the liquid temperature. Bs0 is associated with the solid phase of the elastic constants. DB ¼ Bl0 Bs0 corresponds to temperature expressing the initial driving force of crystallization, which can be enhanced by lowering DB. a represents the solute diffusion factor measuring the size difference between A and B atoms, and is defined as a ¼ ðaB aA Þ=a0 where aA ; aB ; and a0 are the relaxed-state lattice parameters for pure A, pure B and the alloy film. Also, in the spirit of keeping calculations as simple as possible without losing the basic physics contained in the model, c ¼ 0 in the free energy. Minimizing the free energy, h pffiffiffi pffiffiffi i n ¼ A 1=2 cosð2qy= 3Þ cosðqxÞ cosðqy= 3Þ will be used. Substituting this expression into (2) and minimizing with respect to q and A, recalling that u is assumed constant over the scale that n pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffi varies, getting qeq ¼ 3=ð2RÞ and Amin ¼ 4ðt þ t 2 15v DBÞ=ð15v Þ. To simulate microstructure formation in binary alloys, dynamical equations of motion for the fields n and u need to be developed. @n dF dF @nA @nB dF dF þ M B r2 ¼ þ ¼ M A r2 ¼ r M1 r þ r M2 r @t @t @t dn du dnA dnB
2.2. The chosen of parameters The equilibrium diagram is shown in Fig. 1 for the appropriate parameters set as follows: Bs0 ¼ 1:0, Bl2 ¼ 1:8, t ¼ 0:6, L ¼ 2:65, l ¼ 4:0, a ¼ 0:26 and w ¼ 0:088. Our simulation is conducted on a periodic system of size 512Dx 512Dx, the binary alloy film with average density difference u ¼ 0 grows from a liquid phase above the bulk coherent spinodal temperature, the time step Dt ¼ 0:5. According to the equilibrium diagram, we choose DB0 ¼ 0:00886 to implement all of the following simulations, growth at this DB0 above the miscibility gap is typical of experimental conditions and should ensure that phase separation is driven by local stresses and is not due to spinodal decomposition. Initial conditions consisted of a binary unstrained crystalline planar substrate of eight atomic layer thickness in the bottom of the simulation area which is below a symmetric supercooled liquid of components A and B. In what follows the misfit strain is defined as ðafilm asub Þ=asub , for a symmetric mixture of A and B components afilm ¼ ðaA þ aB Þ=2. 3. Results and discussion For the reason that the evident contrast as bright and dark regions can highlight the phase separation in the alloy, grey-scale maps are adopted. The results are presented in Fig. 2 for two components with different atomic sizes and same mobility. Fig. 2 (a) and (c) show the atomic number density difference while Fig. 2(b) and (d) demonstrate the local concentration difference. As shown in Fig. 2(a) and (b), at the initial stage of growth, epitaxial
0.10
¼ ðM A þ M B Þr2 ð^n tn þ mn3 þ cÞ þ ðM B M A Þr2 K r2 u
0.08
þ ðM A M B Þr2 ðn2 Bl2 u þ 2an2 BS0 ðr2 þ r4 Þ þ wu þ lu3 Þ ð3Þ
0.06
@u dF dF @nA @nB dF dF M B r2 ¼ ¼ r M2 r ¼ M A r2 þ r M1 r @t @t @t dn du dnA dnB
0.04
þ ðM A þ M B Þr
2
ðn2 Bl2
2
3
n2 BS0 ð
u þ 2a
ΔΒ0
¼ ðM A M B Þr ð^n tn þ mn þ cÞ ðM B þ M A Þr K r u 2
2
r þ r Þ þ w u þ lu Þ 2
4
0.02 0.00
3
ð4Þ
Liquid
ls ΔΒ 0 =0.0213 c ΔB0 =-0.0139
Solid (0,0.00886)
-0.02 -0.04
where M1 ¼ ðMA þ MB Þ=q2 and M2 ¼ ðM B MA Þ=q2 . To get Eqs. (3) and (4), we used the approximate value of Eq. (1): n ¼ nA þ nB , u ¼ nA nB . Eqs. (3) and (4) couple the dynamics of the fields n and u through a symmetric mobility tensor. The dependence of mobility M A and MB in general depends on local crystal density and the local relative concentration of species A and B. While M A and MB are different, a simple Euler algorithm was used for the time derivative and the spherical Laplacian approximation was used. Assuming a substitutional diffusion between species A and B, that the same M mobility applies for the two species, furthermore, considering that the mobility coefficient is a constant, the
-0.06 -0.08 -0.10 -0.25 -0.20 -0.15 -0.10 -0.05 0.00 ϕ
0.05
0.10
0.15
0.20
0.25
Fig. 1. Phase diagram as a function of DB0 u for two-dimensional triangular system. The blue solid line corresponds to liquid phase, the red solid line corresponds to solid phase, the black square dot corresponds to selected initial condition. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
Y. Peng et al. / Computational Materials Science 123 (2016) 65–69
(a)
(b)
(c)
(d)
67
Fig. 2. Plot of epitaxial growth on a substrate of eight atomic layers in thickness. (a) and (c) show the atomic number density difference. (b) and (d) explain the local concentration difference, in which light regions represent component A which consists of large size atoms and dark regions represent component B which consists of small size atoms. In this figure, M A ¼ M B , misfit strain is 0.088. Times from top to bottom are 80Dt and 1040Dt. Inset: the local concentration difference with 1800-time step. The start and the end of arrows correspond to positions of dislocations nucleation and nucleated dislocations, respectively.
films become unstable and nonplanar. Even without compositional effects, the stress nonuniformity from misfit between substrate and films is sufficient to drive the surface instability as in a single-component film. The addition of the alloy components increases the instability of thin films. Composition variations cause an enhancement of the instability in two ways. First, composition variations increase the stress nonuniformity, causing a larger driving force for the surface instability. Second, the stress nonuniformity results in the enhanced diffusion of a kind of atomic species to surface peaks, causing the instability to grow faster. As presented in Fig. 2(d), it is clear to see that there exists an obvious difference compared to pure epitaxial film. The epitaxial layer exhibits double islands phenomenon: light and dark islands. It fully shows alloy segregation instability couples with the morphological instability, resulting in the simultaneous modulations of the surface profile and alloy composition during epitaxial growth. The strain modulations associated with morphological perturbations may enhance the segregation of depositing alloy species. Reasons for double-island phenomenon are the presence of lateral and vertical separation due to different size components. Although the two components have equal mobility, they own different atomic sizes which result in solute strain. Also, the film is in compressive stress because of the positive misfit strain, and peaks correspond to tensile stress regions, valleys correspond to regions of compressive stress. Therefore, in the presence of the solute strain and the misfit strain, large size atoms are driven toward regions of tensile stress which correspond to peaks and small size atoms are driven toward regions of compressive stress which correspond to valleys in a compressively strained film. With the time going by, large size atoms gathered in the regions of peaks while small atoms gathered in the regions of valleys, that are double islands. What is more, in addition to the vertical separation, there exists lateral phase separation. This coupling creates a lateral phase separation on the length scale of the surface instability and has been predicted as well as verified for binary films [21,22]. The epitaxial system can maintain a lower free energy by transferring
atoms from the valleys to peaks, because the transition leads to a decrease in the contact area between the substrate and a new layer. Thus, atomic transition to upper layers leads to a relaxation in the local strain field. Lateral segregation has two important implications: on the one hand, if epitaxial films thickness had been increased during growth, misfit dislocations would have nucleated at the interface, most probably in valleys where the strain is highest, as shown in the inset of Fig. 2(d); on the other hand, the lateral segregation reduces large size atoms in the regions of the highest strain, thereby reducing the strain compared to a similar morphology without segregation, and hence, retarding dislocation nucleation. Thus, compared with single component, not all islands coalesce with each other and only a few dislocations nucleation in the structure of double islands. What’s more, due to the lateral segregation, light islands grow much quicker than dark islands. By comparing Fig. 2 with Fig. 3, it is clear to find that the sign of misfit strain is the only difference between them. Similarly, the negative misfit strain has some influences on the epitaxial growth. Due to the negative misfit strain, the film is under tensile stress, and peaks correspond to compressive stress regions, valleys correspond to regions of tensile stress. Therefore, in the presence of the solute strain and the misfit strain, large size atoms are driven toward regions of tensile stress which correspond to valleys and small size atoms are driven toward regions of compressive stress which correspond to peaks in a tensile strained film. Accompanied with the epitaxial growth, large size atoms gathered in regions of valleys and small size atoms gathered in regions of peaks, and then double islands are formed: dark and light islands. Surely, in addition to the vertical separation, there exists lateral phase separation. For the same purpose, in order to lower free energy, atoms move from valleys to peaks, the number of small size atoms in valley regions becomes much less. Meanwhile, there are misfit dislocations nucleating in valleys where the strain is highest. Similarly, not all islands coalesce with each other and only a few dislocations nucleation in the structure of double islands as shown in the inset of Fig. 3. What’s more, by comparing the inset of Fig. 2 with that of
68
Y. Peng et al. / Computational Materials Science 123 (2016) 65–69
(a)
(b)
(c)
(d)
Fig. 3. Plot of epitaxial growth on a substrate of eight atomic layers in thickness. (a) and (c) show the atomic number density difference. (b) and (d) explain the local concentration difference, in which light regions represent component A which consists of large size atoms and dark regions represent component B which consists of small size atoms. In this figure, M A ¼ M B , misfit strain is 0.088. Times from top to bottom are 80Dt and 600Dt. Inset: the local concentration difference with 840Dt, dislocations are enclosed by black circles.
Fig. 3, under negative misfit strain, misfit dislocations are earlier to nucleate, therefore, under the condition of positive misfit strain, double islands are much more stable and need much more time to coalesce. Above all, the mechanism that large atoms preferentially segregate to peaks if the misfit is compressive but to valleys if the misfit is tensile is consistent with experimental evidence in the AlAs/InAs system [23]. In order to describe influences of atomic sizes on the phase separation in detail, maximum thickness of atomic layer is shown as a function of time in Fig. 4 by combining with Figs. 2 and 3. In Fig. 4, curve1 and curve 2 exhibit the atomic layer thickness of small size atoms and large size atoms under the positive misfit strain, respectively. Curve 3 and curve 4 show the atomic layer
1
atomic layer thickness
25
2 3 4
20
15
10
5
0
400
800
1200
1600
2000
time/Δt Fig. 4. Maximum thickness of atomic layer is plotting as a function of time. Curve 1 and curve 2 exhibit the atomic layer thickness of small size atoms and large size atoms under the positive misfit strain, respectively. Curve 3 and curve 4 show the atomic layer thickness of large size atoms and small size atoms under the negative misfit strain, respectively.
thickness of large size atoms and small size atoms under the negative misfit strain, respectively. While the misfit strain is positive, large size atoms are driven toward positions of peaks and small size atoms are driven to valleys. Finally, large size atoms gathered in regions of peaks and smaller size atoms gathered in regions of valleys. In this way, double islands are formed. Also, it is shown that the atomic layer thickness of large size atoms is much larger than that of small size atoms. Likely, while the misfit strain is negative, double islands are formed in the way that large size atoms gathered in regions of valleys and small size atoms gathered in regions of peaks, which is on the contrary compared with double islands under positive misfit strain. In addition, it is clear to see that the difference of the atomic layer thickness of double islands in positive misfit strain is larger than that under negative misfit strain. This phenomenon in another way proves that under the condition of positive misfit strain, the double islands are much more stable and need much more time to coalesce. Fig. 5 shows influences of unequal mobilities on the epitaxial growth, the left column shows the atomic number density difference and the right column exhibits the local concentration difference. Dislocations are enclosed by black squares, circles and arrows. Light numbers show positions of phase separation. The misfit strain can induce this morphological instability because the stressed material with perturbed surface has the nonuniform strain energy. The stress in the film with sinusoidally perturbed surface is relaxed at the wave peaks and increased at the troughs by equal magnitudes. Also, the misfit strain is positive, the film is in compressive stress. Therefore, at the initial stage of growth, misfit dislocations nucleated at the interface in troughs where the compressive strain is highest (see regions marked circles and arrows in Fig. 5). And the atoms with greater mobility accumulate at the film surface (see region marked 1 in Fig. 5). It was found that when two components have a significant mobility difference, typically greater than a 2:1 ratio, the effect of mobility is more important than the combined effect of misfit strain and solute
Y. Peng et al. / Computational Materials Science 123 (2016) 65–69
69
Fig. 5. Plot of epitaxial growth on a substrate of eight atomic layers in thickness showing dislocations nucleation and phase separation between peaks and valleys. The left column shows the atomic number density difference. The right column explains the local concentration difference in which dark regions represent the component A, and light regions represent the component B. In this figure, MA ¼ 1; M B ¼ 0:25, misfit strain is 0.04. Time corresponds to 3600Dt. The dislocations are enclosed by black squares, circles and arrows.
strain in determining which component accumulates at the surface. Under compressive stress films, peaks correspond to tensile stress regions, valleys correspond to regions of compressive stress. Therefore, in the presence of the unequal mobilities and the misfit strain, the components with greater mobility are driven toward regions of compressive stress which corresponds to valleys (see regions marked in 2 in Fig. 5). Further, while the strain in valleys increases to a certain number, there exists dislocations nucleation again (see regions marked squares in Fig. 5). 4. Conclusions The binary phase field crystal is employed to study influences of atomic sizes and mobility differences on the liquid phase heteroepitaxial. The results are as follows: (1) Due to the addition of components, different atomic sizes have important impacts on the formation of islands. It was found that large size atoms are driven toward regions of tensile stress which correspond to peaks but small size atoms are driven toward regions of compressive stress which correspond to valleys under positive misfit. (2) As for negative misfit stress, the film is under tensile stress. Large size atoms are driven toward regions of tensile stress which correspond to valleys but small size atoms are driven toward regions of compressive stress which correspond to peaks. (3) Due to the vertical separation and lateral separation, the epitaxial layer exhibits double islands phenomenon: light and dark islands. (4) In the presence of the different mobilities, the component with greater mobility accumulates at the film surface, what’s more, while two components have a significant mobility difference, the effect of mobility is more important than the combined effect of the misfit strain and solute strain in determining which component accumulates at the surface. (5) During the growth process, there are dislocations nucleation in valleys where the strain is highest. Acknowledgements This work was supported by the National Natural Science Foundation of China (Grant Nos. 51575452, 51475378 and 51474176), the Fundamental Research Funds for the Central Universities (Grant No. 3102015ZY025), Research Fund of the State Key Laboratory of Solidification Processing (NWPU) (Grant No. 161-QP-2016), NSFC—Guangdong mutual funds (phase ii) supercomputing science and applied research under special funding and National Supercomputing Center in Guangzhou. Appendix A. Supplementary material Supplementary data associated with this article can be found, in the online version, at http://dx.doi.org/10.1016/j.commatsci.2016. 06.017.
References [1] R. Zhu, E. Pan, P.W. Chung, Fast multiscale kinetic Monte Carlo simulations of three-dimensional self-assembled quantum dot islands, Phys. Rev. B 75 (20) (2007) 205339. [2] B.J. Spencer, M. Blanariu, Shape and composition map of a prepyramid quantum dot, Phys. Rev. Lett. 95 (20) (2005) 206101. [3] W.X. Li, J. Zhang, T.H. Shen, G.A. Jones, Magnetic nanowires fabricated by anodic aluminum oxide template—a brief review, Sci. China Phys. Mech. Astron. 54 (7) (2011) 1181–1189. [4] D. Chen, S. Xiong, S.H. Ran, B. Lin, One-dimensional iron oxides nanostructures, Sci. China Phys. Mech. Astron. 54 (7) (2011) 1190–1199. [5] M.L. Lee, C.W. Leitz, E.A. Fitzgerald, Formation of planar strained layers, 2004, Google Patents. [6] B.J. Spencer, P.W. Voorhees, S.H. Davis, Morphological instability in epitaxially strained dislocation-free solid films, Phys. Rev. Lett. 67 (26) (1991) 3696. [7] D.E. Jesson, K.M. Chen, S.J. Pennycook, T. Thundat, R.J. Warmack, Morphological evolution of strained films by cooperative nucleation, Phys. Rev. Lett. 77 (7) (1996) 1330. [8] R.J. Asaro, W.A. Tiller, Interface morphology development during stress corrosion cracking: Part I. Via surface diffusion, Metall. Trans. 3 (7) (1972) 1789–1796. [9] M.A. Grinfeld, The stress driven instability in elastic crystals: mathematical models and physical manifestations, J. Nonlinear Sci. 3 (1) (1993) 35–83. [10] C. Chen, Z. Chen, J. Zhang, X.J. Du, Phase-field crystal modeling of shape transition of strained islands in heteroepitaxy, Sci. China Phys. Mech. Astron. 55 (11) (2012) 2042–2048. [11] C. Chen, Z. Chen, J. Zhang, T. Yang, X.J. Du, Faceting transitions in crystal growth and heteroepitaxial growth in the anisotropic phase-field crystal model, Chin. Phys. B 21 (11) (2012) 118103. [12] K.R. Elder, M. Katakowski, M. Haataja, M. Grant, Modeling elasticity in crystal growth, Phys. Rev. Lett. 88 (24) (2002) 245701. [13] K.R. Elder, Z.F. Huang, A phase field crystal study of epitaxial island formation on nanomembranes, J. Phys.: Condens. Matter 22 (36) (2010) 364103. [14] M. Herrera, D. Gonzalez, M. Hopkinson, H.Y. Liu, R. Garacia, Influence of the Growth Temperature on the Composition Fluctuations of GaInNAs/GaAs Quantum Wells, Springer, Berlin Heidelberg, 2008, pp. 199–221. [15] A. Cartalade, A. Younsi, E. Regnier, S. Schuller, Simulations of phase-field models for crystal growth and phase separation, Procedia Mater. Sci. 7 (2014) 72–78. [16] J. Novák, S. Hasenohrl, I. Vavra, M. Kucera, Influence of surface strain on the MOVPE growth of InGaP epitaxial layers, Appl. Phys. A 87 (3) (2007) 511–516. [17] K.R. Elder, N. Provatas, J. Berry, P. Stefanovic, M. Grant, Phase-field crystal modeling and classical density functional theory of freezing, Phys. Rev. B 75 (6) (2007) 064107. [18] Z. Huang, K.R. Elder, N. Provatas, Phase-field-crystal dynamics for binary systems: Derivation from dynamical density functional theory, amplitude equation formalism, and applications to alloy heterostructures, Phys. Rev. E 82 (2) (2010) 021605. [19] K.R. Elder, H. Zhi-Feng, P. Nikolas, Amplitude expansion of the binary phasefield-crystal model, Phys. Rev. E: Stat., Nonlin, Soft Matter Phys. 81 (1Pt1) (2010) 132–135. [20] G. Tegze, G. Bansel, Gl. Toth, T. Pusztai, Z. Fanet, Advanced operator splittingbased semi-implicit spectral method to solve the binary phase-field crystal equations with variable coefficients, J. Comput. Phys. 228 (5) (2009) 1612– 1623. [21] B.J. Spencer, P.W. Voorhees, J. Tersoff, Enhanced instability of strained alloy films due to compositional stresses, Phys. Rev. Lett. 84 (11) (2000) 2449–2452. [22] Z.F. Huang, R.C. Desai, Epitaxial growth in dislocation-free strained alloy films: morphological and compositional instabilities, Phys. Rev. B: Condens. Matter 65 (20) (2002) 205419. [23] A.G. Norman, S.P. Ahrenkiel, H. Moutinho, M.M. Al-Jassim, Strain-dependent morphology of spontaneous lateral composition modulations in (AlAs){subm} (InAs){subn} short-period superlattices grown by molecular beam epitaxy, Appl. Phys. Lett. 73 (13) (1998) 1844–1846.
Contents lists available at ScienceDirect
Computational Materials Science journal homepage: www.elsevier.com/locate/commatsci
A binary phase field crystal study for phase segregation of liquid phase heteroepitaxial growth Yingying Peng ⇑, Yanli Lu ⇑, Zheng Chen, Genggeng Yu State Key Lab of Solidification Processing, Northwestern Polytechnical University, 710072 Xi’an, PR China
a r t i c l e
i n f o
Article history: Received 26 March 2016 Received in revised form 16 June 2016 Accepted 18 June 2016
Keywords: Binary PFC Atomic sizes Mobility differences Liquid phase heteroepitaxial
a b s t r a c t The binary phase field crystal (PFC) model is employed to investigate influences of atomic sizes and mobility differences on the liquid phase heteroepitaxial. It was found that large size atoms are driven toward regions of tensile stress which correspond to peaks in a compressively strained film but to valleys in a film with tensile strain. Small size atoms are on the contrary in contrast for large size atoms. Due to the existence of vertical separation and lateral separation resulting from atomic size differences, the epitaxial layer exhibits double-island phenomenon: light and dark islands. In the presence of different mobilities, atoms with greater mobility accumulate at the film surface. And in the process of epitaxial growth, there are misfit dislocations nucleating in valleys where the strain is highest. Ó 2016 Elsevier B.V. All rights reserved.
1. Introduction Due to unique optical and electronic properties, with potential applications in optoelectronics and semiconductor devices, self-assembled of nanostructure in semiconductor heteroepitaxy has been intensively studied [1–4]. These structures include junctions, quantum dot islands and superlattices in which planar interfaces are highly desire [5–7]. However, the difference of lattice constant between the deposited material and substrate results in mismatch strain. Morphological instability induced by the strain leads to nonplanar surfaces and defects in the film. This stressdriven morphological instability [8,9] has been studied extensively, and it is well established that the morphological instability can change the growth mode of a thin film from a layer-by-layer mechanism to a three-dimensional island mechanism. However, previous works mainly focus on the study of pure materials [10–12], with the development of alloy thin films, more and more investigations about alloy heteroepitaxy have been reported. A latest report on alloy heteroepitaxy from Elder [13] et al has shown influences of substrate thickness on the growth of the single island and multiple islands. Compared with pure materials, alloy thin films have more complicated compositions. As for epitaxial growth of alloy thin films, alloy components may be prone to phase separation and according to Vegard’s law, the lattice constant of the film is generally a function of alloy compositions. Therefore, the growth of alloy thin films needs further study. Recently, alloy phase ⇑ Corresponding authors. E-mail addresses: [email protected] (Y. Peng), [email protected] (Y. Lu). http://dx.doi.org/10.1016/j.commatsci.2016.06.017 0927-0256/Ó 2016 Elsevier B.V. All rights reserved.
separation [14–16] has captured many scholars’ attentions. The main issues motivated us are that how atomic sizes and mobility differences affect processes of phase separation and further influence the stability of epitaxial thin films. During the actual production process of thin films, it is very difficult to control the thickness of thin films because of the fast growth rate. In order to get nanoscale films, we find another way known as phase field crystal model to predict the growth of thin films. Liquid phase epitaxy as a mature technology has been used in the production of semiconductor optoelectronic devices, as well as magnetic garnets, superconductors, ferroelectrics and other optical materials. Therefore, the purpose of this paper is to illustrate how the two-dimensional binary PFC model [17,18] addresses compositional effects in alloy liquid phase heteroepitaxy, focusing on influences of atom size differences and mobility differences on phase separation. 2. Simulation method 2.1. Model The binary phase field crystal model proposed by Elder [19] recently has not only incorporated elasticity and plasticity on atomic length and diffusive timescales and multiple crystal orientations which are caused by the periodic oscillations of the crystal lattice, but also involves the solidification, phase segregation, solute diffusion and other important features. For a binary alloy made up of A and B atoms, to define the total number density q ¼ qA þ qB . To simplify calculations it is convenient to fist
66
Y. Peng et al. / Computational Materials Science 123 (2016) 65–69
, introduce the following dimensionless fields: nA ¼ ðqA qA Þ=q . Also, it indicates that the atomic number density nB ¼ ðqB qB Þ=q field n and the concentration field u (To show the atomic density difference, it can be approximated as the concentration field), which gives
dynamics of the fields n and u decouple, the respective equations of motions have the form:
@n dF ¼ M r2 @t dn
ð5Þ
Þ=q n ¼ nA þ nB ¼ ðq q u ¼ ðnB nA Þ þ ðq B q A Þ ¼ ðqB qA Þ=q
@u dF ¼ M r2 du @t
ð6Þ
ð1Þ
The numerical discretization of Eqs. (5) and (6) is implemented using the mixed explicit-implicit Fourier transformations and operator splitting by Tegze et al. [20].
where F is the free energy functional given by
Z F¼
n t v w u K ! dx ^ n n3 þ n4 þ cu þ u2 þ u4 þ j r u2 j 4 2 3 2 4 2 ð2Þ
where ^ Bl0 þ Bl2 u2 þ Bs0 ð2r2 þ r4 Þ þ 4auBs0 ðr2 þ r4 Þ. Briefly outlining the physical meaning of the model parameters is as follows: the Bl0 , Bl2 , Bs0 , t, v , w, K, c, a and u are constants, whose values depend on the material parameters. t, c and v are model parameters associated with the Taylor coefficient. The length scale of phase separation is decided by the interplay of u, K and w. Bl2 determines liquid volume modulus with the changing of atomic concentration field. Bl0 is the liquid temperature. Bs0 is associated with the solid phase of the elastic constants. DB ¼ Bl0 Bs0 corresponds to temperature expressing the initial driving force of crystallization, which can be enhanced by lowering DB. a represents the solute diffusion factor measuring the size difference between A and B atoms, and is defined as a ¼ ðaB aA Þ=a0 where aA ; aB ; and a0 are the relaxed-state lattice parameters for pure A, pure B and the alloy film. Also, in the spirit of keeping calculations as simple as possible without losing the basic physics contained in the model, c ¼ 0 in the free energy. Minimizing the free energy, h pffiffiffi pffiffiffi i n ¼ A 1=2 cosð2qy= 3Þ cosðqxÞ cosðqy= 3Þ will be used. Substituting this expression into (2) and minimizing with respect to q and A, recalling that u is assumed constant over the scale that n pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffi varies, getting qeq ¼ 3=ð2RÞ and Amin ¼ 4ðt þ t 2 15v DBÞ=ð15v Þ. To simulate microstructure formation in binary alloys, dynamical equations of motion for the fields n and u need to be developed. @n dF dF @nA @nB dF dF þ M B r2 ¼ þ ¼ M A r2 ¼ r M1 r þ r M2 r @t @t @t dn du dnA dnB
2.2. The chosen of parameters The equilibrium diagram is shown in Fig. 1 for the appropriate parameters set as follows: Bs0 ¼ 1:0, Bl2 ¼ 1:8, t ¼ 0:6, L ¼ 2:65, l ¼ 4:0, a ¼ 0:26 and w ¼ 0:088. Our simulation is conducted on a periodic system of size 512Dx 512Dx, the binary alloy film with average density difference u ¼ 0 grows from a liquid phase above the bulk coherent spinodal temperature, the time step Dt ¼ 0:5. According to the equilibrium diagram, we choose DB0 ¼ 0:00886 to implement all of the following simulations, growth at this DB0 above the miscibility gap is typical of experimental conditions and should ensure that phase separation is driven by local stresses and is not due to spinodal decomposition. Initial conditions consisted of a binary unstrained crystalline planar substrate of eight atomic layer thickness in the bottom of the simulation area which is below a symmetric supercooled liquid of components A and B. In what follows the misfit strain is defined as ðafilm asub Þ=asub , for a symmetric mixture of A and B components afilm ¼ ðaA þ aB Þ=2. 3. Results and discussion For the reason that the evident contrast as bright and dark regions can highlight the phase separation in the alloy, grey-scale maps are adopted. The results are presented in Fig. 2 for two components with different atomic sizes and same mobility. Fig. 2 (a) and (c) show the atomic number density difference while Fig. 2(b) and (d) demonstrate the local concentration difference. As shown in Fig. 2(a) and (b), at the initial stage of growth, epitaxial
0.10
¼ ðM A þ M B Þr2 ð^n tn þ mn3 þ cÞ þ ðM B M A Þr2 K r2 u
0.08
þ ðM A M B Þr2 ðn2 Bl2 u þ 2an2 BS0 ðr2 þ r4 Þ þ wu þ lu3 Þ ð3Þ
0.06
@u dF dF @nA @nB dF dF M B r2 ¼ ¼ r M2 r ¼ M A r2 þ r M1 r @t @t @t dn du dnA dnB
0.04
þ ðM A þ M B Þr
2
ðn2 Bl2
2
3
n2 BS0 ð
u þ 2a
ΔΒ0
¼ ðM A M B Þr ð^n tn þ mn þ cÞ ðM B þ M A Þr K r u 2
2
r þ r Þ þ w u þ lu Þ 2
4
0.02 0.00
3
ð4Þ
Liquid
ls ΔΒ 0 =0.0213 c ΔB0 =-0.0139
Solid (0,0.00886)
-0.02 -0.04
where M1 ¼ ðMA þ MB Þ=q2 and M2 ¼ ðM B MA Þ=q2 . To get Eqs. (3) and (4), we used the approximate value of Eq. (1): n ¼ nA þ nB , u ¼ nA nB . Eqs. (3) and (4) couple the dynamics of the fields n and u through a symmetric mobility tensor. The dependence of mobility M A and MB in general depends on local crystal density and the local relative concentration of species A and B. While M A and MB are different, a simple Euler algorithm was used for the time derivative and the spherical Laplacian approximation was used. Assuming a substitutional diffusion between species A and B, that the same M mobility applies for the two species, furthermore, considering that the mobility coefficient is a constant, the
-0.06 -0.08 -0.10 -0.25 -0.20 -0.15 -0.10 -0.05 0.00 ϕ
0.05
0.10
0.15
0.20
0.25
Fig. 1. Phase diagram as a function of DB0 u for two-dimensional triangular system. The blue solid line corresponds to liquid phase, the red solid line corresponds to solid phase, the black square dot corresponds to selected initial condition. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
Y. Peng et al. / Computational Materials Science 123 (2016) 65–69
(a)
(b)
(c)
(d)
67
Fig. 2. Plot of epitaxial growth on a substrate of eight atomic layers in thickness. (a) and (c) show the atomic number density difference. (b) and (d) explain the local concentration difference, in which light regions represent component A which consists of large size atoms and dark regions represent component B which consists of small size atoms. In this figure, M A ¼ M B , misfit strain is 0.088. Times from top to bottom are 80Dt and 1040Dt. Inset: the local concentration difference with 1800-time step. The start and the end of arrows correspond to positions of dislocations nucleation and nucleated dislocations, respectively.
films become unstable and nonplanar. Even without compositional effects, the stress nonuniformity from misfit between substrate and films is sufficient to drive the surface instability as in a single-component film. The addition of the alloy components increases the instability of thin films. Composition variations cause an enhancement of the instability in two ways. First, composition variations increase the stress nonuniformity, causing a larger driving force for the surface instability. Second, the stress nonuniformity results in the enhanced diffusion of a kind of atomic species to surface peaks, causing the instability to grow faster. As presented in Fig. 2(d), it is clear to see that there exists an obvious difference compared to pure epitaxial film. The epitaxial layer exhibits double islands phenomenon: light and dark islands. It fully shows alloy segregation instability couples with the morphological instability, resulting in the simultaneous modulations of the surface profile and alloy composition during epitaxial growth. The strain modulations associated with morphological perturbations may enhance the segregation of depositing alloy species. Reasons for double-island phenomenon are the presence of lateral and vertical separation due to different size components. Although the two components have equal mobility, they own different atomic sizes which result in solute strain. Also, the film is in compressive stress because of the positive misfit strain, and peaks correspond to tensile stress regions, valleys correspond to regions of compressive stress. Therefore, in the presence of the solute strain and the misfit strain, large size atoms are driven toward regions of tensile stress which correspond to peaks and small size atoms are driven toward regions of compressive stress which correspond to valleys in a compressively strained film. With the time going by, large size atoms gathered in the regions of peaks while small atoms gathered in the regions of valleys, that are double islands. What is more, in addition to the vertical separation, there exists lateral phase separation. This coupling creates a lateral phase separation on the length scale of the surface instability and has been predicted as well as verified for binary films [21,22]. The epitaxial system can maintain a lower free energy by transferring
atoms from the valleys to peaks, because the transition leads to a decrease in the contact area between the substrate and a new layer. Thus, atomic transition to upper layers leads to a relaxation in the local strain field. Lateral segregation has two important implications: on the one hand, if epitaxial films thickness had been increased during growth, misfit dislocations would have nucleated at the interface, most probably in valleys where the strain is highest, as shown in the inset of Fig. 2(d); on the other hand, the lateral segregation reduces large size atoms in the regions of the highest strain, thereby reducing the strain compared to a similar morphology without segregation, and hence, retarding dislocation nucleation. Thus, compared with single component, not all islands coalesce with each other and only a few dislocations nucleation in the structure of double islands. What’s more, due to the lateral segregation, light islands grow much quicker than dark islands. By comparing Fig. 2 with Fig. 3, it is clear to find that the sign of misfit strain is the only difference between them. Similarly, the negative misfit strain has some influences on the epitaxial growth. Due to the negative misfit strain, the film is under tensile stress, and peaks correspond to compressive stress regions, valleys correspond to regions of tensile stress. Therefore, in the presence of the solute strain and the misfit strain, large size atoms are driven toward regions of tensile stress which correspond to valleys and small size atoms are driven toward regions of compressive stress which correspond to peaks in a tensile strained film. Accompanied with the epitaxial growth, large size atoms gathered in regions of valleys and small size atoms gathered in regions of peaks, and then double islands are formed: dark and light islands. Surely, in addition to the vertical separation, there exists lateral phase separation. For the same purpose, in order to lower free energy, atoms move from valleys to peaks, the number of small size atoms in valley regions becomes much less. Meanwhile, there are misfit dislocations nucleating in valleys where the strain is highest. Similarly, not all islands coalesce with each other and only a few dislocations nucleation in the structure of double islands as shown in the inset of Fig. 3. What’s more, by comparing the inset of Fig. 2 with that of
68
Y. Peng et al. / Computational Materials Science 123 (2016) 65–69
(a)
(b)
(c)
(d)
Fig. 3. Plot of epitaxial growth on a substrate of eight atomic layers in thickness. (a) and (c) show the atomic number density difference. (b) and (d) explain the local concentration difference, in which light regions represent component A which consists of large size atoms and dark regions represent component B which consists of small size atoms. In this figure, M A ¼ M B , misfit strain is 0.088. Times from top to bottom are 80Dt and 600Dt. Inset: the local concentration difference with 840Dt, dislocations are enclosed by black circles.
Fig. 3, under negative misfit strain, misfit dislocations are earlier to nucleate, therefore, under the condition of positive misfit strain, double islands are much more stable and need much more time to coalesce. Above all, the mechanism that large atoms preferentially segregate to peaks if the misfit is compressive but to valleys if the misfit is tensile is consistent with experimental evidence in the AlAs/InAs system [23]. In order to describe influences of atomic sizes on the phase separation in detail, maximum thickness of atomic layer is shown as a function of time in Fig. 4 by combining with Figs. 2 and 3. In Fig. 4, curve1 and curve 2 exhibit the atomic layer thickness of small size atoms and large size atoms under the positive misfit strain, respectively. Curve 3 and curve 4 show the atomic layer
1
atomic layer thickness
25
2 3 4
20
15
10
5
0
400
800
1200
1600
2000
time/Δt Fig. 4. Maximum thickness of atomic layer is plotting as a function of time. Curve 1 and curve 2 exhibit the atomic layer thickness of small size atoms and large size atoms under the positive misfit strain, respectively. Curve 3 and curve 4 show the atomic layer thickness of large size atoms and small size atoms under the negative misfit strain, respectively.
thickness of large size atoms and small size atoms under the negative misfit strain, respectively. While the misfit strain is positive, large size atoms are driven toward positions of peaks and small size atoms are driven to valleys. Finally, large size atoms gathered in regions of peaks and smaller size atoms gathered in regions of valleys. In this way, double islands are formed. Also, it is shown that the atomic layer thickness of large size atoms is much larger than that of small size atoms. Likely, while the misfit strain is negative, double islands are formed in the way that large size atoms gathered in regions of valleys and small size atoms gathered in regions of peaks, which is on the contrary compared with double islands under positive misfit strain. In addition, it is clear to see that the difference of the atomic layer thickness of double islands in positive misfit strain is larger than that under negative misfit strain. This phenomenon in another way proves that under the condition of positive misfit strain, the double islands are much more stable and need much more time to coalesce. Fig. 5 shows influences of unequal mobilities on the epitaxial growth, the left column shows the atomic number density difference and the right column exhibits the local concentration difference. Dislocations are enclosed by black squares, circles and arrows. Light numbers show positions of phase separation. The misfit strain can induce this morphological instability because the stressed material with perturbed surface has the nonuniform strain energy. The stress in the film with sinusoidally perturbed surface is relaxed at the wave peaks and increased at the troughs by equal magnitudes. Also, the misfit strain is positive, the film is in compressive stress. Therefore, at the initial stage of growth, misfit dislocations nucleated at the interface in troughs where the compressive strain is highest (see regions marked circles and arrows in Fig. 5). And the atoms with greater mobility accumulate at the film surface (see region marked 1 in Fig. 5). It was found that when two components have a significant mobility difference, typically greater than a 2:1 ratio, the effect of mobility is more important than the combined effect of misfit strain and solute
Y. Peng et al. / Computational Materials Science 123 (2016) 65–69
69
Fig. 5. Plot of epitaxial growth on a substrate of eight atomic layers in thickness showing dislocations nucleation and phase separation between peaks and valleys. The left column shows the atomic number density difference. The right column explains the local concentration difference in which dark regions represent the component A, and light regions represent the component B. In this figure, MA ¼ 1; M B ¼ 0:25, misfit strain is 0.04. Time corresponds to 3600Dt. The dislocations are enclosed by black squares, circles and arrows.
strain in determining which component accumulates at the surface. Under compressive stress films, peaks correspond to tensile stress regions, valleys correspond to regions of compressive stress. Therefore, in the presence of the unequal mobilities and the misfit strain, the components with greater mobility are driven toward regions of compressive stress which corresponds to valleys (see regions marked in 2 in Fig. 5). Further, while the strain in valleys increases to a certain number, there exists dislocations nucleation again (see regions marked squares in Fig. 5). 4. Conclusions The binary phase field crystal is employed to study influences of atomic sizes and mobility differences on the liquid phase heteroepitaxial. The results are as follows: (1) Due to the addition of components, different atomic sizes have important impacts on the formation of islands. It was found that large size atoms are driven toward regions of tensile stress which correspond to peaks but small size atoms are driven toward regions of compressive stress which correspond to valleys under positive misfit. (2) As for negative misfit stress, the film is under tensile stress. Large size atoms are driven toward regions of tensile stress which correspond to valleys but small size atoms are driven toward regions of compressive stress which correspond to peaks. (3) Due to the vertical separation and lateral separation, the epitaxial layer exhibits double islands phenomenon: light and dark islands. (4) In the presence of the different mobilities, the component with greater mobility accumulates at the film surface, what’s more, while two components have a significant mobility difference, the effect of mobility is more important than the combined effect of the misfit strain and solute strain in determining which component accumulates at the surface. (5) During the growth process, there are dislocations nucleation in valleys where the strain is highest. Acknowledgements This work was supported by the National Natural Science Foundation of China (Grant Nos. 51575452, 51475378 and 51474176), the Fundamental Research Funds for the Central Universities (Grant No. 3102015ZY025), Research Fund of the State Key Laboratory of Solidification Processing (NWPU) (Grant No. 161-QP-2016), NSFC—Guangdong mutual funds (phase ii) supercomputing science and applied research under special funding and National Supercomputing Center in Guangzhou. Appendix A. Supplementary material Supplementary data associated with this article can be found, in the online version, at http://dx.doi.org/10.1016/j.commatsci.2016. 06.017.
References [1] R. Zhu, E. Pan, P.W. Chung, Fast multiscale kinetic Monte Carlo simulations of three-dimensional self-assembled quantum dot islands, Phys. Rev. B 75 (20) (2007) 205339. [2] B.J. Spencer, M. Blanariu, Shape and composition map of a prepyramid quantum dot, Phys. Rev. Lett. 95 (20) (2005) 206101. [3] W.X. Li, J. Zhang, T.H. Shen, G.A. Jones, Magnetic nanowires fabricated by anodic aluminum oxide template—a brief review, Sci. China Phys. Mech. Astron. 54 (7) (2011) 1181–1189. [4] D. Chen, S. Xiong, S.H. Ran, B. Lin, One-dimensional iron oxides nanostructures, Sci. China Phys. Mech. Astron. 54 (7) (2011) 1190–1199. [5] M.L. Lee, C.W. Leitz, E.A. Fitzgerald, Formation of planar strained layers, 2004, Google Patents. [6] B.J. Spencer, P.W. Voorhees, S.H. Davis, Morphological instability in epitaxially strained dislocation-free solid films, Phys. Rev. Lett. 67 (26) (1991) 3696. [7] D.E. Jesson, K.M. Chen, S.J. Pennycook, T. Thundat, R.J. Warmack, Morphological evolution of strained films by cooperative nucleation, Phys. Rev. Lett. 77 (7) (1996) 1330. [8] R.J. Asaro, W.A. Tiller, Interface morphology development during stress corrosion cracking: Part I. Via surface diffusion, Metall. Trans. 3 (7) (1972) 1789–1796. [9] M.A. Grinfeld, The stress driven instability in elastic crystals: mathematical models and physical manifestations, J. Nonlinear Sci. 3 (1) (1993) 35–83. [10] C. Chen, Z. Chen, J. Zhang, X.J. Du, Phase-field crystal modeling of shape transition of strained islands in heteroepitaxy, Sci. China Phys. Mech. Astron. 55 (11) (2012) 2042–2048. [11] C. Chen, Z. Chen, J. Zhang, T. Yang, X.J. Du, Faceting transitions in crystal growth and heteroepitaxial growth in the anisotropic phase-field crystal model, Chin. Phys. B 21 (11) (2012) 118103. [12] K.R. Elder, M. Katakowski, M. Haataja, M. Grant, Modeling elasticity in crystal growth, Phys. Rev. Lett. 88 (24) (2002) 245701. [13] K.R. Elder, Z.F. Huang, A phase field crystal study of epitaxial island formation on nanomembranes, J. Phys.: Condens. Matter 22 (36) (2010) 364103. [14] M. Herrera, D. Gonzalez, M. Hopkinson, H.Y. Liu, R. Garacia, Influence of the Growth Temperature on the Composition Fluctuations of GaInNAs/GaAs Quantum Wells, Springer, Berlin Heidelberg, 2008, pp. 199–221. [15] A. Cartalade, A. Younsi, E. Regnier, S. Schuller, Simulations of phase-field models for crystal growth and phase separation, Procedia Mater. Sci. 7 (2014) 72–78. [16] J. Novák, S. Hasenohrl, I. Vavra, M. Kucera, Influence of surface strain on the MOVPE growth of InGaP epitaxial layers, Appl. Phys. A 87 (3) (2007) 511–516. [17] K.R. Elder, N. Provatas, J. Berry, P. Stefanovic, M. Grant, Phase-field crystal modeling and classical density functional theory of freezing, Phys. Rev. B 75 (6) (2007) 064107. [18] Z. Huang, K.R. Elder, N. Provatas, Phase-field-crystal dynamics for binary systems: Derivation from dynamical density functional theory, amplitude equation formalism, and applications to alloy heterostructures, Phys. Rev. E 82 (2) (2010) 021605. [19] K.R. Elder, H. Zhi-Feng, P. Nikolas, Amplitude expansion of the binary phasefield-crystal model, Phys. Rev. E: Stat., Nonlin, Soft Matter Phys. 81 (1Pt1) (2010) 132–135. [20] G. Tegze, G. Bansel, Gl. Toth, T. Pusztai, Z. Fanet, Advanced operator splittingbased semi-implicit spectral method to solve the binary phase-field crystal equations with variable coefficients, J. Comput. Phys. 228 (5) (2009) 1612– 1623. [21] B.J. Spencer, P.W. Voorhees, J. Tersoff, Enhanced instability of strained alloy films due to compositional stresses, Phys. Rev. Lett. 84 (11) (2000) 2449–2452. [22] Z.F. Huang, R.C. Desai, Epitaxial growth in dislocation-free strained alloy films: morphological and compositional instabilities, Phys. Rev. B: Condens. Matter 65 (20) (2002) 205419. [23] A.G. Norman, S.P. Ahrenkiel, H. Moutinho, M.M. Al-Jassim, Strain-dependent morphology of spontaneous lateral composition modulations in (AlAs){subm} (InAs){subn} short-period superlattices grown by molecular beam epitaxy, Appl. Phys. Lett. 73 (13) (1998) 1844–1846.