Langmuir 1999, 15, 5883-5892
5883
Diffusion of Interacting Adsorbates on a Square Lattice† A. A. Tarasenko Institute of Physics, National Academy of Sciences of Ukraine, Prospect Nauki 46, 252028, Kijv 28, Ukraine
L. Jastrabı´k Institute of Physics, Academy of Sciences of the Czech Republic, Na Slovance 2, 180 40 Praha 8, Czech Republic
C. Uebing* Max-Planck-Institut fu¨ r Eisenforschung, Postfach 140444, D-40074 Du¨ sseldorf, Germany, and Lehrstuhl fu¨ r Physikalische Chemie II, Universita¨ t Dortmund, D-44227 Dortmund, Germany Received August 24, 1998. In Final Form: March 22, 1999 A two-dimensional lattice gas model with square symmetry is investigated by using the real-space renormalization group (RSRG) approach. Using the RSRG transformation with a cluster of 20 spins, we have investigated the influence of pairwise nearest and next nearest neighbor pair interactions as well as three- and four-particle interactions on the chemical adatom diffusion coefficient. In addition we have analyzed adsorption isotherms and adatom density fluctuations. We have also studied how interactions in the activated state influence surface diffusion.
by Fick’s law
1. Introduction The migration of adsorbates on solid surfaces plays an essential role in many physical and chemical processes such as adsorption, desorption, melting, roughening, crystal and film growth, catalysis, and corrosion. Understanding surface diffusion is one of the keys to controlling these processes. For a quantitative description of diffusion, two important concepts describing the general diffusion problem in three dimensions have been developed, each of them provides a different view of diffusion and allows to define its own “coefficient of diffusion”. In the following we will briefly recall these two concepts for the specific case of surface diffusion. The conceptually most simplest surface diffusion coefficient is the tracer surface diffusion coefficient Dt, which addresses the random walk of individual tagged particles on a two-dimensional lattice, i.e.,
Dt ) lim tf∞
1 〈|r b(t) - b r (0)|2〉 4t
(1)
Here b r(t) denotes the total displacement of the tagged adparticle as a function of time t. The tracer surface diffusion coefficient is a single-particle diffusion coefficient. There are several experimental techniques available to directly measure Dt for adsorbed particles. Principal among them is field ion microscopy (FIM). In contrast, the chemical surface diffusion coefficient Dc is a manyparticle diffusion coefficient phenomenologically defined * To whom correspondence may be addressed at Max-PlanckInstitut fu¨r Eisenforschung Max-Planck-Strasse 1, D-40237 Du¨sseldorf, Germany. Phone: 49-211-6792-290. FAX: 49-211-6792268. Email:
[email protected]. † Presented at the Third International Symposium on Effects of Surface Heterogeneity in Adsorption and Catalysis on Solids, held in Poland, August 9-16, 1998.
J ) - Dc∇n
(2)
Here n is the density of particles on the surface. The diffusion flux of particles J is driven by the gradient ∇n. It is important to note that the chemical surface diffusion coefficient is of interest for describing mass-transfer processes along the surface. Therefore, in the remaining part of this paper we will focus on the chemical diffusion coefficient. In recent years the effects of lateral interactions on the chemical surface diffusion coefficient of adsorbed particles have been intensively investigated using many different theoretical methods: mean-field1 and Bethe-Peierls approximations,2 the real-space renormalization group treatment,3,4 transfer matrix techniques5 and Monte Carlo modeling.6 It was found that adparticle interactions can strongly influence surface diffusion, especially at low temperatures and in the close vicinity of surface phase transitions. From simple physical considerations, it is intuitively expected that attractive interactions between adsorbed species inhibit the adparticle migration and thus slow surface diffusion. In contrast, repulsive interactions are expected to accelerate surface diffusion. Despite their simplicity, these rules describe the qualitative behavior of surface diffusion processes for many systems, at least in the limit of high temperatures, quite well.6 (1) Tarasenko, A. A.; Chumak, A. A. Sov. Phys. Solid State 1982, 24, 1683. (2) Tarasenko, A. A.; Chumak, A. A. Sov. Phys. Solid State 1980, 22, 1716. (3) Tarasenko, A. A.; Chumak, A. A. Poverkhnost′ Fizika, Khimija, Mekhanika 1989, 11, 98 (in Russian). (4) Tarasenko, A. A.; Jastrabik, L.; Uebing, C. Phys. Rev. B 1998, 57, 10166. (5) Myshlyatev, A. V.; Stepanov, A. A.; Uebing, C.; Zhdanov, V. P. Phys. Rev. B 1995, 52, 5977. (6) Uebing, C.; Gomer, R. J. Chem. Phys. 1991, 95, 7626, 7636, 7641, 7648.
10.1021/la981097u CCC: $18.00 © 1999 American Chemical Society Published on Web 08/07/1999
5884 Langmuir, Vol. 15, No. 18, 1999
Tarasenko et al.
a two-dimensional square lattice with lattice constant a (as shown in Figure 1). Foreign particles, adsorbed on the surface, are assumed to exclusively occupy the minima of the potential, i.e., the adsorption sites. If the depth of these potential minima, , is much larger than the thermal energy, . kBT, the adparticles will stay within the minima and from time to time perform jumps to nearest neighbor empty sites. The duration of such a jump τj is much shorter than the mean waiting time between jumps. In this case we can define a set of site occupation numbers {ni} by
{
1, if site i is occupied 0, if site i is empty.
ni )
(3)
A given set {ni} specifies a configuration of the whole adparticle system which in thermodynamic equilibrium is described by the statistical operator F,
F ) Q-1 exp[β(µNa - Ha)] Figure 1. Square lattice gas showing the c(2 × 2) ordered lattice gas phase. Blocks of sites for 5 × 4 RSRG transformation investigated in this paper are shown.
In general, the determination of the chemical surface diffusion coefficient requires the solution of a kinetic equation for a system of many strongly interacting particles. However, if interactions are restricted to affect only the ground state and not the activated state of the diffusing adparticles, the problem can be reduced to the calculation of purely thermodynamic quantities.7-11 In this case, the problem of determining the adparticle diffusion coefficient is reduced to the calculation of the free energy of the system. The model was used widely in most papers concerned with the theory of surface diffusion and computer simulation of the process due to the simplicity of the expression for the diffusion coefficient. In this project we have investigated interacting lattice gases on a square lattice. Besides pairwise nearest neighbor interactions we have also considered the influence of pairwise next nearest neighbor interactions and of three- and four-particle interactions affecting the ground-state of adsorbed particles. In addition, interactions in the activated state at the saddle point of the surface potential relief are taken into account. We have determined thermodynamic properties, such as adsorption isotherms and adparticle density fluctuations, as well as the chemical diffusion coefficient. The mean square density fluctuations of adparticles are found to determine the critical behavior of the chemical diffusion coefficient. The outline of this paper is as follows: the derivation of the diffusion equation and expressions for the chemical diffusion coefficient are described in section 2. The realspace renormalized group (RSRG) approach employed in this work is described in section 3. The results of our calculations are finally presented in section 4. 2. Diffusion of Adparticles on a Square Lattice Let us consider an idealized solid surface of square symmetry. The potential relief minima of the surface form (7) Bowker, M.; King, D. A. Surf. Sci. 1978, 71, 583. (8) Reed, D. A.; Ehrlich, G. Surf. Sci. 1981, 105, 603. (9) Reed, D. A.; Ehrlich, G. Surf. Sci. 1981, 102, 588. (10) Chumak, A. A.; Tarasenko, A. A. Surf. Sci. 1980, 91, 694. (11) Tarasenko, A. A.; Tomchuk, P. M.; Chumak, A. A. Fluctuations in Bulk and on Surfaces of Solids; Naukova Dumka: Kijv, 1992, in Russian.
(4)
Here β ≡ 1/kBT; µ, Na ) ∑ini, and Ha denote the chemical potential, the number of adparticles, and the Hamiltonianof the system, respectively. The latter has the following form
Ha ) - Na + φN
∑ ninj + φ2N ∑ ninj + 〈N〉
φ3
∑ 〈3〉
〈NN〉
ninjnk + φ4
∑ ninjnknl
(5)
〈4〉
Here φN, φ2N, φ3, and φ4 denote pairwise nearest and next nearest neighbor interactions as well as three- and fourparticle interactions. The summations are carried out over all corresponding pairs of nearest (〈N〉) and next nearest neighbors (〈NN〉) and over all elementary triangles (〈3〉) and squares (〈4〉). Q is the grand partition function
Q)
∑ exp[β(µNa - Ha)]
(6)
{n}
Here the summation is carried out over all 2N configurations of the adparticle system. The occupation numbers are time dependent due to the jumps of adparticles. Here we restrict the considerations to nearest neighbor jumps only. An adparticle in the ith site can jump to one of its nearest neighbor sites if the destination is empty. The equation, describing the balance of adparticles for any given site i can be written as10 4
ni (t + ∆t) - ni(t) )
(Jji - Jij) ∑ j)1
(7)
where the sum goes over all nearest neighbors of site i (labeled 1, 2, 3, and 4 in Figure 1). Jij represents the number of jumps from site i to site j during the time interval ∆t. For small values of ∆t, the regular and fluctuating parts of the adparticle flux can be defined by
Jij ) νijnihj∆t + δJij
(8)
Here we have introduced a “hole” occupation number hj ≡ 1 - nj as abbreviation. νij is the frequency of jumps from the site i to site j. The multiplier nihj ensures that the
Diffusion of Interacting Adsorbates
Langmuir, Vol. 15, No. 18, 1999 5885
jump proceeds from occupied to empty sites only. δJij describes the fluctuations of the adparticle flux during the small time interval ∆t. The balance equation describes the relaxation kinetics of adparticle occupation numbers. The characteristic frequency of this process is equal to the mean value of the adparticle jump frequency 〈νijnihj〉 (the angular brackets 〈‚‚‚〉 denote averaging with the statistical operator F) and the characteristic length is equal to the lattice constant a. In general, the particle density fluctuations are large, i.e., the occupation numbers vary abruptly from 1 to 0 and vice versa and it is impossible to solve the nonlinear balance equation or to linearize it over small adparticle density fluctuations. However, in many cases the detailed kinetic description is not required when one investigates the relaxation of macroscopically large surface density inhomogeneities or long-wavelength, low-frequency density fluctuations. The characteristic length of such density disturbances is much greater than the lattice constant a, and its decay proceeds as a result of a large number of jumps of many adparticles. Therefore, the characteristic frequency of this process must be much less than the mean frequency of adparticle jumps. Let us consider the states of the adparticle system averaged over a large time interval τ (〈νijnihj〉τ . 1). The averaging smoothes out all processes with characteristic frequencies ω > τ-1. We assume that the averaging over τ is equivalent to averaging with a statistical operator F˜ which has the form of the statistical operator F but contains the chemical potential of adparticles µi gradually varying in space and time. The statistical operator F˜ is the localequilibrium operator of the system. This approximation works well if the system has two relaxation scales. The first (fast) scale of relaxation is determined by establishing local equilibrium in different small regions of the whole surface. The characteristic frequency and length of this stage are equal to 〈νijnihj〉 and a, respectively. During the time 〈νijnihj〉-1 local equilibrium is established over the distances a. In the second (slow) stage different regions of the surface come into equilibrium with each other via adparticle diffusion. The characteristic frequency ω and the characteristic length l must satisfy the conditions l . a and ω . τ-1. During the time interval τ many adparticles visit any given site i. Thus, the fluctuations of the occupation numbers and other physical quantities, averaged over τ, will be small as averages are taken over a large number of adparticles. One can obtain an explicit form of the local-equilibrium operator F˜ by developing the exponential term into a series of small deviations of the chemical potential δµi ) µi - µ,
F˜ ≈ F + δF ) F[1 + β
∑i δµi (ni - θ)]
(9)
Here θ ≡ 〈ni〉 is the mean surface coverage. The averaging of the balance equation with the localequilibrium statistical operator describes the evolution of the fluctuation the occupation number of the ith site
∂δni(t) ∂t
and N
δ(νijnjhj) ) β
∑ δµk 〈νijnihj(nk- θ)〉 k)1
(12)
JiL(t) is the Langevin source of density fluctuations. Its correlation function was derived in ref 10. To exclude fluctuations of the chemical potential δµk from eqs 10 and 12, an expression, relating fluctuations of adparticle fluxes δ(νijnjhj) with adparticle density fluctuations δni must be obtained. For this purpose details of the adparticle jump mechanism must be known. Here we restrict the considerations to nearest neighbor jumps only. An adparticle at site i can jump to one of its nearest neighbor sites if the destination is empty. This situation is schematically depicted in Figure 1 for the particle on site i jumping to site 1. In the activated state the jumping adparticle increases its distance to the sites 2 and 4 and approaches the sites 5 and 7. The diffusing adparticle must surmount the saddle point of the potential between the initial site and the final site. In the case of interacting lattice gases, the activation energies of jumps are modified in the presence of adjacent adparticles affecting (i) the ground state of the potential (i.e., the adsorption energies) and (ii) the activated state (i.e., the saddle point). However, the interactions between an activated particle at the saddle point of the potential barrier and its nearest neighbors are frequently neglected. This simple model can be justified if one supposes that the adsorbed particle considerably changes its charge during the jump. In a chemisorbed state, for instance, the particle can (partially) lose an electron from its outer shell (for electropositive adsorbates) or attract an additional electron from the substrate (for electronegative adsorbates), thus obtaining a positive or negative electrical charge. The charge could cause strong disturbances of the electron density around the adparticle and, in principle, could result in anisotropic nonadditive interactions between adparticles separated by some lattice constants a. At the saddle point the interaction between the activated particle and the substrate atoms decreases and the particle could restore its neutral state, thus weakening the lateral interaction with its neighbors. In this case it would be reasonable to neglect the interaction of an activated particle with its surrounding. Assuming that the jump frequency is influenced by the presence of nearest and next nearest neighbors affecting the ground state of the potential only, one can write the following expression
νi1 ) ν exp[βni(φN
∑ nj + φ2N ∑ nj + φ3 ∑ njnk + 〈N〉
〈2N〉
φ4
∑ 〈4〉
〈3〉
njnknl)] (13)
4
)
[δ(νjinjhi) - δ(νijnihj)] + JLi (t) ∑ j)1
(10)
where N
δni ) β
() ∂θ
δµj 〈ni(nj - θ)〉 δµi ∑ ∂µ j)1
(11)
As already mentioned the summations are carried out over all pairs of nearest and next nearest neighbors of site i and over all elementary triangles and squares having one vertex at site i. Similar expressions can be obtained for all the other nearest neighbor jumps, affecting the occupation of site i. With these expressions one can obtain a rather simple expression for the chemical diffusion coefficient which has first been given by Zhdanov (for a
5886 Langmuir, Vol. 15, No. 18, 1999
Tarasenko et al.
Using eqs 17 or 18, one can obtain the following expression for the chemical diffusion coefficient
detailed derivation see ref 12)
Dc ) D(0) exp(βµ) P00 θ-1
∂(βµ) ∂ ln θ
(14)
Here D(0) is the chemical diffusion coefficient in the limit of θ f 0 (which corresponds to the noninteracting lattice gas, i.e., the Langmuir gas), P00 is the probability of finding two holes in nearest neighbor sites, θ-1∂(βµ)/∂ ln θ is the thermodynamical factor. It is quite useful to introduce the free energy of the system F as
F ) kBT N-1 ln Q
(15)
The first and second derivatives of the free energy over µ and φN allow calculation of all quantities of eq 14 according to
θ)
∂F ∂µ
P00 ) 1 - 2θ + 〈nin1〉 ) 1 - 2θ -
1 ∂F 2 ∂φN
( )
∂µ ∂2F )θ ∂ ln θ ∂µ2
-1
(16)
Thus, the problem of calculating the chemical diffusion coefficient is reduced to the determination of the free energy F of the system of adparticles. In the following we would like to analyze the more general situation that lateral interactions also influence the saddle point of the potential for jumping adparticles. If interactions in the activated state are taken into account, a much more complex expression for the diffusion coefficient arises which includes many-particle correlation functions. We consider here a model introduced in ref 13 where only nearest neighbor interactions between the jumping adparticle and its surrounding are considered. Next nearest neighbor as well as three- and four-particle interactions in the ground state and in the activated state are neglected. In this case an activated adparticle at the saddle point interacts with its four nearest neighbors (labeled 2, 4, 5, and 7 in Figure 1). As the distances between activated adparticles and their neighbors are close to the lattice constant a (actually the relative distance increases by about 12% only), it is reasonable to assume that the interaction energy in the activated state is numerically equal to φN. Thus one obtains the following expression for the jump frequency for adparticle, which differs considerably from the previous case (eq 13)
νi1′ ) ν exp[βniφN(
∑ nj - ∑ nj)] 〈N〉
(17)
〈N*〉
Here the symbol 〈N*〉 represents the summation over all nearest neighbor sites of a given saddle point. By use of the actual labels for the nearest neighbor sites given in Figure 1, eq 17 simplifies to
νi1′ ) ν exp[βφ (n3 - n7 - n5)]
(17)
(12) Zhdanov, V. P. Elementary Physicochemical Processes on Solid Surfaces; Plenum: New York, 1991. (13) Uebing, C.; Gomer, R. Surf. Sci. 1997, 381, 33.
Dc′ ) D(0)θ-1
∂(βµ) exp[β(φN - µ - )] {〈nin1〉 + ∂ ln θ 2〈n3nin1〉 ∆ + 〈n3nin1n6〉 ∆2} (19)
where ∆ ≡ exp βφN-1. Equation 19 contains two manyparticle correlation functions 〈n3nin1〉 and 〈n3nin1n6〉 which cannot be calculated as derivatives of the free energy F. Therefore, approximate multiplicative expressions for the correlation functions must be used to calculate the chemical diffusion coefficient Dc′. 3. Real-Space Renormalization Group Treatment of the Square Lattice In the preceding section we have obtained expressions for the chemical diffusion coefficients Dc and Dc′ by using simple models for adparticle jumps. In the absence of interactions in the activated state the problem is reduced to the determination of the free energy of the system F. For this purpose approximate methods must be used. Even for the simplest model the problem remains too complex to be solved exactly. The well-known Onsager’s solution for a 2D Ising spin model was obtained in zero magnetic field, which is equivalent to half coverage, θ ) 0.5. In this section we briefly outline the RSRG method used for calculation of the free energy and its desired derivatives (see eqs 14 and 19). It is well-known that the lattice gas model with nearest neighbor interactions is equivalent to the Ising spin model in an external magnetic field. For the description of such a spin model a spin variable s ) (1 is needed. Spin variables and occupation numbers are simply related via
ni ) (1 + si)/2
(20)
Thus, empty lattice sites correspond to s ) -1; occupied ones are represented by s ) 1. It is easily possible to obtain the equivalent reduced Hamiltonian of the Ising model in the following form
HI ) Nc + h
∑i si + k ∑ sisj + d ∑ sisj + 〈N〉
〈2N〉
t
∑ sisjsk + f ∑ sisjsksl 〈3〉
(21)
〈4〉
where
c)
1 φ4 (µ + - φN - φ2N - φ3) 2 16
h)
φ4 1 (µ + - 3φ3) - φN - φ2N 2 4 k)d)-
φN φ 3 φ 4 4 2 8
φ4 1 (φ + φ3) 4 2N 16
t)-
φ3 φ4 8 16
f)-
φ4 16
(22)
Diffusion of Interacting Adsorbates
Langmuir, Vol. 15, No. 18, 1999 5887
The free energies of lattice gas and Ising spin systems are identical apart from a constant c,
F(µ,φN,φ2N,φ3,φ4) ) c + F(h,k,d,t,f)
(23)
Strong lateral interactions cause phase transitions in the system. The exact critical value of the nearest neighbor pair interaction parameter k* in the absence of an external magnetic field and other interactions was obtained by Kramers and Wannier14 as k* ) (0.5 ln(1 + 21/2) ≈ (0.4406868. However, to investigate the critical behavior of the system in an external magnetic field or with other pair and multiparticle interactions approximate methods are required. In the RSRG method of Niemeyer and van Leeuwen15 and Nauenberg and Nienhuis16,17 the whole lattice is divided into blocks with L sites as shown in Figure 1. The blocks must have the symmetry of the lattice and form the same lattice as the original one with the lattice constant L1/2a. A “block spin” SR is assigned to each block. For odd L the block spin SR is determined by the so-called “majority rule”18 L
sj) ∑ j)1
SR ) sign(
where sign(x) )
{
+1, if x > 0 (24) -1, if x < 0
Each value of the spin block SR ) (1 corresponds to 2L-1 site spin configurations. The RSRG transformation of the spin system reduces the number of independent variables and constitutes a transition from the set of N site spins {si} to N/L block spins {SR}. It can be performed as a partial summation in the grand partition function over all site spin configurations, which leaves all block spins unchanged. The RSRG transformation can be written as
original lattice to the lattice composed of the block spins will be exact. Unfortunately the problem cannot be simplified by changing the order of summation in eq 25. Therefore, some approximations must be used in order to carry out the partial summation in eqs 25 or 26. In the framework of the RSRG approach, one usually employs periodic boundary conditions. It is assumed that the whole lattice is given by the periodic continuation of a small cluster of blocks. In the following we will consider the 5 × 4 cluster shown in Figure 1. We suppose that the whole lattice is obtained as the periodic continuation of a single 5 × 4 cluster consisting of four different blocks of 20 different site spins. In this case one obtains the following relations between original and renormalized values of the interaction parameters
h1(S1 + S2 + S3 + S4) + 2k1 (S1 + S3) (S2 + S4) + 4d1 (S1S3 + S2S4) + 4t1 (S1S2S3 + S1S2S4 + S1S3S4 + S2S3S4) + 4f1S1S2S3S4 + g ) ΦS1S2S3S4 (27) The renormalized interaction parameters are labeled by the index 1. If definite values of the block spins S1,2,3,4 ) (1 are introduced into eq 27, one can easily obtain the system of equations for the renormalized interaction parameters
4h1 + 8k1 + 8d1 + 16t1 + 4f1 + g ) Φ1111 -4h1 + 8k1 + 8d1 - 16t1 + 4f1 + g ) Φ-1-1-1-1 2h1 - 8t1 - 4f1 + g ) Φ111-1 -2h1 + 8t1 - 4f1 + g ) Φ1-1-1-1 -8k1 + 8d1 + 4f1 + g ) Φ1-11-1
SR)Const
Q)
∑ exp[HI ({SR}) +g] ) ∑ ∑ {S}
{S}
{s}
-8d1 + 4f1 + g ) Φ11-1-1
exp [HI({si})], (25)
where the second summation on the right-hand side is carried out over all configurations {si} that do not change the values of the block spins {SR}. The key supposition of the RSRG method is that the result of the partial summation HI({SR}) would in essence have the same form as the original Hamiltonian eq 21, but probably with other interaction parameters and some additional terms that weakly influence the critical behavior of the system. We rewrite the RSRG transformation in the following form
Solving the system eq 28, one obtains the following relations between original and renormalized values of the interaction parameters
h1 )
∑ {s}
exp[HI({si})]} ≡ Φ({SR})
(14) Kramers, H. A.; Wannier, G. H. Phys. Rev. 1941, 60, 252. (15) Niemeyer, Th.; van Leeuwen, J. M. J. Physica 1974, 71, 17. (16) Nauenberg, M.; Nienhuis, B. Phys. Rev. Lett. 1974, 33, 1598. (17) Nienhuis, B.; Nauenberg, M. Phys. Rev. Lett. 1975, 35, 477. (18) The interested reader is referred to Th. Niemeijer and J. M. J. van Leeuwen (Niemeijer, Th.; van Leeuwen, J. M. J. Renormalization Theory for Ising-like Spin Systems in Phase Transitions and Critical Phenomena; Domb, C., Green, M. S., Eds.; Academic Press: London, 1972; Vol. 6, Chapter 7.) for a detailed description of the RSRG method and to G. D. Mahan and F. H. Claro (Mahan, G. D.; Claro, F. H. Phys. Rev. B 1977, 16, 1168) for applications of the RSRG method to lattice gas models.
1 (Φ + Φ-1-1-1-1 - 2Φ1-11-1) 32 1111
d1 )
1 (Φ + Φ-1-1-1-1 + 2Φ1-11-1 - 4Φ11-1-1) 64 1111
t1 )
1 (Φ - Φ-1-1-1-1 - 2Φ111-1 + 2Φ1-1-1-1) 64 1111
(26)
If one can carry out the partial summation in eq 25, then the real-space renormalization transformation of the
1 (Φ - Φ-1-1-1-1 + 2Φ111-1 - 2Φ1-1-1-1) 16 1111 k1 )
SR ) Const
HI({SR}) + g ) ln{
(28)
f1 )
1 + Φ-1-1-1-1 + 2Φ1-11-1 + 4Φ11-1-1 (Φ 64 1111 4Φ111-1 - 4Φ1-1-1-1)
g)
1 + Φ-1-1-1-1 + 2Φ1-11-1 + 4Φ11-1-1 + (Φ 64 1111 4Φ111-1 + 4Φ1-1-1-1) (29)
Any thermodynamical state of the system can be unequivocally represented as a point p(h,k,d,t,f) in the fivedimensional space of the spin Hamiltonian interaction parameters. An RSRG transformation may be considered as a jump from the initial point P to a point P1 )
5888 Langmuir, Vol. 15, No. 18, 1999
Tarasenko et al.
(h1,k1,d1,t1,f1). Starting from the point P1 as initial, the next RSRG transformation moves the system into point P2. This process can be repeated by an infinite number of times. The sequence of points PP1, P2, ..., Pm, ..., P∞ forms a “trajectory” in the phase hyperspace. As was shown by Nauenberg and Nienhuis,16,19 the free energy of the system F can be evaluated in the series of sequential RSRG transformations of the original Hamiltonian eq 21 ∞
F ) kBT
∑ L-mg(Pm) m)0
(30)
Here Pm ) (hm,km,dm,tm,fm) is the point representing the interaction parameters after the mth RSRG transformation; h0 ) h, k0 ) k, d0 ) d, t0 ) t, f0 ) f. The most important property of any RSRG transformation is the existence of fixed points of the system of Eqs 29. The fixed points are determined by the conditions h1 ) h, k1 ) k, d1 ) d, t1 ) t, f1 ) f. The unstable fixed points of the system eq 29 correspond to the critical points of the Hamiltonian eqs 21. Investigating the behavior of the trajectories in the vicinities of the fixed points, one can construct the global phase diagram of the spin system or the equivalent lattice gas system. In this work we consider only some particular critical points and phase transitions in order to investigate how the different types of interactions between adparticles influence surface diffusion. A fixed point corresponding to the ferromagnetic phase transition due to the pair nearest neighbor attraction is located at Pfc = (0.0, 0.2951, 0.0707, 0.0, 0.0096). The critical surface intersects the k-axis at kc = 0.3955. The relative error compared to the exact value k* ) 0.440 686 8 is about 10%. Analogously, the antiferromagnetic critical point is located at Pafc = (0.0, -0.2906, 0.09648, 0.0, 0.0012). The critical surface intersects the k-axis at kc = -0.4055, and the relative error is about 8%. 4. Results and Discussion 4.1. Nearest Neighbor Interactions in the Ground and Activated State. Using the RSRG transformation described in the preceding section (5 × 4 cluster) we have calculated thermodynamic properties (i.e., adsorption isotherms and adparticle density fluctuations) and the coverage dependences of the chemical diffusion coefficients Dc and Dc′ for different values of φN. We will first consider the influence of nearest neighbor interactions in the ground and activated state. Next nearest neighbor as well as three- and four-particle interactions are neglected. We have calculated the coverage dependence of the chemical diffusion coefficient Dc′ calculated from eq 19, which accounts for the interactions of an activated adparticle with its nearest neighbors. As already mentioned we assume that the interaction energy in the activated state is numerically equal to φN. Coverage dependences ln Dc′(θ)/D(0) are plotted in Figure 2 for attractive (φN < 0, Figure 2a) and repulsive (φN > 0, Figure 2b) interactions, respectively. For comparison the coverage dependences ln Dc(θ)/D(0) obtained from eq 14 in the absence of interactions in the activated state are also plotted as dashed lines. It is quite obvious that attractive interactions in the ground state cause a significant decrease of Dc relative to θ ) 0 where Dc ≡ D(0) (dashed lines of Figure 2a). In contrast, repulsive interactions in the ground state cause a substantial increase of Dc (dashed lines of Figure 2b). Strong repulsion between adparticles causes second-order phase transitions from disordered (19) Nienhuis, B.; Nauenberg, M. Phys. Rev. B 1975, 11, 4152.
Figure 2. Coverage dependences of the normalized chemical diffusion coefficients ln[Dc′/D(0)] (continuous lines) and ln[Dc/ D(0)] (dotted lines) for attractive (a) and repulsive (b) nearest neighbor interactions φN as indicated (in kJ/mol) and φ2N ) φ3 ) φ4 ) 0, T ) 300 K. Arrows indicate the critical points in (b).
lattice gas phase to the c(2 × 2) structure around half coverage. This is clearly seen in the coverage dependences of the mean square density fluctuations plotted in Figure 3. At the critical points the chemical diffusion coefficient has sharp minima (at θ ≈ 0.42 and 0.58) and the coverage dependences are nonanalytical. In ref 5 the same behavior was found to be due to the critical behavior of the thermodynamic factor. At half coverage the chemical diffusion coefficient shows a sharp maximum which is due to the c(2 × 2) ordering.20 In the presence of interactions in the activated state the changes induced by ground state interactions are largely compensated. This is clearly seen for θ f 1, where Dc′ > D(0) in the case of attractive interactions (Figure 2a) and Dc′ < D(0) in case of repulsive interactions (Figure 2b). However, the peculiarities of the chemical diffusion coefficient related to the order-disorder phase transition are still present (Figure 2b); see the tiny minima at the critical points and the small maximum at half coverage. In the limits of θ f 0, 1, jumping adparticles have none or three nearest neighbors, respectively. Therefore, in the absence of interactions in the activated state the limiting values of Dc are given by
lim Dc ) D(0) θf0
lim Dc ) D(0) exp(3βφN) θf1
(31)
If interactions in the activated state are considered, the limiting value of Dc′ for zero coverage is the same, but for (20) Tarasenko, A. A.; Jastrabik, L.; Nieto, F.; Uebing, B. Submitted to Phys. Rev. B.
Diffusion of Interacting Adsorbates
Langmuir, Vol. 15, No. 18, 1999 5889
Figure 4. Adsorption isotherms θ(µ + ) for different values of the next nearest neighbor interactions φ2N (in kJ/mol) as indicated and φN ) φ3 ) φ4 ) 0, T ) 300 K. Note that negative values of φ2N correspond to attractive interactions. The dashed line (φ2N ) 0) represent the noninteracting case (Langmuir gas).
Figure 3. Coverage dependences of the mean density fluctuations ln[ωT/4] for attractive (a) and repulsive (b) nearest neighbor interactions φN as indicated (in kJ/mol) and φ2N ) φ3 ) φ4 ) 0, T ) 300 K. The dashed line (φN ) 0) represents the noninteracting case (Langmuir gas).
θ f 1 the limit is given by
lim Dc′ ) D(0) θf0
lim Dc′ ) D(0) exp(-βφN) θf1
(32)
The differences between eqs 31 and 32 are due to the fact that the ground state of the jumping adatom is modified by mutual interactions with three adjacent adatoms (for θ f 1) while the saddle point is modified by four adjacent adatoms. As a consequence diffusion can be accelerated (inhibited) in the case of attractive (repulsive) interactions. It is certainly interesting to note that interactions in the activated state may change diffusion coefficients by several orders of magnitude without notably influencing the ordering-induced peculiarities at the critical points and at half coverage (Figure 2b). We have used a rather simple model, accounting for the interaction of the activated particle with its surrounding and obtained qualitatively different coverage dependences of the chemical diffusion coefficient. It is intuitively obvious that the model of jumps strongly influences diffusion as the interactions of activated particles strongly influence the activation energy of diffusion Ea. The magnitude of the activation energy is approximately equal to the mean value of the height of the potential barrier, which an activated particle must surmount during its jump from one site to another. The barrier height depends on the interactions of an activated particle with its nearest neighbors; roughly speaking, it depends on the local, shortrange order in the system of adparticles. The short-range order depends on the adparticle density and correlations between nearest neighboring adparticles. Therefore, one should expect that the actual model of jumps strongly
influence the coverage dependence of the diffusion coefficient. However, the diffusion coefficient peculiarities, which are related to the order-disorder phase transitions, should remain unchanged. The latter are affected by the long-range ordering in the system and must be universal for all models of jumps (i.e., these peculiarities should depended on the system’s Hamiltonian and dimensionality only). 4.2. Next Nearest Neighbor Interactions in the Ground State. In this section we will investigate the influence of the next nearest neighbor interactions φ2N on the chemical diffusion coefficient. For the sake of simplicity we neglect all other interactions, i.e., φN ) φ3 ) φ4 ) 0. We also assume that φ2N affects the ground state only. In this case the Hamiltonian eq 5 of the adatom system reduces to the Ising Hamiltonian for two uncoupled sublattices (the open and filled circles in Figure 1). Adsorption isotherms calculated from
θ(µ) )
∂F 1 ∂F ≡ 1+β ∂µ 2 ∂h
(
)
(33)
are shown in Figure 4. For weak interaction (high temperatures) the dependences are close to the Langmuir isotherm (the isotherm for lattice gas without lateral interaction, shown as a dashed line)
θ(µ) ) exp(µ + ) [1 + exp(µ + )]-1
(34)
As interaction increases (temperature decreases) a peculiarity appears at θ ) 0.5. The type of the peculiarity depends on the sign of the interaction parameter φ2N. For repulsive interaction the peculiarity turns out to be a broad, almost horizontal plateau at half monolayer coverage. The plateau corresponds to an ordered p(2 × 1) structure, which consists of vertical or horizontal alternating empty and filled rows. The order-disorder phase transition is of second order and correspdonds to the phase transition in an Ising model with antiferromagnetic interaction between spins. For attractive interaction (φ2N < -3 kJ/mol) the peculiarity turns into vertical jump of adatom density at zero magnetic field h. It is easy to see that the mean value of adatom density θ related linearly to the spontaneous magnetization 〈si〉. Therefore, the jump corresponds directly to the first-order transition in a ferromagnetic Ising model.
5890 Langmuir, Vol. 15, No. 18, 1999
Tarasenko et al.
Figure 5. Coverage dependences of the normalized chemical diffusion coefficient ln[Dc/D(0)] for different values of φ2N (in kJ/mol) as indicated and φN ) φ3 ) φ4 ) 0, T)300 K. The dashed line (φ2N)0) represent the noninteracting case (Langmuir gas).
Figure 7. Adsorption isotherms (a) for different values of φ3 (in kJ/mol) and (b) φ4 (in kJ/mol) as indicated. All other interactions are ignored. Positive values correspond to repulsive interactions. Figure 6. Coverage dependences of the mean density fluctuations ln[χT/4] for different values of φ2N (in kJ/mol) as indicated and φN ) φ3 ) φ4 ) 0, T ) 300 K. The dashed line (φ2N ) 0) represents the noninteracting case (Langmuir gas).
The coverage dependences of the chemical diffusion coefficient (eq 14) are shown in Figure 5. It is quite obvious that in the noninteracting limit, i.e., for θ f 0, Dc is not affected by the interactions. However, for nonzero values of θ the general trend already discussed in section 4.1 is clearly visible. Attractive interaction tends to decrease Dc while repulsive interaction increases Dc. From that point of view nearest and next nearest neighbor interactions behave at least qualitatively similar. However, the most interesting observation is the absence of peculiarities of Dc induced by the order-disorder phase transition. Even in case of strong repulsive interactions the Dc(θ) graphs exhibit neither minima at the critical coverages nor maxima at half coverage. It may be connected with the specific type of the ordered structure p(2 × 1). Figure 6 shows the coverage dependence of the mean square adparticle density fluctuations
〈(ni - θ) (nj - θ)〉 ) χT/4 ≡
β ∂2F 4 ∂h2
(35)
It is quite obvious that the coverage dependences of density fluctuations (and the system phase diagram) are symmetrical around θ ) 0.5 due to the particle-hole symmetry of the Hamiltonian eq 5 (in the absence of three and four particle interactions, φ3 ) φ4 ) 0). At low temperatures χT(θ) exhibits a deep and narrow minimum at half
monolayer coverage and low temperatures but remains analytical at this point. This minimum, which is not present in case of attractive next nearest neighbor interactions, is due to the suppression of density fluctuations in the well-ordered p(2 × 1) lattice gas phase at low temperatures. Density fluctuations (i.e., the displacement of an adatom from its “right” position in the filled adatom rows to a site in the empty rows) increases the free energy of the system and is thermodynamically unfavorable. At coverages q * 0.5, there are density fluctuations which are stabilized by the nonstoichiometry. These fluctuations do not require additional energy for their existence and cannot be removed from the system by adatom jumps. Therefore χT increases when θ slightly deviates from the half monolayer coverage. The coverage dependences of density fluctuations (and the system phase diagram) are symmetrical around θ ) 0.5 due to the particle-hole symmetry of the Hamiltonian eq 5. It is valid for pair interactions between adparticles only. The three-, four-, and other multiparticle interactions do not maintain this symmetry, and coverage dependences of physical quantities (and phase diagrams also) would be asymmetrical as is demonstrated in the following section. 4.3. Multiparticle Interactions in the Ground State. In this section we will investigate how three- and four-particle ground-state interactions influence adsorption isotherms, chemical diffusion, and adparticle density fluctuations. To differentiate between contributions from different interactions we “turn off” all other interactions except φ3 or φ4. The results of our calculations are shown in Figures 6a, 7a, and 8a for three-particle interactions (φ4 ) 0) and in Figures 6b, 7b, and 8b for four-particle interactions (φ3 ) 0), respectively.
Diffusion of Interacting Adsorbates
Figure 8. Coverage dependences of the normalized chemical diffusion coefficient ln[Dc/D(0)] (a) for different values of φ3 (in kJ/mol) and (b) φ4 (in kJ/mol) as indicated. All other interactions are ignored. Positive values correspond to repulsive interactions. The arrows indicate tiny maxima at θ ) 3/4.
Adsorption isotherms are shown in Figure 7. In case of repulsive three- or four-particle interactions, φ3,4 > 0, the isotherms exhibit large plateaus at θ ) 3/4, which are attributable to a p(2 × 2) ordered arrangement of holes. In addition, for strong repulsive three-particle interactions, φ3 > 0, a small plateau is also visible at θ ) 1/2 (Figure 7a). This plateau corresponds to a c(2 × 2) ordered lattice gas phase. For attractive interactions, φ3,4 < 0, the isotherms became steeper and at low temperatures vertical jumps may appear, which correspond to first-order phase transition similar to the case of strong pair attraction between adparticles. It is probably interesting to note that the adsorption isotherms are not symmetrical. This is clearly seen in Figure 7 especially for strongly repulsive interactions. The asymmetry is due to the multiparticle interactions which do not maintain the symmetry of particle and holes. In Figure 8 we show the coverage dependence of the chemical diffusion coefficient Dc calculated from eq 14. The diffusion coefficient grows rapidly as θ f 3/4 and exhibits a tiny maximum at the stoichiometric density corresponding to the well-ordered p(2 × 2) lattice gas phase. This behavior is quite similar to the case of pairwise nearest neighbor repulsions at half coverage. It should be noted that strong three-particle repulsions cause significant increases of Dc around half coverage, i.e., when the c(2 × 2) lattice gas phase is formed. Four-particle repulsion causes very slow growth of the diffusion coefficient as density increases from zero to θ ) 0.75 as the number of elementary squares composed of adparticles is very small for strong lateral repulsion. As θ f 0.75 the number of the squares increases rapidly causing steep growth of the diffusion coefficient. Three-
Langmuir, Vol. 15, No. 18, 1999 5891
Figure 9. Coverage dependences of the mean density fluctuations ln[χT/4] (a) for different values of φ3 (in kJ/mol) and (b) φ4 (in kJ/mol) as indicated. All other interactions are ignored. Positive values correspond to repulsive interactions.
and four-particle attraction inhibits surface diffusion. The dependences for the three- and four-particle interactions are very similar with respect to each other. The formation of p(2 × 2) and c(2 × 2) ordered lattice gas phases is also reflected by the corresponding coverage dependences of the adparticle density fluctuations (Figure 9), which exhibit minima at θ ) 3/4 (corresponding to the p(2 × 2) phase) and θ ) 1/2 (corresponding to the c(2 × 2) phase for φ3 > 0, Figure 9a). As already mentioned, these minima are due to the strong suppression of adparticle density fluctuations at low temperatures (or strong interactions), when the ordering is almost perfect. The coverage dependences of the adparticle density fluctuations for three- and four-particle interactions are also similar to each other. They are asymmetrical in the absence of particle-hole symmetry of the corresponding Hamiltonian. Three- and four-particle attraction increases considerably the density fluctuations similar to the case of pair attraction between nearest neighbor adparticles. 5. Summary We have investigated how pairwise nearest (φN) and next nearest neighbor interactions (φ2N) as well as threeand four-particle interactions (φ3,4) influence the thermodynamical properties and the surface diffusion of adatoms on a square lattice. In addition, we have analyzed the influence of interactions between activated adparticles in the saddle point of the potential and their surrounding. Using a RSRG transformation with a cluster of four blocks containing five site spins each, we have calculated the free energy and other physical quantities for various sets of interactions. We have obtained adsorption isotherms as well as coverage dependences of the chemical diffusion
5892 Langmuir, Vol. 15, No. 18, 1999
Tarasenko et al.
coefficient and of adparticle density fluctuations. The multiparticle interactions significantly change the phase diagram of the system. The formation of an ordered p(2 × 2) structure of holes is possible for strong repulsive three- or four-particle interactions at θ ) 3/4.
meinschaft, by the International Association for the promotion of cooperation with scientists from the New Independent States of the former Soviet Union INTAS-96-0533, and by the Fonds der chemischen Industrie (FCI).
Acknowledgment. This work has been supported by the Heisenberg program of the Deutsche Forschungsge-
LA981097U
