PHYSICAL REVIEW B VOLUME 56, NUMBER 9 1 SEPTEMBER 1997-I Effect of the vacancy interaction on antiphase domain growth in a two-dimensional binary alloy Marcel Porta, Carlos Frontera, Eduard Vives, and Teresa CastaŽn Departament d'Estructura i Constituents de la Mate ria, Facultat de FiŽsica, Universitat de Barcelona, Diagonal 647, E-08028 Barcelona, Catalonia, Spain Received 26 March 1997; revised manuscript received 16 May 1997 We have performed a Monte Carlo simulation of the influence of diffusing vacancies on the antiphase domain growth process in a binary alloy after a quench through an order-disorder transition. The problem has been modeled by means of a Blume-Emery-Griffiths Hamiltonian whose biquadratic coupling parameter K controls the microscopic interactions between vacancies. The asymmetric term L has been taken as L 0 and the ordering dynamics has been studied at very low temperature as a function of K inside the range 0.5 K/J 1.40 with J 0 being the ordering energy . The system evolves according to the Kawasaki dynamics so that the alloy concentration is conserved while the order parameter is not. The simulations have been performed on a two-dimensional square lattice and the concentration has been taken so that the system corresponds to a stoichiometric alloy with a small concentration of vacancies. We find that, independently of K, the vacancies exhibit a tendency to concentrate at the antiphase boundaries. This effect gives rise, via the vacancy-vacancy interaction described by K), to an effective interaction between bulk diffusing vacancies and moving interfaces that turns out to strongly influence the domain growth process. One distinguishes three different behaviors: i For K/J 1 the growth process of ordered domains is anisotropic and can be described by algebraic laws with effective exponents lower than 1/2; ii K/J 1 corresponds to the standard Allen-Cahn growth; iii for K/J 1 we found that, although the motion of the interface is curvature driven, the repulsive effective interaction between both the vacancies in the bulk and those at the interfaces slows down the growth. S0163-1829 97 06933-6 I. INTRODUCTION merical simulation studies12­15 provide evidences for slow growth, either logarithmic growth laws or algebraic laws The dynamical evolution of a binary alloy after a quench with small exponents. from a high-temperature disordered phase has been one of A general feature commonly observed during the early- the prototypes in the study of relaxational processes towards time evolution in systems with annealed defects is the ten- equilibrium. It has been found that the late stages of this dency to concentrate the disorder at the domain walls. In- process obey dynamical scaling and the typical domain size deed, experiments10,16 and numerical simulations8,12­15,17 R(t) dominates all other lengths.1,2 In this regime, the do- reveal that the vacancies and the excess particles tend to mains grow in time according to a power law R(t) tx with accumulate at the domain walls. This is accompanied by a the growth exponent x satisfying a remarkable universality. depletion in the bulk defect concentration, which renders the For the case of a pure, ideal system, the exponent x 1/2 is excess internal energy unsuitable to measure the total associated with cases where the order parameter is not amount of interfaces. Only at late times, as the interfaces conserved,3 whereas x 1/3 describes systems with con- disappear and the system approaches equilibrium, the an- served order parameter.4 For the nonconserved case a typical nealed defects may dissolve again into the bulk, provided example is a binary alloy undergoing an order-disorder phase that they display no cooperative phenomena. Simultaneously transition. Most theoretical studies are limited to ideal con- an overshooting in the bulk order parameter is observed.12,18 ditions, that is, to a pure stoichiometric binary alloy. In such Very recently it has been suggested12,18,19 that this is a ge- conditions, theory3 and numerical simulations5­8 definitively neric effect in ordering dynamics, coming from a subtle agree about the value of the kinetic exponent x 1/2. Much competition between nonequilibrium internal energy and less unanimity is obtained from the experiments. This is cer- nonequilibrium entropy. tainly due to the imperfections always present in real mate- In previous13,14 studies of the diluted square Ising model rials. Some examples are vacancies, third-component with nearest-neighbor interactions, it was obtained that the impurities,9 nonstoichiometry,10 dislocations, etc. It is of effect of a small concentration of vacancies is dramatic, lead- great interest to elucidate in which manner and to what ex- ing to an extremely slow growth described by a logarithmic tent the presence of these imperfections modifies the ideal growth-law14 or even to a complete pinning of the process.13 asymptotic growth law. This difference, obtained on the same model in the limit of We shall not consider here the problem of quenched low vacancy concentration, comes from special features in- disorder,11 but concentrate on mobile punctual defects such troduced in the coupled dynamics used; that is, the system as vacancies, third-component impurities, and excess par- evolves according to the nonconserved Glauber dynamics ticles in off-stoichiometric binary alloys. These belong to the but the vacancy concentration is forced to be constant. We category commonly named annealed disorder. Despite the notice that in these studies, the vacancies exhibit a natural theoretical suggestion that such a kind of disorder should not tendency to cluster. The results obtained in the present work modify the asymptotic growth law,12 experiments,10 and nu- show that the true asymptotic growth behavior is definitively 0163-1829/97/56 9 /5261 10 /$10.00 56 5261 © 1997 The American Physical Society 5262 PORTA, FRONTERA, VIVES, AND CASTAŽN 56 algebraic with effective exponents smaller than the Allen- TABLE I. Bond energies and the same measured with respect to Cahn value. the A-B bond for the BEG model as a function of the parameters J, Concerning the present scenario for the effect of excess K, and L. particles in a nonstoichiometric binary alloy without vacan- cies, very recently it was suggested8 that the existence of Bond Energy Excess energy effective interactions between diffusing excess particles and A-B J K 0 those localized at the antiphase domain boundaries is crucial A-A J K 2L 2J 2L in determining the essential time dependence of the growth law. It was shown that when these specific interactions are B-B J K 2L 2J 2L not present the main assumptions underlying the Allen-Cahn A-V 0 J K theory are fulfilled. On the other hand there is experimental B-V 0 J K evidence that small deviations from the stoichiometric com- V-V 0 J K position may provoke drastic modifications in the growth law. It has been reported that the ordering kinetics in pair wise and restricted to nearest neighbors n.n. only. The Cu 10 0.79Au 0.21 shows a crossover from the standard Allen- spin-1 BEG Hamiltonian is Cahn growth law, for stoichiometric Cu3Au,10,20 to a loga- rithmic growth law. In Ref. 8 it was suggested that to ac- n.n. n.n. n.n. count for such behavior, additional interactions to the ones H J S 2 2 2 2 iS j K Si Sj L Si Sj SiSj , 1 present in the nearest neighbors Ising model are needed. i, j i, j i, j Indeed, the authors showed how this can be accomplished by where J stands for the atom-atom ordering interaction, K is a simply extending the interactions to next nearest neighbors. biquadratic coupling parameter accounting for the energy Our main interest here will be to incorporate the effect of difference between atom-atom pairs and those involving va- such interactions between annealed impurities in a more gen- cancies, and L is an asymmetry term accounting for the en- eral framework. ergy difference between A-A and B-B pairs. When the pa- The three-state Blume-Emery-Griffiths BEG model is rameter K promotes the formation of atom-atom pairs, the especially adequate for our purposes.21 For given values of vacancies tend to cluster, whereas if the vacancy-atom pairs the model parameters (K and L), the third value of the spin are preferred, cooperative effects for the vacancies are not variable may represent either a vacancy or an impurity in expected, at least in the limit of low vacancy concentration. particular an excess particle . In the present work we neglect Although the hypothesis of a rigid lattice is crude, since the asymmetric term (L 0) and restrict ourselves to the lattice deformations have a strong influence on the atom dy- case of vacancies in a stoichiometric binary alloy, so that the namics around the vacancies, this effect is partially taken additional coupling coming from the interplay between the into account in the phenomenological character of the con- diffusive motion of both vacancies and excess particles is not stant K. We also expect that the assumption of pair interac- considered here. tions only, disregarding many-body effects, only introduces The influence of mobile vacancies on the kinetics of or- quantitative but not qualitative changes concerning the order- dering arises from the coupling between their diffusive dy- ing dynamics. namics and the motion of the domain walls in this case, We have restricted the present study to concentrations so antiphase-domain boundaries APB . This intercoupling de- that the system corresponds to a stoichiometric binary alloy pends on the two following facts concerning the behavior of with small concentration of vacancies. Since we have taken the vacancies: their tendency to precipitate at the APBs and J 0, the ordered phase will be antiferromagneticlike, with their tendency to cluster. Whereas the former is encountered almost all the bonds of the A-B kind. for all values of the model parameters studied in this work, As we have mentioned in the Introduction our goal here is the interaction among vacancies, and furthermore their ten- to study the influence of annealed vacancies on the kinetics dency to cluster, depends on K. The combination of both of domain growth. This influence originates in the interplay effects gives rise to an effective interaction controlled by K) between the two following specific interactions: i the between bulk diffusing vacancies and those localized at the vacancy-APB interaction here the APB is thought as a se- interfaces that turns out to be crucial in determining the es- quence of A-A and B-B bonds and ii the vacancy-vacancy sential time dependence of the growth law. interaction. When the the vacancy-APB interaction is attrac- The organization of this paper is the following. We start tive, the vacancies tend to concentrate at the APBs and there- by defining the model and the region of parameters of inter- fore the vacancy-vacancy interaction introduces an effective est here Sec. II . In Sec. III we provide the details of the interaction between bulk vacancies and APBs. simulations and describe the algorithms used. In Sec. IV we The next step is to calculate these specific interactions for present the results and discuss them in Sec. V. Finally, in the different values of the model parameters. They are ob- Sec. VI we summarize our main conclusions. tained from the bond energies, expressed with respect to the energy of the A-B bond, summarized in Table I. For notation II. THE MODEL AND PARAMETERS we introduce the following reduced parameters K* K/J, L* L/J, with J 0. It follows that K* 1 separates the We assume an underlying rigid square lattice with tendency for the vacancies to cluster (K* 1) or not i 1, . . . ,N l l sites that can be occupied either by A at- (K* 1). Moreover, for L* 1 the specific interaction be- oms (Si 1), B atoms (Si 1), or vacancies (Si 0). Fol- tween vacancies and APBs is attractive, favoring the absorp- lowing standard procedures, the interactions are taken to be tion of vacancies at the interfaces. For values of L* 1 the 56 EFFECT OF THE VACANCY INTERACTION ON . . . 5263 FIG. 1. Regions in the space of the parameters K* and L* where different dynamics are expected. The solid line K* 1 sepa- rates the region of vacancy attraction (K* 1) from the region of vacancy repulsion (K* 1). The solid lines L* 1 separate the region where the vacancies tend to precipitate at the antiphase boundaries ( L* 1) from the region of vacancy-antiphase bound- ary repulsion ( L* 1). Solid circles indicate the points studied in FIG. 2. a kBT/J vs cV phase diagram of the BEG model for the present work. The square corresponds to the case K* L* 0 K* L* 0. b kBT/J vs K* phase diagram of the BEG model for and the diamonds are the points where the model is formally L* 0 and cV 0.06. The points correspond to the Monte Carlo equivalent to a nonstoichiometric binary alloy. results whereas the solid lines are obtained from mean-field calcu- lations. Dashed lines are just guides to the eyes. In the inset, we behavior is more complex, in particular, if K* 1 a compe- show the region of interest. The arrows indicate the working tem- perature and vacancy concentration tition between K* and the asymmetric term L* appears. The results of the above analysis are illustrated in Fig. 1. factor width. It is known that the structure factor provides an We have indicated, in white, the region where the vacancies overall description of the ordering process. In particular, the exhibit tendency to accumulate at the APBs. The line K* 1 study of its time evolution will provide information about the separates the regions with vacancy attraction and vacancy dynamical scaling properties. This second group of simula- repulsion. Black circles indicate the points studied in the tions constitutes the major part of results presented and, con- present work, all sitting along the line L* 0. The point at trarily to the equilibrium simulations, is very time consum- K* 0 and L* 0 indicated by a square has been previ- ing. ously studied in14 and corresponds to a diluted Ising model. The points with K* 1 and L* 1, indicated by dia- monds, correspond to a nonstoichiometric binary alloy with- A. Equilibrium simulations: Phase diagram out vacancies.8 In order to obtain the phase diagram of model 1 in the It is worth mentioning that the BEG model, with K 0, particular case of L* 0, we have performed Monte Carlo L 0, and J 0, has also been used for the study of the simulations in the grand canonical ensemble using the Leg- ordering dynamics via vacancies with cV 4 10 4. Such a endre transformation restricted dynamics, only allowing exchanges between atoms and vacancies, may strongly modify the dynamics.22 This mechanism is not considered in the present study. H 2 GC H Si , 2 III. MONTE CARLO SIMULATIONS where stands for the chemical potential difference between atoms either A or B) and vacancies. This is because we We have performed different Monte Carlo simulations of restrict our study to the case of stoichiometric composition the model defined in Sec. II with L* 0, J 0, and K* in- (NA NB) and the chemical potentials of A and B are then side the range 0.5 K* 1.4. First we need to know, in the equal. The simulations have been performed on a system of region of vacancy attraction (K* 1), the temperature at linear size l 128 using the Glauber dynamics implemented which the system separates into two phases. This is impor- into the Metropolis algorithm. The different runs are ex- tant in order to characterize the equilibrium state at the point tended up to 1500 Monte Carlo steps per site MCs , our unit to which quenches have been performed. Next, in order to of time. characterize the time evolution of the ordering process sub- Figures 2 a and 2 b show two different sections of the sequent to the quench, we focus on the study of the structure phase diagram obtained by Monte Carlo simulations. In the 5264 PORTA, FRONTERA, VIVES, AND CASTAŽN 56 same figure we also show the results obtained by using stan- have been averaged over equivalent directions. The size and dard mean-field techniques self-consistent field method .23 shape of the ordered domains has been obtained by fitting the This has been done for completeness and in order to reduce averaged profiles to a Lorentzian function powered to 3/2 in the range of model parameters to be explored numerically. order to reproduce the Porod's law for the decay of the tail at Figure 2 a temperature versus vacancy concentration cor- long q's:24 responds to a section with fixed K* 0, and Fig. 2 b tem- perature versus K*) to a section with fixed vacancy concen- 3/2 tration c S q,t a t 4 V 0.06. Phase I corresponds to an atomic disordered phase with randomly diluted vacancies; phase II 1 q/ t 2 corresponds to an atomic ordered AB phase with randomly where q is the distance to the superstructure peak diluted vacancies, which exhibit only short-range order; q k ( 1 phase III corresponds to a phase separation region with co- 2 , 12 ) , a ( t ) is the height of the peak, and ( t ) its existing ordered AB domains with low concentration of va- width. Only data with S(q,t) S0 2.5 10 5 has been con- cancies and vacancy clusters with a low concentration of sidered for the fits. S0 has been obtained from a completely disordered atoms. The inset in Fig 2 a shows an enlarged disordered system with the same concentration of particles portion of the coexistence region in the limit of very low and vacancies. The quantities a(t) and (t) provide infor- vacancy concentration. In this limit both phase transitions mation about the square order parameter growth and the in- (I II and II III) are very separated, so that they occur verse domains size, respectively. almost independently. This is reflected in the straight-line- 2. N-fold way algorithm shaped boundaries in Fig 2 b . One obtains that the order- disorder transition temperature (I II) is independent of K*, For the particular case of K* 0 and L* 0 we have used while the temperature for the phase separation transition the N-fold algorithm25 in order to reach very long times (107 (II III) depends linearly on K*. We expect that for higher MCs in the evolution of a system of linear size l 200. The values of c possible exchanges have been classified into 11 different V the coupling between both phase transitions makes the boundaries become curved. In any case notice classes, according to their energy change. A class including that, when K* 1, the coexistence region III disappears, those exchanges not modifying the system configuration has and that the AB ordered phase II extends down to T 0. also been taken into account. This is done in order to com- pare with the standard dynamical simulations in which such B. Nonequilibrium simulations exchanges are considered. The time elapsed after each ex- change has been taken as the average time needed, in a stan- Although most of the results, presented in the next sec- dard Monte Carlo simulation, for the acceptance of a useful tion, have been obtained following the standard Kawasaki proposed exchange. In this case, the structure factor evolu- dynamics, alternative optimized algorithms have been used tion has also been studied. The averages have been per- when specially long simulations were needed. This subsec- formed over 30 independent runs. tion is devoted to a description of the different algorithms used in the study of the time evolution of the process that 3. Optimized multigrid algorithm follows a thermal quench from very high temperature disor- dered phase to T 0.1J/k Due to the very large number of classes for K* 0 the B performed on a stoichiometric binary alloy with a small concentration of vacancies fixed at N-fold algorithm becomes difficult to construct. Therefore, a c simpler but less optimized algorithm has been constructed in V 0.06 being cA cB 0.47). The different values of K* studied correspond to final states into the ordered phases ei- order to reach very long times. Starting from a standard mul- ther II and III. tigrid algorithm8,14 we have made, for each sublattice, a list of the exchanges whose probability of being accepted is not 1. Standard dynamical simulations negligible. When a sublattice is chosen, only the exchanges present in the list are attempted. The method turns out to be These are simulations performed using the standard Me- very efficient when the number of attempted exchanges is tropolis algorithm together with the Kawasaki dynamics. The low. In particular, for the case K* 0.6 and L* 0, we have linear system size is l 200 even though some initial studies performed simulations up to 107 MCs with no major diffi- were performed on systems of l 100. Starting from an ini- culties the system size was l 200). tial disordered configuration, the runs have been typically extended up to 20 000 MCs. Moreover, averages over about 20 independent realizations have been performed. From each IV. DOMAIN GROWTH RESULTS simulation we have extracted the time evolution of the struc- Figure 3 shows snapshots of the microconfigurations as ture factor defined as they evolve with time after the quench for three different 2 values of the parameter K* 0.6, 1.0, and 1.4). The system S k 1 2 size is l 200 and the quench temperature T 0.1J/k N Siexp i , 3 B . In i 1,N a k r i this case we have used the standard dynamical algorithm described before. The vacancies are indicated in black while where k are the reciprocal space vectors, a is the lattice ordered regions are indicated in white. Due to the tendency parameter, and r i is the vector position of site i. We have exhibited by the vacancies to concentrate at the interfaces, focused on the profiles along the 10 and 11 directions the antiphase domain structure naturally shows up. The pro- around the superstructure peak at k ( 1 1 2 2 ) . Moreover, they cess of absorption of vacancies at the interfaces starts at very 56 EFFECT OF THE VACANCY INTERACTION ON . . . 5265 FIG. 3. Snapshots of the evolving domain structure for K* 0.6, K* 1.0, and K* 1.4. Vacancies (A and B particles are painted in black white . The simulations are performed in a 200 200 square lattice with cV 0.06 at T 0.1J/kB . early times and follows until the whole interface network For the case K* 1 the algebraic growth breaks down before saturates. This initial regime was studied previously8 for the 104 MCs due to finite-size effects. This is consistent with case of a nonstoichiometric binary alloy. Simultaneously, it the snapshots shown in Fig. 3. Indeed, K* 1 corresponds to was discussed in a more general context and suggested to be the fastest evolution and, thus, finite-size effects appear be- a generic effect in ordering dynamics.12,18,19 In any case, this fore. For the other two cases the algebraic regime extends for is a transient, prior to the long-time domain growth regime of much longer times, and finite-size effects are not found. We interest here. Nevertheless, we remark that this phenomenon obtain that while for K* 1 and K* 1.4 both (10) and of vacancy precipitation at the interfaces results in crucial (11) evolve with the same exponent, for K* 0.6 the two importance in the subsequent evolution dictated by interface exponents are clearly different. This is indicative of the ex- reduction specially when the interaction between vacancies istence of anisotropic growth and, as we shall discuss below, is switched on (K* 1). We now come back to Fig. 3. Clear it is related to the local accumulation of vacancies at the differences can be observed in relation to both the orientation of the the antiphase boundaries and the speed of the evolu- tion towards equilibrium. For K* 0.6 vacancy attraction , the domains appear squarelike with the interfaces preferably directed along the 10 direction. In the case K* 1.4 va- cancy repulsion , although the domains are also square shaped, the interfaces are directed along the 11 direction. Moreover, in both cases, the interfaces tend to be flat or almost flat , at least in the regime depicted in Fig. 3. As we shall discuss below, this is a consequence of the vacancy- vacancy interaction that introduces energy barriers for the motion of the vacancies localized at the interfaces favoring, in each case, the different orientation of the APBs. Contrar- ily, no preferred orientation for the boundaries is observed when K* 1. Remember that in this case there is no specific interaction between vacancies. Concerning the speed of the different evolutions, the fastest process occurs for K* 1, whereas for the other two cases it is clearly slower, appar- ently even more for K* 1.4. The introduction of specific interactions between vacancies seems to be behind the slower evolutions, although the underlying physics is, in both cases, different. Before proceeding further it is interesting to look at the quantitative results obtained from structure factor calcula- tions. In Fig. 4 we show the time evolution of the peak width (t) of the structure factor along the two relevant directions, (10) open circles and (11) filled circles for the same three FIG. 4. Width of the structure factor vs time for K* 0.6, selected values of K* as in Fig. 3. Dashed lines indicate the K* 1.0, and K* 1.4. Open circles correspond to the (10) direc- regimes of algebraic domain growth and the numbers on top tion and filled circles to the (11) direction. Dashed lines indicate the are the corresponding fitted values of the kinetic exponents. regions where the growth exponents, written on top, are fitted. 5266 PORTA, FRONTERA, VIVES, AND CASTAŽN 56 FIG. 6. Structure factor width vs time for K* 0. Open filled symbols correspond to the 10 11 direction. The end of the plateau is taken to be at the inflection point, as indicated by the FIG. 5. Time needed to reach the algebraic regime vs K* black arrow. The results shown by circles have been obtained following circles . Straight lines are the predicted slopes of these curves from standard simulations whereas the long time results, indicated by the energy barriers obtained from the model. squares, have been obtained by using the N-fold algorithm. The system size is 200 200, cV 0.06, and T 0.1J/kB . vertices of the interfaces see Fig. 3 . Moreover, while the values of the exponents for K* 1.4 and K* 0.6 are defini- lengths, (10) and (11) , evolve with the same power law. tively smaller than 1/2, for K* 1 it is perfectly consistent This happens to be the case for K* 1 and K* 1.4, with the standard Allen-Cahn value. Notice the long plateau whereas for K* 0.6 both sets of profiles scale indepen- obtained in the case K* 1.4. One needs to perform very dently according to widths evolving with different power long simulations before reaching the algebraic growth re- laws. gime. In fact, we have also obtained such behavior for some values of K* 1 inside the region 0.5 K* 0. In addi- tion, the extension of the plateau depends on K*, suggesting that it is related to the existence of an activated process with energy barriers depending on K*. These energies increase but not symmetrically as one varies K* from K* 1, either to K* 1 or to K* 1, and in both cases hinder the motion of the vacancies at the vertices of the interfaces. The expres- sions for the associated energy barriers (Eb* Eb /J) are Eb* 1 K* for K* 1 and Eb* 3(K* 1) for K* 1. In Fig. 5 black dots are the times needed to reach the algebraic regime for the different values of K*. These have been esti- mated from the simulations as the ending points of the pla- teau. Simultaneously we have plotted straight lines the bar- rier passing time, defined as exp(Eb*/kBT). As an example, to illustrate how black dots in Fig. 5 have been obtained, we show the case of K* 0 Fig. 6 . The evolution up to 107 MCs has been obtained by following first the standard dy- namical simulations circles and next, by using the N-fold way algorithm squares . The arrow indicates the estimated time for the ending point of the plateau. We have tested the existence of dynamical scaling. Fig- ures 7 a , 7 b , and 7 c show the scaled structure factor profiles. Those along the (10) direction have been shifted downwards four decades in order to clarify the picture. Pro- files along each direction have been conveniently scaled with the corresponding . The overlap of the data is satisfactory FIG. 7. log-log plot of the scaled structure factor profiles in the except for the tails at large q's. Nevertheless, to prove the (10) and (11) directions for K* 0.6, K* 1.0, and K* 1.4. The existence of dynamical scaling it is necessary to have not profiles in the 10 direction have been shifted four decades below only the collapse of the curves but also one requires that both in order to clarify the picture. 56 EFFECT OF THE VACANCY INTERACTION ON . . . 5267 FIG. 9. Growth exponent vs K* obtained from the structure factor width . Open filled circles correspond to the evolution of (10) ( (11)). All the simulations are performed in a 200 200 square lattice with cV 0.06 at T 0.1J/kB . FIG. 8. Ratio (10) / (11) vs time for different values of the parameter K*. All these simulations are performed in a 200 200 of the intercoupling between bulk diffusing vacancies and square lattice with cV 0.06 at T 0.1J/kB . the interface motion. For K* 1 the intercoupling proceeds via an effective interaction originated from the vacancy- In the following two figures we present a complete study vacancy specific interaction. Moreover, this introduces dif- of both the anisotropic character of the growth and the ki- ferences in the internal structure of the interface. In the par- netic growth exponent s for a wide range of values of the ticular case of K* 1, this intercoupling reduces to a simple interaction parameter K*. Figure 8 shows the time evolution encounter between curvature-driven interfaces and mobile of the ratio (10) / (11) for different values of K*. For vacancies that mutually cross their respective trajectories. K* 0.8, definitively increases with time, showing that the This does not make the curvature ineffective but may slow shape of the ordered domains becomes more and more spike- down the domain growth. It has been shown that8 the effect like. For 0.8 K* 1.0 the ratio remains constant around of this simple intercoupling does not modify the essential 1, indicating that the domains are circular during all the time dependence of the growth law but modifies the growth evolution. For K* 1 after an initial decrease, reaches the rate prefactor that decreases as the mobile impurity concen- value 1/ 2, indicating that, at long times, the domains are tration increases. squarelike and grow isotropically. Figure 9 shows the growth We next discuss separately the other two cases. We start exponents obtained by fitting an algebraic growth law with the case of vacancy attraction (K* 0.6) and point out (t) tn to the evolution of both (10) open circles and some other relevant features present in Fig. 3. Notice the (11) filled circles . Note that the Allen-Cahn value increasing concentration of vacancies at the interfaces as (n 0.5) is only reached for values of K* 1. they evolve with time. This is premonitory of the phase sepa- We have also studied the behavior of the structure factor ration process, eventually reachable at longer times. During at large q's. It is affected by two different phenomena. i the regime shown in Fig. 3 one is mainly concerned with a Along the direction (1,1) for q 0.5 the structure factor is nonhomogeneous distribution of vacancies along the inter- distorted due to the existence of the nonscaling fundamental faces. The local accumulation at the vertices deserves special peak at k (00) strictly speaking the value of the structure mention. This is the signature of a previous fast process of factor at k (00) is always zero, but the peak exhibits some interface reduction due to the high curvature of the vertices . finite width due to the existence of disorder in the system . The further evolution is hindered until the vacancies diffuse ii For the cases in which the vacancies dissolve into the along the interfaces. This involves activated processes. In- bulk, there is a homogeneous nonscaling background that, in deed, in our simulations we have observed how the temporal turn, may evolve in time. pinning of the high-curvature portion of the interfaces pro- vokes that the further evolution of the interconnecting inter- V. DISCUSSION faces with lower curvature does not fulfill the main as- sumptions underlying the Allen-Cahn theory. The The differences in the values obtained for the kinetic ex- importance of this temporal pinning depends on K*, since ponent of the growth law lie on the different characteristics the energy barriers (Eb* 1 K*) hindering the motion of a 5268 PORTA, FRONTERA, VIVES, AND CASTAŽN 56 FIG. 10. Snapshots of some evolving domains directly ex- tracted from our simulations for K* 0. Vacancies (A and B par- ticles are shown as black white . The simulations are performed in a 200 200 square lattice with cV 0.06 at T 0.1J/kB . vacancy at the corner of the interface increase as K* de- in the limit of low vacancy concentration . Other studies13 creases from K* 1. In Fig. 10 we show this effect for the on the diluted ferromagnetic Ising model, encountered that, case K* 0. In particular, one observes how the single rect- by coupling both dynamics differently the simultaneous angular domains evolve so that they become more and more vacancy-spin exchange and spin flip is not allowed so that platelike. This is because the longer interfaces with a lower the vacancy at the corner is effectively pinned, the growth concentration of vacancies are the only ones able to evolve. stops. In the case of the alloy, the dynamics implemented Moreover, they remain almost flat and parallel to the 10 follows directly from the requirement of the conservation direction. This feature is not observed in previous studies of law for the number of particles and we found that the order- the diluted antiferromagnetic Ising model14 notice that it ing process is definitively described by a growth law that, corresponds to the BEG with K* 0) so that it cannot be although slower than the Allen-Cahn law, is definitively al- attributed exclusively to the interaction which favors the gebraic. For some values of K* this algebraic regime is pre- vacancies to be nearest neighbors at the interface . In Fig. 11 ceded by a plateau whose extension in time depends on K*. we schematically illustrate this mechanism of evolution. This More precisely, it is longer the bigger the energy barriers are. represents a typical domain in the regime under discussion In particular, we found that for K* 0 this algebraic regime here, that is, the regime characterized because the width of is visible only after 105 MCs see Fig. 6 . We notice that the interface is small and the only relevant length is the size the BEG model, in the particular case of K L 0, corre- of the AB ordered domains. Two facts have to be taken into sponds to a diluted Ising model. In view of this, we believe account: the attractive vacancy-vacancy interaction defined that the growth law for the diluted antiferromagnet is alge- by the Hamiltonian and the conserved character of the imple- braic. In fact, the authors of 14 did not exclude this possi- mented dynamics Kawasaki . The combination of both in- bility in their discussion. Concerning the complete pinning of troduces energy barriers hindering the motion of the vacan- the process reported by 13 , it is certainly due to the ex- cies circles at the interfaces, which hinders the shrinkage of tremely low value of the temperature. the domain. In Fig. 11 we have indicated in black white the Later on, as the width of the interface increases, this pin- vacancies with associated energy barriers so that they are ning effect, localized at the corner, becomes less important induced to go outward inward . The vacancies in gray do and one expects the system to cross over to a completely not have any preference. The crucial point is that the inwards different regime. In this asymptotic regime, not reached in motion of the vacancies at the corner is strongly hindered our simulations, discussing the interface, as formed by the whereas for its neighbors, it is favored. This provokes that increasing accumulation of vacancies, is meaningless. the shrinkage of the domains proceeds by displacing the flat Rather, one would deal with a phase-separation process, interfaces and accumulating the excess vacancies at the which is not studied in the present work. vertices.26 The domains then become spikelike along the 11 We now focus on the discussion for K* 1.4. In this case, direction, breaking down the single-length dynamical scaling the repulsive interaction between vacancies favors them to be so that the growth becomes anisotropic. In Ref. 14 , the next nearest neighbors at the interfaces. As previously, the authors coupled the Glauber dynamics for the spins to the driving force is contained in the corner but its motion is conserved dynamics for the number of vacancies in such a hindered by a barrier of energy 2(K* 1). The interface way that the vacancy at the corner is not pinned. Neverthe- connecting two vertices is directed along the 11 direction less they found that the dynamical evolution of the ordering and evolves so that it displaces in a parallel manner. The process is effectively described by a logarithmic growth law energy barrier associated to this mechanism is FIG. 11. Schematic represen- tation of the evolution of a square domain typically observed in the K* 1 simulations. The circles represent vacancies whereas the white black squares represent A (B) particles. The mechanism of evolution is emphasized by show- ing in black white those vacan- cies that have the tendency to move outward inwards . The va- cancies painted in grey are the ones that are, in this sense, indif- ferent. 56 EFFECT OF THE VACANCY INTERACTION ON . . . 5269 order-disorder transition. The study is performed by Monte Carlo simulations in the limit of low vacancy concentration for a wide range of values of the biquadratic coupling pa- rameter K*, which controls the specific vacancy-vacancy in- teraction. For all values of K* inside the range 0.5 K* 1.4 we found that the vacancies tend to concen- trate at the interfaces. This feature introduces, via the param- eter K*, an intercoupling between diffusing bulk vacancies and moving interfaces. In the particular case of K* 1 this intercoupling does not take place via any specific interaction and the ordering process is consistent with the Allen-Cahn law. In fact we find that n 1/2 for K* 1 whereas the ex- ponent is clearly smaller when such interaction is present, no matter if it is attractive (K* 1) or repulsive (K* 1). Nev- ertheless, our results clearly show that the growth is defini- tively algebraic. When K* 1 the attractive vacancy interac- tion favors the increasing accumulation of vacancies at the interfaces as the system evolves. The regime of interest here FIG. 12. Domain area vs time of a single domain for K* 1.2 corresponds to the very initial stages of a phase separation solid and dotted lines and K* 1.4 dashed line . The solid dot- process when the width of the interfaces is small and the ted line corresponds to cV 0.06 (cV 0.04 . We simultaneously only relevant length is the size of the AB ordered domains. show, in the inset, the time evolution of the total amount of bulk Our main finding is that for quenches inside the coexistence vacancies in the system. The simulations are performed in a 200 region the growth for the binary alloy is, in this regime, 200 square lattice at T 0.1J/kB . anisotropic. This is related to the existence of energy barriers depending on K*) that hinder the motion of the vacancies at E the vertices of the 10 squarelike domains and simulta- b* 3(K* 1). The vacancies, in excess as a consequence of the interface reduction, are in this case vacancy repul- neously provoke local accumulations of vacancies along the sion expelled to the bulk. Moreover, the migrating interface 11 directions. Furthermore, these barriers may delay the has to cope with the effective repulsive interaction due to the apparition of the algebraic regime of the growth process to bulk vacancies, whose concentration increases as the system very late times. Concerning the underlying mechanism for evolves. Indeed, one expects this should interfere with the the motion of the interfaces, it is not purely curvature driven dynamics. At this point, it is interesting to discriminate and the two effective exponents needed to describe the pro- whether or not this interference makes the curvature-driven cess are lower than the Allen-Cahn value. Nevertheless, both mechanism become ineffective. In this sense, we have veri- tend to n 1/2 as K* 1. For K* 1 the specific vacancy- fied that the interface evolves covering a domain area con- vacancy interaction is repulsive and the interfaces of the stant in time. This is shown if Fig. 12 for single domains square-like domains are in this case directed along the 11 directly extracted from our simulations. These calculations directions. As in the previous case, the motion of the vacan- have been performed using the optimized multigrid algo- cies at the vertices is an activated process. The associated rithm discussed previously. A linear time dependence for the energy barriers depend on K* and the algebraic growth re- domain area evolution is obtained. This is commonly ac- gime shows up only after the time needed for surpassing the cepted as indicative that the motion of the interface is curva- barrier. Since the interface motion has to cope with the re- ture driven.27 It then follows that the effect of the intercou- pulsive interaction with diffusing vacancies, we found that pling between mobile bulk vacancies and the evolving the ordering process clearly slows down. Nevertheless, this interfaces is, in this case, to slow down the global process as repulsion does not make the curvature ineffective. The effec- it is revealed by a decreasing in the effective growth expo- tive exponent is lower than the Allen-Cahn value but ap- nent but it does not make ineffective the curvature driven proaches n 1/2 as K* 1, as expected. mechanism, which remains as the underlying mechanism for the motion of the interfaces. Moreover the effective exponent ACKNOWLEDGMENTS tends to the ideal 1/2 Allen-Cahn value as K* 1. We acknowledge financial support from the ComisioŽn In- VI. SUMMARY terministerial de Ciencia y TecnologiŽa CICyT, project num- ber MAT95-504 and supercomputing support from Fun- We use a Blume-Emery-Griffiths model with L 0) to dacioŽ Catalana per a la Recerca F.C.R. and Centre de study the influence of mobile vacancies on the kinetics of SupercomputacioŽ de Catalunya CESCA . M.P. and C.F. domain growth in a stoichiometric binary alloy after also acknowledge financial support from the Comissionat per quenches to a very low temperature (T 0.1J/kB) through an a Universitats i Recerca Generalitat de Catalunya . 5270 PORTA, FRONTERA, VIVES, AND CASTAŽN 56 1 Dynamics of Ordering Processes in Condensed Matter, edited by 16 K. Oki, H. Sagane, and T. Eguchi, J. Phys. France Colloq. 38, S. Komura and H. Furukawa Plenum, New York, 1988 . C7-414 1977 . 2 J.D. Gunton, M. San Miguel, and P.S. Sahni, in Phase Transitions 17 T. Ohta, K. Kawasaki, and A. Sato, Phys. Lett. A 126, 93 1987 . and Critical Phenomena, edited by C. Domb and J.L. Lebowitz 18 H. Gilho"j, C. Jeppesen, and O.G. Mouritsen, Phys. Rev. E 53, Academic, New York, 1983 , Vol. 8. 5491 1996 . 3 S.M. Allen and J.W. Cahn, Acta Metall. 27, 1085 1979 . 19 H. Gilho"j, C. Jeppesen, and O.G. Mouritsen, Phys. Rev. Lett. 75, 4 I.M. Lifshitz and V.V. Slyozov, J. Phys. Chem. Solids 19, 35 3305 1995 . 1961 . 20 S.E. Nagler, R.F. Shannon, C.R. Harkless, M.A. Singh, and R.M. 5 M.K. Phani, J.L. Lebowitz, M.H. Kalos, and O. Penrose, Phys. Nicklow, Phys. Rev. Lett. 61, 718 1988 . Rev. Lett. 45 366 1980. 21 M. Blume, V.J. Emery, and R.B. Griffiths, Phys. Rev. A 4, 1071 6 P.S. Sahni, G. Dee, J.D. Gunton, M. Phani, J.L. Lebowitz, and M. 1971 . Kalos, Phys. Rev. B 24, 410 1981 . 22 E. Vives and A. Planes, Phys. Rev. Lett. 68, 812 1992 ; Phys. 7 H.C. Fogedby and O.G. Mouritsen, Phys. Rev. B 37, 5962 Rev. B 47, 2557 1993 ; C. Frontera, E. Vives, and A. Planes, 1988 , and references therein. ibid. 48, 9321 1993 . 8 M. Porta and T. CastaŽn, Phys. Rev. B 54, 166 1996 . 23 Details of the mean-field calculation, using the self-consistent 9 F. Bley and M. Fayard, Acta Metall. 24, 575 1976 . field method, will be given elsewhere. 10 R.F. Shannon, C.R. Harkless, and S.E. Nagler, Phys. Rev. B 38, 24 G. Porod, in Small Angle X-Ray Scattering, edited by O. Glatter 9327 1988 . and L. Kratky Academic, New York, 1982 . 11 D.A. Huse and C.L. Henley, Phys. Rev. Lett. 54, 2708 1985 . 25 A.B. Bortz, M.H. Kalos, and J.L. Lebowitz, J. Comput. Phys. 17, 12 H. Gilho"j, C. Jeppesen, and O.G. Mouritsen, Phys. Rev. E 52, 10 1975 . 1465 1995 . 26 If this analysis is done in a circular domain, it can be seen that it 13 D.J. Srolovitz and G.N. Hassold, Phys. Rev. B 35, 6902 1987 . first evolves to a square domain and then proceeds as discussed. 14 O.G. Mouritsen and P.J. Shah, Phys. Rev. B 40, 11 445 1989 ; 27 E. Domany and D. Kandel, in Cellular Automata and Modeling of P.J. Shah and O.G. Mouritsen, ibid. 41, 7003 1990 . Complex Physical Systems, edited by P. Manneville, N. Boccara, 15 O.G. Mouritsen, P.J. Shah, and J.V. Andersen, Phys. Rev. B 42, G.Y. Vichniac, and R. Bidaux Springer-Verlag, Berlin, 1990 , 4506 1990 . p. 98.