- Email: [email protected]
Cellular automaton for the order-disorder transition
Pergamon 0960-0779(94)00289-4
Chaos, Solitons & Fractals Vol. 6, pp. 439M-45, 1995 Copyright © 1995 Elsevier Science Ltd Printed in Great Britain. All rights reserved 0960--0779/95 $9.50 + .00
C e l l u l a r A u t o m a t o n for t h e O r d e r - D i s o r d e r T r a n s i t i o n
J.PALANDI
la, R . M . C . d e
A L M E I D A 1, J . R . I G L E S I A S 1 a n d M . K I W I 3
1 Instituto de F/sica, Universidade Federal do Rio Grande do Sul C.P.15051, 91501-970 - Porto Alegre- RS - Brazil 2 Departamento de F/sica, Universidade Federal de Santa Maria Campus de Camobi, 97119-900 - Santa Maria- RS - Brazil 3 Facultad de Fisica, Pontificia Universidad Catdlica Casilla 306, C,o r r e o 22 - Santiago - Chile
A b s t r a c t - We propose a probabilistic cellular automaton model with a rule that interchanges neighboring atoms to simulate the diffnsion driven dynamics of binary alloys. In particular, we consider the ordering transition of two dimensional binary alloys on a square lattice by assuming antiferromagnetic couplings and measure the energy, specific heat, and correlation functions as functions of temperature. We perform damage spreading simulations and the phase space structure that emerges consists of a single paramagnetic thermodynamic state at high temperatures and a two-basin phase space at low temperatures. The results are compatible with an order-disorder phase transition at a temperature close to the Onsager exact value.
INTRODUCTION Alloys like A1Zn, for example, undergo an order-disorder phase transition from a disordered state where the two atomic species are randomly distributed over the lattice sites to an ordered state where the lattice separates into two pure phases. On tile other hand, alloys like CuZn(fl-brass) and CuAu, for example, undergo an order-disorder phase transition to an ordered mixed state where each atomic species preferentially occupy sites on two cubic sublattices. For these substitutional binary alloys, the model Hamiltonian for nearest-neighbor interactions in terms of local occupation variables can be reduced to the Ising model Hamiltonian when the only degree of t'reedom is the thermally activated diffusion process[l]. The equilibrium properties of these systeuls may then be described by ferro- and antiferromagneticlike couplings, respectively, and the diffusion process can be modeled by interchanging two different neighboring atoms. The main approaches to investigate complex dynamical systems through computer simulation are molecular dynamics[2], Monte Carlo [3,4] and cellular automata[5] techniques. In particular, several cellular automata models have been applied to simulate Ising models[6]. For binary alloys, realistic simulations of the order-disorder problem should take into account not only the conservation of the number of atoms of each species but also the possibility of diffusion. All early attempt to match these conditions was performed [7,8], where the number of atoms of each species is conserved and time evolution is driven by thermally induced diffusion. Here we develop further this Ising 2D probabilistic cellular automaton model with nearest neighbor interchange dynamics to simulate the ordering of an anti-ferromagnetic 2D binary alloy. 439
440
J. PALANDI et al.
CELLULAR
AUTOMATON
MODEL
We consider N x N square lattices with periodic boundary conditions. The sites are occupied by either an A or B atom represented by Ising spins ,g'i,j = -1-1, where i and j are the coordinates of the lattice sites. Taking antiferromagnetic, nearest neighbors interactions, the energy of the system reads E = -J
~
,9'i,.i,~'k,,
(J < 0)
(1)
where < ij, kl > stands for the sum over neighboring sites (i,j) and (k, l). The evolution of the system is ruled by a dynamics whose basic step is the interchange of an A and a B atom on nearest neighboring sites with a probability given by
p =
exp(-/~zXE)
(2)
(1 + e x p ( - / J A E ) ) where fl is the inverse of temperature and A E is the energy difference due to the exchange of the two neighboring atoms. This rule differs from standard Kawasaki dynamics[9] because it considers interchange of only neighboring atoms in order to simulate a diffusion driven evolution. In each simulation step (ss) (i) the lattice is divided into two sublattices of horizontal alternate pairs; (ii) all pairs on each sublattice are simultaneously updated according to the rule based on Eq.2; (iii) steps (i) and (ii) are repeated with one lattice parameter horizontal translation; and (iv) the three steps are repeated in the vertical direction[10]. Fig.1 displays snapshots of the lattice evolution at different temperatures, suggesting the following scenario: at very high temperature there is no correlation between the sites, and as the temperature is lowered, short range correlations appear (Fig. l(a)) until a critical transition temperature is reached, when long range order is established. In the low temperature phase the system relaxes to a one-grain structure (Fig.l(b)), where local defects may surviv(, for very long times, with the positions of the defects continuously changing in time.
initial
•
.
500 ss
~
2 5 0 0 ss
~
(b)
initial
100 ss
200 ss
Fig.1. Snapshots of the evolution of the system from a random initial state at (a) T = 2.5 and (b) T = 1.5.
Probabilistic cellular automaton model
441
RESULTS The rth--~neighbor equal-time correlation function at time step t after relaxation is given by N
1
N
=
1
. =
(k,I)=~
where N~ is the n u m b e r of r ~h neighbors of a given site and (k, l) = < i,j > , stands for the summation over the (k,l) sites that are r Lh neighbors to site (i,j). We adopt sample sizes of 32 x 32, 64 x 64, and 96 x 96 sites, and measure the correlation functions up to the 22-shell for the 32 x 32 lattice, and up to the 43-shell for the others. The relaxation times are chosen large enough to ensure that the energy is stabilized around its mean value. At high temperatnres, the correlation function (Fig.2) shows a short range order that increases as the temperature approaches the phase transition from above. At T = 1.5, the correlation function is practically independent of distance, indicating a long range ordered phase. 1.0
,~ " £ k ' k "
k 'k'a'h'/='X'£'
m
[]
o
o
o
o
o
o
o
o z~ o
0.5
O o
O o
<>
/~" £ ' ,~ ' []
o
o
D
T=1.5 2.0 2.5 3.0
<>
o v'-
c)
0.0 o
o
<>
<>
o
-0.5
'O 0
-I .0
[]
I"1
[]
D
[]
O
[]
D
D
n
[]
[]
[]
D
4.4..A..~..~,4.4.4.4..&,4.4.4.4.4 10
30
20
r
Fig.2. Time average (over 200 ss after a relaxation of 400 ss) of the correlation functions up to tile 30 ~ shell for the 64 × 64 lattice at T = 1.5, T = 2.0, T = 2.5, and T = 3.0. Tile error bars are smaller than the size of tile symbols used.
1.0
o []
"<>'~
r=8
' o n
96x96 64X64
0.5
o
0.0
-0.5
-1.0
,i 1.5
2.0
,
2.5
3.0
T
Fig.3. Time average (over 50 ss after a relaxation of 100 ss) of the correlation functions as functions of temperature for the 1~, 8~ , and 15~ shells for 32 x 32, 64 x 64, and 96 x 96 lattices.
J. PALANDIet al.
442
To calculate the correlation flmction as function of temperature for different shells we simulate an annealing process by fixing a high initial temperature, relaxing the system for I00 s,~, and taking time averages of the correlation function (Eq. 3) over 50 ss. The temperature is then lowered and the same procedure is repeated. In Fig.a the vertical dotted line indicates the transition temperature T~ = 2.27 obtained by Onsager for the square lattice Ising model[11]: the data are eo, sistent with the correlation function going to zero at T = T~ in the infinite sample size limit (r --+ oo). To calculate the energy of the system as a function of the temperature we simulate an annealing process analogous to the above mentioned, but now letting the system relax for 10 ss for the 32 × 32 and 50 ss for the 64 × 64 samples, and taking time averages of the energy per bond over 50 ss and 100 ss, respectively. Ensemble averages over 20 samples are also taken. The specific heat per bond is then calculated as the mean square energy fluctuation. Fig,4 shows simulation data that agree fairly well with theoretical results and are consistent with a phase transition at T = 2,27. 2.0 []
1.6
o
64x64 32x32 Onsager
1.2 ¢0
0.4
0.0
1.0
"
3.0
2.0
4.0
T
Fig.4. Specific heat pet" bond averaged over time and samples as function of temperature for 32 × 32 and 64 × 64 lattices. The dotted line stands for the exact Onsager result and the error bars correspond to the ensemble average.
Finally, to investigate how the stochasticity of the initial conditions determines the final state of the system, independent of the inherent stochasticity of the thermal noise, we performed damage spreading simulations. We consider two initial configurations. ,~0) and ,~(2), which differ in a number D(0) of misoriented spins. The damage spreading idea is to let them evolve subject to the same thermal noise in such a way that they should be brought together when, below T~, both are in the basin of attraction of the same thermodynamic state. Hence, the asymptotic number of misoriented spins, D ( t --+ oo), may signal the transition temperature (see, for example, [12,13,14]). We measured D ( t ---+ oo) for two different protocols: /) Glauber-like protocol (as it is called in [14]) where a same random number 1 > z] > 0 is compared with the interchange probability P (Eq. 2) of each configuration and ii) Direction protocol[15J where an interchange direction is first randomly chosen, associated with a random number z2 and then whenever an interchange in this direction is possible, the associated probability P is compared with Zl. Again, the same random nmnber sequences, for zl and z2 respectively, are used in both samples evolutions. As the two configurations evolve, following either one of the two protocols, we monitor the damage given by
D(t)
=
~
=,
,,--,
,
,
(4)
Probabilistic cellular automaton model
443
up to enough long times to assure that the asymptotic value D(t ~ oo) has been attained. We have adopted samples of size 32 x 32 and 64 x 64 and for" each value of the t e m p e r a t u r e and each protocol, a sample is relaxed (about 200 ss), a copy with given initial damage (D(0) = 0.1N 2 and D(0) = 0.9N 2) is made, the two lattices are allowed to evolve until the damage is stabilized (over 300 ss) and a time average is performed (over 100 ss). For low temperatures and initial damages, the asymptotic value D(t --+ oc) is low, indicating that both configurations remains in the same basin of attraction. As the t e m p e r a t u r e increases, still below the transition, D(t --~ oc) increases due to thermal fluctuations around the energy minimum. After transition, either D(t ~ ec) = 0.5N 2 for Clauber-like protocol, implying that the configurations are unrelated, or D(t --+ cx~) decreases slowly to zero for direction protocol. In the case where the initial damage is high the situation is analogous, with the configurations in opposite basins at low temperatures. Finally, above transition temperature, the final damage is independent of the initial conditions. The results for a single pair of samples are plotted in Fig.5 and for an ensemble of 20 pairs of samples in Fig.6. Both figures show that there are only two basins of attraction in phase space at low temperatures, and only one high t e m p e r a t u r e stable thermodynanfic state.
1.0
II~l
0.2 I ~-_~o i i G l a u b e r T / ~ , ~
t
_.
J
0.0 1.6
2.0
1.8
iJ
~
-I
2.2
o,o,.o,
2.4
2.6
2.8
T Fig.5. Time average damage as function of t e m p e r a t u r e for Glauber-like and direction protocols for :12 x 32 and 64 x 64 lattices.
1.0 0.8 0.6
!
0.4
i
3cl
0.2
~/;~,~, ~ D{O)=0.9 . ~ o ) r = D . ( g ) =0.1 "
0.0 1.8
~*K
I
I
2.0
a
I
2.2
I
i
I
2.4 T
t
I
2.6
a
I
*
2.8
Fig.6. Ensemble average da,Illage as ['unct, ioll of t e m p e r a t u r e for Glauber-like and direction protocols for the 32 × 32 lattice.
444
J. PALANDt et al.
A manifest difference between the two protocols in Fig.5 is observed in the vicinity of the phase transition, where the direction protocol yields a more precise deterinination of the critical temperature. By allowing forbidden interclianges of atoms, that is, by neglecting the direction of the stochastic electrochemical potential which may counteract the thermal noise, the Glauber-like protocol leads to an overestimation of the transition telnperature, as it was suggested in [15]. On the other hand, in Fig.6 finite size effects hinder the difference between the two protocols: the energy barriers between the two thermodynamic states are not infinite and the transition probability from one basin of attraction to the other is not zero, even below T = 2.27. CONCLUSIONS We propose a cellular automaton model to simu]ate the diffusion driven time evolution of binary alloys where an order-disorder phase transition oCcllrs. We considered a square lattice with first-neighbor antiferromagnetic couplings and measured the energy, specific heat and correlation functions as functions of the temperature. The model presented is tractable within reasonable computing times in a standard workstation for sample sizes up to 96 x 96. The obtained order-disorder phase transition temperature is consistent with the Onsager value for the Ising square lattice and an analysis of damage spreading confirms that the phase space consists of a single thermodynamic state equivalent to a paranmgnetic one at high temperature and a two basin phase space, each basin corresponding to one of the possible antiferromagnetic ordering punctuated by local defects randomly walking through the system at low temperature. The damage spreading data of direction protocol provides a more precise determination of the critical temperature than Glauber-like one.
Acknowledgement -This work has been partially supported by brazilian agencies CNPq (Conselho Nacional de Desenvolvimento Cient/fico e Tecnoldgico), and (',APES (Coordena~;£o de Aperfei~.oamento de Pessoal de N/vel Superior). Partial support has also been granted by the chilean Consejo Nacional de Investigaciones Cientificas y Tecnoldgicas (('ONI(!YT) and the Fundacidn Andes.
References [1] K. Huang, Statistical Mechanics, John Wiley, New York (1987). [2] A. Rahman, Correlations in the motion of atoms in liquid argon, Phys.Rev. A 136, 405 (1964). [3] N. Metropolis, A. W. Rosenbluth, M. N. Itosenbhlth, A. H. Teller, and E. Teller, Equation of state calculations by fast computing machines, .1.(~heln. Phys 21, 1087 (1953). [4] Applications of the Monte Carlo Method in Statistical Physics, edited by K. Binder, Springer, Berlin (1984); D. W. Heermann, An Introduction to Computer Simulation Methods in Theoretical Physics, Springer, Berlin (.1986). [5] S. Wolfram, Statistical mechanics of cellular automata. Rev.Mod.Phys. 55, 601 (1983); . . . . . . Theory and Applications of Cellular Automata, World Scientific, Singapore (1986).
,
[6] See for instance: Y. Pomeau, Invariant in cellular automata, &Phys. A: Math.Gen. 17, L415 (1984); M. Creutz, Deterministic lsing dynamics, Ann.Phys. 167, 62 (1986); T. Toffoli and N. Margolus, Cellular Automata Machines, MIT Press, ('mnbridge, Mass. (1987). [7] M. C. L. Adams and J. R. Iglesias, A cellular antolnata for the order-disorder transition, in Workshop in Computational Physics and Cellular Automata, edited by A. Pires, D. P. Landau, and H. J. Herrmann, World Scientific, Singapore (1990), p. 163.
Probabilistic cellularautomatonmodel
445
[8] J. Palandi, R. M. C. de Almeida, .J.R. Iglesias, and M. Kiwi, submitted for publication (1993). [9] K. Kawasaki, Kinetics of Ising models, in Pha.~e Tvansition.s and Critical Phenomena, edited by C. Domb and M. S. Green, Academic, New York (1972), Vol. 2, p.443. [10] P. M. C. de Oliveira, T. J. P. Penna, S. M. Moss de Oliveira, and R. M. Zorzenon, Cellular automata as microcanonical simulators, J. Phys. A: M ath.Gen. 24, 219 (l 991 ). [11] L.Onsager, Crystal statistics.l. A two-dimensional model with an order-disorder transition, Phys.Rev. 65, 117 (1944); and B. Kaufmann, Crystal statistics.lI. Partition function evaluated by spinor analysis, Phys. Rev. 76, 1232 (1949). [12] B. Derrida, and G. Weisbuch, Dynamical phase transitions in 3- dimensional spin glasses, Europhys. Lett. 4, 657 (1987). [13] H.E. Stanley, D. Stauffer, .J. Kertdsz,and H..J. Herrmann, Dynamics of spreading phenomena in two-dimensional Ising models, Phys. Rev. Lett. 59, 2326 (1987). [14] A.M. Mariz, H.J. Herrmann, and L. de Arcangelis, Coml)arative study of damage spreading in the Ising model using heat-bath, Glauber, and Metrot)olis dynaulics, .l. Stat. Phys. 59, 1043 (1990). [15] R.M.C. de Almeida, Equivalence between Glauber and heat-bath dynamics in damage spreading simulations, J.Phys. France 3, 951 (1993).
Chaos, Solitons & Fractals Vol. 6, pp. 439M-45, 1995 Copyright © 1995 Elsevier Science Ltd Printed in Great Britain. All rights reserved 0960--0779/95 $9.50 + .00
C e l l u l a r A u t o m a t o n for t h e O r d e r - D i s o r d e r T r a n s i t i o n
J.PALANDI
la, R . M . C . d e
A L M E I D A 1, J . R . I G L E S I A S 1 a n d M . K I W I 3
1 Instituto de F/sica, Universidade Federal do Rio Grande do Sul C.P.15051, 91501-970 - Porto Alegre- RS - Brazil 2 Departamento de F/sica, Universidade Federal de Santa Maria Campus de Camobi, 97119-900 - Santa Maria- RS - Brazil 3 Facultad de Fisica, Pontificia Universidad Catdlica Casilla 306, C,o r r e o 22 - Santiago - Chile
A b s t r a c t - We propose a probabilistic cellular automaton model with a rule that interchanges neighboring atoms to simulate the diffnsion driven dynamics of binary alloys. In particular, we consider the ordering transition of two dimensional binary alloys on a square lattice by assuming antiferromagnetic couplings and measure the energy, specific heat, and correlation functions as functions of temperature. We perform damage spreading simulations and the phase space structure that emerges consists of a single paramagnetic thermodynamic state at high temperatures and a two-basin phase space at low temperatures. The results are compatible with an order-disorder phase transition at a temperature close to the Onsager exact value.
INTRODUCTION Alloys like A1Zn, for example, undergo an order-disorder phase transition from a disordered state where the two atomic species are randomly distributed over the lattice sites to an ordered state where the lattice separates into two pure phases. On tile other hand, alloys like CuZn(fl-brass) and CuAu, for example, undergo an order-disorder phase transition to an ordered mixed state where each atomic species preferentially occupy sites on two cubic sublattices. For these substitutional binary alloys, the model Hamiltonian for nearest-neighbor interactions in terms of local occupation variables can be reduced to the Ising model Hamiltonian when the only degree of t'reedom is the thermally activated diffusion process[l]. The equilibrium properties of these systeuls may then be described by ferro- and antiferromagneticlike couplings, respectively, and the diffusion process can be modeled by interchanging two different neighboring atoms. The main approaches to investigate complex dynamical systems through computer simulation are molecular dynamics[2], Monte Carlo [3,4] and cellular automata[5] techniques. In particular, several cellular automata models have been applied to simulate Ising models[6]. For binary alloys, realistic simulations of the order-disorder problem should take into account not only the conservation of the number of atoms of each species but also the possibility of diffusion. All early attempt to match these conditions was performed [7,8], where the number of atoms of each species is conserved and time evolution is driven by thermally induced diffusion. Here we develop further this Ising 2D probabilistic cellular automaton model with nearest neighbor interchange dynamics to simulate the ordering of an anti-ferromagnetic 2D binary alloy. 439
440
J. PALANDI et al.
CELLULAR
AUTOMATON
MODEL
We consider N x N square lattices with periodic boundary conditions. The sites are occupied by either an A or B atom represented by Ising spins ,g'i,j = -1-1, where i and j are the coordinates of the lattice sites. Taking antiferromagnetic, nearest neighbors interactions, the energy of the system reads E = -J
~
,9'i,.i,~'k,,
(J < 0)
(1)
where < ij, kl > stands for the sum over neighboring sites (i,j) and (k, l). The evolution of the system is ruled by a dynamics whose basic step is the interchange of an A and a B atom on nearest neighboring sites with a probability given by
p =
exp(-/~zXE)
(2)
(1 + e x p ( - / J A E ) ) where fl is the inverse of temperature and A E is the energy difference due to the exchange of the two neighboring atoms. This rule differs from standard Kawasaki dynamics[9] because it considers interchange of only neighboring atoms in order to simulate a diffusion driven evolution. In each simulation step (ss) (i) the lattice is divided into two sublattices of horizontal alternate pairs; (ii) all pairs on each sublattice are simultaneously updated according to the rule based on Eq.2; (iii) steps (i) and (ii) are repeated with one lattice parameter horizontal translation; and (iv) the three steps are repeated in the vertical direction[10]. Fig.1 displays snapshots of the lattice evolution at different temperatures, suggesting the following scenario: at very high temperature there is no correlation between the sites, and as the temperature is lowered, short range correlations appear (Fig. l(a)) until a critical transition temperature is reached, when long range order is established. In the low temperature phase the system relaxes to a one-grain structure (Fig.l(b)), where local defects may surviv(, for very long times, with the positions of the defects continuously changing in time.
initial
•
.
500 ss
~
2 5 0 0 ss
~
(b)
initial
100 ss
200 ss
Fig.1. Snapshots of the evolution of the system from a random initial state at (a) T = 2.5 and (b) T = 1.5.
Probabilistic cellular automaton model
441
RESULTS The rth--~neighbor equal-time correlation function at time step t after relaxation is given by N
1
N
=
1
. =
(k,I)=
where N~ is the n u m b e r of r ~h neighbors of a given site and (k, l) = < i,j > , stands for the summation over the (k,l) sites that are r Lh neighbors to site (i,j). We adopt sample sizes of 32 x 32, 64 x 64, and 96 x 96 sites, and measure the correlation functions up to the 22-shell for the 32 x 32 lattice, and up to the 43-shell for the others. The relaxation times are chosen large enough to ensure that the energy is stabilized around its mean value. At high temperatnres, the correlation function (Fig.2) shows a short range order that increases as the temperature approaches the phase transition from above. At T = 1.5, the correlation function is practically independent of distance, indicating a long range ordered phase. 1.0
,~ " £ k ' k "
k 'k'a'h'/='X'£'
m
[]
o
o
o
o
o
o
o
o z~ o
0.5
O o
O o
<>
/~" £ ' ,~ ' []
o
o
D
T=1.5 2.0 2.5 3.0
<>
o v'-
c)
0.0 o
o
<>
<>
o
-0.5
'O 0
-I .0
[]
I"1
[]
D
[]
O
[]
D
D
n
[]
[]
[]
D
4.4..A..~..~,4.4.4.4..&,4.4.4.4.4 10
30
20
r
Fig.2. Time average (over 200 ss after a relaxation of 400 ss) of the correlation functions up to tile 30 ~ shell for the 64 × 64 lattice at T = 1.5, T = 2.0, T = 2.5, and T = 3.0. Tile error bars are smaller than the size of tile symbols used.
1.0
o []
"<>'~
r=8
' o n
96x96 64X64
0.5
o
0.0
-0.5
-1.0
,i 1.5
2.0
,
2.5
3.0
T
Fig.3. Time average (over 50 ss after a relaxation of 100 ss) of the correlation functions as functions of temperature for the 1~, 8~ , and 15~ shells for 32 x 32, 64 x 64, and 96 x 96 lattices.
J. PALANDIet al.
442
To calculate the correlation flmction as function of temperature for different shells we simulate an annealing process by fixing a high initial temperature, relaxing the system for I00 s,~, and taking time averages of the correlation function (Eq. 3) over 50 ss. The temperature is then lowered and the same procedure is repeated. In Fig.a the vertical dotted line indicates the transition temperature T~ = 2.27 obtained by Onsager for the square lattice Ising model[11]: the data are eo, sistent with the correlation function going to zero at T = T~ in the infinite sample size limit (r --+ oo). To calculate the energy of the system as a function of the temperature we simulate an annealing process analogous to the above mentioned, but now letting the system relax for 10 ss for the 32 × 32 and 50 ss for the 64 × 64 samples, and taking time averages of the energy per bond over 50 ss and 100 ss, respectively. Ensemble averages over 20 samples are also taken. The specific heat per bond is then calculated as the mean square energy fluctuation. Fig,4 shows simulation data that agree fairly well with theoretical results and are consistent with a phase transition at T = 2,27. 2.0 []
1.6
o
64x64 32x32 Onsager
1.2 ¢0
0.4
0.0
1.0
"
3.0
2.0
4.0
T
Fig.4. Specific heat pet" bond averaged over time and samples as function of temperature for 32 × 32 and 64 × 64 lattices. The dotted line stands for the exact Onsager result and the error bars correspond to the ensemble average.
Finally, to investigate how the stochasticity of the initial conditions determines the final state of the system, independent of the inherent stochasticity of the thermal noise, we performed damage spreading simulations. We consider two initial configurations. ,~0) and ,~(2), which differ in a number D(0) of misoriented spins. The damage spreading idea is to let them evolve subject to the same thermal noise in such a way that they should be brought together when, below T~, both are in the basin of attraction of the same thermodynamic state. Hence, the asymptotic number of misoriented spins, D ( t --+ oo), may signal the transition temperature (see, for example, [12,13,14]). We measured D ( t ---+ oo) for two different protocols: /) Glauber-like protocol (as it is called in [14]) where a same random number 1 > z] > 0 is compared with the interchange probability P (Eq. 2) of each configuration and ii) Direction protocol[15J where an interchange direction is first randomly chosen, associated with a random number z2 and then whenever an interchange in this direction is possible, the associated probability P is compared with Zl. Again, the same random nmnber sequences, for zl and z2 respectively, are used in both samples evolutions. As the two configurations evolve, following either one of the two protocols, we monitor the damage given by
D(t)
=
~
=,
,,--,
,
,
(4)
Probabilistic cellular automaton model
443
up to enough long times to assure that the asymptotic value D(t ~ oo) has been attained. We have adopted samples of size 32 x 32 and 64 x 64 and for" each value of the t e m p e r a t u r e and each protocol, a sample is relaxed (about 200 ss), a copy with given initial damage (D(0) = 0.1N 2 and D(0) = 0.9N 2) is made, the two lattices are allowed to evolve until the damage is stabilized (over 300 ss) and a time average is performed (over 100 ss). For low temperatures and initial damages, the asymptotic value D(t --+ oc) is low, indicating that both configurations remains in the same basin of attraction. As the t e m p e r a t u r e increases, still below the transition, D(t --~ oc) increases due to thermal fluctuations around the energy minimum. After transition, either D(t ~ ec) = 0.5N 2 for Clauber-like protocol, implying that the configurations are unrelated, or D(t --+ cx~) decreases slowly to zero for direction protocol. In the case where the initial damage is high the situation is analogous, with the configurations in opposite basins at low temperatures. Finally, above transition temperature, the final damage is independent of the initial conditions. The results for a single pair of samples are plotted in Fig.5 and for an ensemble of 20 pairs of samples in Fig.6. Both figures show that there are only two basins of attraction in phase space at low temperatures, and only one high t e m p e r a t u r e stable thermodynanfic state.
1.0
II~l
0.2 I ~-_~o i i G l a u b e r T / ~ , ~
t
_.
J
0.0 1.6
2.0
1.8
iJ
~
-I
2.2
o,o,.o,
2.4
2.6
2.8
T Fig.5. Time average damage as function of t e m p e r a t u r e for Glauber-like and direction protocols for :12 x 32 and 64 x 64 lattices.
1.0 0.8 0.6
!
0.4
i
3cl
0.2
~/;~,~, ~ D{O)=0.9 . ~ o ) r = D . ( g ) =0.1 "
0.0 1.8
~*K
I
I
2.0
a
I
2.2
I
i
I
2.4 T
t
I
2.6
a
I
*
2.8
Fig.6. Ensemble average da,Illage as ['unct, ioll of t e m p e r a t u r e for Glauber-like and direction protocols for the 32 × 32 lattice.
444
J. PALANDt et al.
A manifest difference between the two protocols in Fig.5 is observed in the vicinity of the phase transition, where the direction protocol yields a more precise deterinination of the critical temperature. By allowing forbidden interclianges of atoms, that is, by neglecting the direction of the stochastic electrochemical potential which may counteract the thermal noise, the Glauber-like protocol leads to an overestimation of the transition telnperature, as it was suggested in [15]. On the other hand, in Fig.6 finite size effects hinder the difference between the two protocols: the energy barriers between the two thermodynamic states are not infinite and the transition probability from one basin of attraction to the other is not zero, even below T = 2.27. CONCLUSIONS We propose a cellular automaton model to simu]ate the diffusion driven time evolution of binary alloys where an order-disorder phase transition oCcllrs. We considered a square lattice with first-neighbor antiferromagnetic couplings and measured the energy, specific heat and correlation functions as functions of the temperature. The model presented is tractable within reasonable computing times in a standard workstation for sample sizes up to 96 x 96. The obtained order-disorder phase transition temperature is consistent with the Onsager value for the Ising square lattice and an analysis of damage spreading confirms that the phase space consists of a single thermodynamic state equivalent to a paranmgnetic one at high temperature and a two basin phase space, each basin corresponding to one of the possible antiferromagnetic ordering punctuated by local defects randomly walking through the system at low temperature. The damage spreading data of direction protocol provides a more precise determination of the critical temperature than Glauber-like one.
Acknowledgement -This work has been partially supported by brazilian agencies CNPq (Conselho Nacional de Desenvolvimento Cient/fico e Tecnoldgico), and (',APES (Coordena~;£o de Aperfei~.oamento de Pessoal de N/vel Superior). Partial support has also been granted by the chilean Consejo Nacional de Investigaciones Cientificas y Tecnoldgicas (('ONI(!YT) and the Fundacidn Andes.
References [1] K. Huang, Statistical Mechanics, John Wiley, New York (1987). [2] A. Rahman, Correlations in the motion of atoms in liquid argon, Phys.Rev. A 136, 405 (1964). [3] N. Metropolis, A. W. Rosenbluth, M. N. Itosenbhlth, A. H. Teller, and E. Teller, Equation of state calculations by fast computing machines, .1.(~heln. Phys 21, 1087 (1953). [4] Applications of the Monte Carlo Method in Statistical Physics, edited by K. Binder, Springer, Berlin (1984); D. W. Heermann, An Introduction to Computer Simulation Methods in Theoretical Physics, Springer, Berlin (.1986). [5] S. Wolfram, Statistical mechanics of cellular automata. Rev.Mod.Phys. 55, 601 (1983); . . . . . . Theory and Applications of Cellular Automata, World Scientific, Singapore (1986).
,
[6] See for instance: Y. Pomeau, Invariant in cellular automata, &Phys. A: Math.Gen. 17, L415 (1984); M. Creutz, Deterministic lsing dynamics, Ann.Phys. 167, 62 (1986); T. Toffoli and N. Margolus, Cellular Automata Machines, MIT Press, ('mnbridge, Mass. (1987). [7] M. C. L. Adams and J. R. Iglesias, A cellular antolnata for the order-disorder transition, in Workshop in Computational Physics and Cellular Automata, edited by A. Pires, D. P. Landau, and H. J. Herrmann, World Scientific, Singapore (1990), p. 163.
Probabilistic cellularautomatonmodel
445
[8] J. Palandi, R. M. C. de Almeida, .J.R. Iglesias, and M. Kiwi, submitted for publication (1993). [9] K. Kawasaki, Kinetics of Ising models, in Pha.~e Tvansition.s and Critical Phenomena, edited by C. Domb and M. S. Green, Academic, New York (1972), Vol. 2, p.443. [10] P. M. C. de Oliveira, T. J. P. Penna, S. M. Moss de Oliveira, and R. M. Zorzenon, Cellular automata as microcanonical simulators, J. Phys. A: M ath.Gen. 24, 219 (l 991 ). [11] L.Onsager, Crystal statistics.l. A two-dimensional model with an order-disorder transition, Phys.Rev. 65, 117 (1944); and B. Kaufmann, Crystal statistics.lI. Partition function evaluated by spinor analysis, Phys. Rev. 76, 1232 (1949). [12] B. Derrida, and G. Weisbuch, Dynamical phase transitions in 3- dimensional spin glasses, Europhys. Lett. 4, 657 (1987). [13] H.E. Stanley, D. Stauffer, .J. Kertdsz,and H..J. Herrmann, Dynamics of spreading phenomena in two-dimensional Ising models, Phys. Rev. Lett. 59, 2326 (1987). [14] A.M. Mariz, H.J. Herrmann, and L. de Arcangelis, Coml)arative study of damage spreading in the Ising model using heat-bath, Glauber, and Metrot)olis dynaulics, .l. Stat. Phys. 59, 1043 (1990). [15] R.M.C. de Almeida, Equivalence between Glauber and heat-bath dynamics in damage spreading simulations, J.Phys. France 3, 951 (1993).