Annihilation and Coagulation Reactions in Low-Dimensional

We also analyze the steady-state regime of annihilation and coagulation ..... time t, divided by the total number n(t) of particles present on the per...
0 downloads 0 Views 798KB Size
Download PDF
Langmuir 1996, 12, 61-69

61

Annihilation and Coagulation Reactions in Low-Dimensional Substrata: Effects of Probability of Reaction and Short Range Interactions† M. Hoyuelos and H. O. Ma´rtin* Departamento de Fı´sica, Facultad de Ciencias Exactas y Naturales, Universidad Nacional de Mar del Plata, Funes 3350, 7600 Mar del Plata, Argentina Received October 25, 1994. In Final Form: February 27, 1995X We study the bimolecular coagulation A + A f A and annihilation A + A f 0 reactions with diffusion and probability p of reaction in one dimension and annihilation reaction with diffusion and short range interactions in two-dimensional fractal and nonfractal percolation clusters. In these cases there is no input of particles and the particle density ρ decays as a function of time t. The effect of p and interactions become relevant at short times. In one dimension analytical approximations of the density for coagulation and annihilation reactions are found which agree very well with Monte Carlo results of ρ(t) for all times and for small values of p. For percolation clusters and large values of the nearest-neighbor repulsive interaction U between particles, we find that a simple mean-field approximation works at short times. The length of the interval where this approximation holds increases as U increases. For the case of repulsive nearest-neighbor with attractive next nearest-neighbor interactions we find an exponential decay of the density at short times. We also analyze the steady-state regime of annihilation and coagulation reactions of immobile reactants with reaction probability and input of particles in a one-dimensional lattice. Analytical approximations of the particle density and of the n-particle correlations are obtained. The limit of zero input rate is solved exactly. The analytical results were confirmed by Monte Carlo simulations.

I. Introduction Recently much effort has been dedicated to the study of diffusion-reaction systems; see, e.g., refs 1-3. This is mainly due to the anomalous behavior that appears when the diffusion occurs on nonhomogeneous substrata (e.g., fractals4 and multifractals5 or on one-dimensional systems).6 This behavior cannot be predicted by mean-field approximations. We consider the annihilation (A + A f 0) and coagulation (A + A f A) reactions. The particles diffuse independently and react instantaneously and irreversibly when two of them meet. In the mean-field approach, one assumes that the law of mass-action holds, i.e., that the variation of the particle density ρ with time t is proportional to ρ2. The process is then described by the classical textbook second-order † Presented at the symposium on Advances in the Measurement and Modeling of Surface Phenomena, San Luis, Argentina, August 24-30, 1994. X Abstract published in Advance ACS Abstracts, January 1, 1996.

(1) Blumen, A.; Klafter, J.; Zumofen, G. Optical Spectroscopy of Glasses; Zschokke, I., Ed.; Reidel: Dordrecht, 1981; p 199. Kuzovkov, V.; Kotomin, E. Rep. Prog. Phys. 1988, 51, 1479. Ben-Avraham, D.; Burschaka, M. A.; Doering, C. R. J. Stat. Phys. 1990, 60, 695. Toussaint, D.; Wilczek, F. J. Chem. Phys. 1983, 78, 2642. Ra´cz, Z. Phys. Rev. Lett. 1985, 55, 1707. Doering, C. R.; Ben-Avraham, D. Phys. Rev. Lett. 1989, 62, 2563. Jiang, Z.; Ebner, C. Phys. Rev. A 1990, 41, 5333. Argyrakis, P.; Kopelman, R. Phys. Rev. A 1990, 41, 2114. Argyrakis, P.; Kopelman, R. Phys. Rev. A 1990, 41, 2121. Wio, H. S.; Rodrı´guez, M. A.; Briozzo, C. B.; Pesquera, L. Phys. Rev. A 1991, 44, R813. Larralde, H.; Araujo, M.; Havlin, S.; Stanley, H. E. Phys. Rev. A 1992, 46, 855. Privman, V. Phys. Rev. E 1994, 50, 50. (2) Anacker, L. W.; Kopelman, R. J. Chem. Phys. 1984, 81, 6402. Kopelman, R. J. Stat. Phys. 1986, 42, 185; Philos. Mag. B 1987, 56, 717; Science 1988, 241, 1620. Newhouse, J. S.; Kopelman, R. J. Chem. Phys. 1988, 92, 1538. (3) Kang, K.; Redner, S. Phys. Rev. A 1985, 32, 435. (4) Meakin, P.; Stanley, H. E. J. Phys. A 1984, 17, L173. Newhouse, J. S.; Argyrakis, P.; Kopelman, R. Chem. Phys. Lett. 1984, 107, 48. Anacker, L. W.; Kopelman, R. Phys. Rev. Lett. 1987, 58, 289. Cle´ment, E.; Sander, L. M.; Kopelman, R. Phys. Rev. A 1989, 39, 6472. Argyrakis, P.; Kopelman, R. Phys. Rev. A 1990, 41, 2114. Sheu, W.; Lindenberg, K.; Kopelman, R. Phys. Rev. A 1990, 42, 2279. (5) Albano, E.; Ma´rtin, H. O. Phys. Rev. A 1989, 39, 6003. Ma´rtin, H. O.; Albano, E. Z. Phys. B 1990, 80, 147. (6) Albano, E.; Ma´rtin, H. O. J. Phys. Chem. 1988, 92, 3594.

0743-7463/96/2412-0061$12.00/0

reaction dρ/dt ) -kρ2, where k is a constant. This approximation is valid in the reaction-limited case, when the reaction rate is sufficiently low and many collisions occur before a particle reacts. One has ρ(t) ∼ t-1. On the other hand, if the diffusion time is much larger than the reaction time, then diffusion governs the process and we have the diffusion-limited case. In this case the density behaves as3,4,7

ρ ∼ t-γ,

γ ) min(ds/2,1)

(1)

where ds is the spectral dimension8 of the substratum. For some fractal structures ds < 2 and for d dimensional Euclidean lattices, ds ) d. For ds < 2, the reaction order x ) 1 + 2/ds (dρ/dt ∝ ρx) is anomalous (x > 2). We modify the model so that when two particles collide they react with probability p (0 < p e 1)9 (for a related model see ref 10) or the particles diffuse coupled by nearestneighbor (NN) and next-nearest-neighbor (NNN) adimensional interactions U (U ) U′/kBT) and W (W ) W′/ kBT), respectively. In many physical and chemical processes, the reaction does not take place instantaneously and the species can collide many times before the reaction occurs. This can be the case when the reaction depends on the orientation of the molecules or when the species must overcome an effective energy barrier to react. The motivation for the present work is not to reproduce or predict any specific experiment but to analyze the general effects of the introduction of the interactions or probability of reaction present in almost all catalytic reactions. Nevertheless there are very few experimental results concerning the influence of short range potentials. (7) Zumofen, G.; Blumen, A.; Klafter, J. J. Chem. Phys. 1985, 82, 3198. Havlin, S.; Ben-Avraham, D. Adv. Phys. 1987, 36, 695. Zumofen, G.; Klafter, J.; Blumen, A. Phys. Rev. A 1991, 43, 7068. (8) Alexander, S.; Orbach, R. J. Phys. Lett. (Paris) 1982, 43, L625. Rammal, R.; Toulouse, G. J. Phys. Lett. (Paris) 1983, 44, L13. (9) Braunstein, L.; Ma´rtin, H. O.; Grynberg, M.; Roman, H. E. J. Phys. A 1992, 25, L255. (10) Privman, V.; Nielaba, P. Europhys. Lett. 1992, 18 (8), 673.

© 1996 American Chemical Society

62

Langmuir, Vol. 12, No. 1, 1996

Hoyuelos and Ma´ rtin

In the present work we analyze the behavior of ρ(t) for the reactions A + A f A and A + A f 0 with probability p of reaction in a one-dimensional lattice and for the reaction A + A f 0 on two-dimensional percolation clusters with interactions U and W between particles. We use the percolation cluster on a square lattice with the lattice covering probability q > qc, where qc = 0.593 is the critical probability (see, e.g., ref 11). If q ) qc, the percolation cluster is fractal and ds = 1.32 (see ref 12 and references cited therein), and if qc < q < 1, it is not fractal, but it is still a disordered structure. All these cases are interesting from the experimental point of view, because they are related to actual situations in heterogeneous structures (see, e.g., ref 2 and references cited therein). On the other hand, probably the surfaces of most materials are fractal on a molecular scale (see ref 13 and references cited therein) and percolation is a simple model to describe fractal substrata. At short times the density is large and the decay of ρ(t), as a function of time t, strongly depends on p or on the short range interactions. At long times the density is low, the particles are spread out, and we recover the behavior of eq 1 regardless of p or interactions. Specifically we expect that in the asymptotic regime (t f ∞)

ρ(t) ∼ t-1/2

(2)

for one dimension and

ρ(t) ∼ t-0.66, ρ(t) ∼ (t/ln t)-1,

q ) qc q > qc

(3)

for percolation clusters. For q > qc the percolation cluster has dimension two, where a logarithmical correction to eq 1 appears.14 In the previous cases there is no input of particles and the aim is to analyze the influence of probability of reaction and interactions and to obtain the approximate behavior of ρ(t) at short times. Another interesting case is when the system is surrounded by a gas of particles at constant pressure, so, there is a constant input rate of particles adsorbed into the system. We study the steady-state of annihilation and coagulation reactions between immobile reactants with a constant input rate of particles in a one-dimensional lattice. It is considered that when two particles are at first neighbors they react with probability p (related works can be found in, e.g., refs 15-17). Very recently17 analytical results have been obtained for annihilation and coagulation reactions between immobile reactants with probability p of reaction. Now we introduce the input of particles into the model and study the steady-state regime. The input changes the behavior of the system, even when the input rate per lattice site tends to zero. As we (11) Stauffer, D. Phys. Rep. 1979, 54, 1; Introduction to the Percolation Theory; Taylor and Francis: London, 1985. Stauffer, D.; Coniglio, A.; Adam, M. Adv. Polym. Sci. 1982, 44, 103. Herrmann, H. J. Phys. Rep. 1986, 136, 153. Feder, J. Fractals; Plenum Press: New York, 1988. Stauffer, D.; Aharony, A. Introduction to Percolation Theory, 2nd ed.; Taylor and Francis: London, 1992. (12) Bunde, A., Havlin, S., Eds. Fractals and Disordered Systems; Springer Verlag: Berlin, 1991; p 130. (13) Farin, D.; Avnir, D. J. Phys. Chem. 1987, 91, 5517. (14) See, e.g., Meakin and Stanley in ref 4. (15) Kenkre, V. M.; Van Horn, H. M. Phys. Rev. A 1981, 23, 3200. (16) Schno¨rer, H.; Kuzovkov, V.; Blumen, A. Phys. Rev. Lett. 1989, 63, 805. Oshanin, G.; Burlatsky, S.; Cle´ment, E.; Graff, D. S.; Sander, L. M. J. Phys. Chem. 1994, 98, 7390. (17) Majumdar, S. N.; Privman, V. J. Phys. A: Math. Gen. 1993, 26, L743.

will see later, the results obtained for this case are quite different from the case with no input. In this article we summarize our previous results of annihilation and coagulation reactions (see refs 18-20). Some new results are presented in sections II-D and IIID. The paper is organized as follows. Section II is devoted to coagulation and annihilation reactions with diffusion and probability p of reaction in one dimension. In section II-A we present the model and the Monte Carlo simulation. In section II-B we obtain the exact rate equation. In section II-C we analyze the behavior of this rate equation and find that, for p , 1, the reaction order x is x ) 2 at short times and x ) 3 at long times. From the analysis of the crossover between these two regimes, a scaling function is obtained. By use of this scaling, analytical approximations of ρ(t) for coagulation and annihilation reactions and for all values of p are presented in section II-D. The agreement between these approximations and the Monte Carlo data is good for all p. In section III we study the annihilation reaction A + A f 0 with diffusion and short range lateral interactions in fractals and nonfractal heterogeneous percolation clusters. The model and the Monte Carlo simulation are described in section III-A. The results and discussion are presented in sections III-B, III-C, and III-D. Specifically, in section III-B we present the general comments, in section III-C the case with repulsive NN interaction (U > 0 and W ) 0), and in section II-C the case with repulsive NN and attractive NNN interactions (U > 0 and W < 0). In section III-C an approximate rate equation of second order is proposed. From this equation one obtains ρ(t) which works reasonably well at short times for large values of U. From the crossover from this mean-field approximation of ρ(t) at short times to the asymptotic behavior at very long times an approximate scaling function is found. The collapse of the Monte Carlo data for different values of U gives support to this approximation. Let us consider the case where the attractive NNN interaction is present. Due to this interaction groups of particles separated by NNN distances appear at short times. In section III-C we analyze the formation of these groups and a crude approximation of the rate equation of the first order is proposed. From this approximation and exponential decay of ρ(t) is obtained which is confirmed by Monte Carlo data. Then, the influence of attractive NNN interaction is to shift the reaction order from x ) 2 (for W ) 0) to x ) 1 (for W < 0) at short times. In section IV we study the steady-state regime of annihilation and coagulation reactions between immobile reactants in a one-dimensional lattice with probability of reaction and input of particles. The model and the equation system for the density and the n-particle correlation Cn are deduced in section IV-A. The approximate solutions of these equations in the steady-state regime are obtained in section IV-B. The numerical data support these analytical results. Finally, in section V we state our general conclusions. II. Coagulation and Annihilation Reactions in One Dimension A. The Model and the Monte Carlo Simulation. We start discussing the case of coagulation reaction A + A f A. The case A + A f 0 will be discussed at the end of section II-D. (18) Hoyuelos, M.; Ma´rtin, H. O. Phys. Rev. E 1993, 48, 71; Phys. Rev. E 1993, 48, 3309. (19) Hoyuelos, M.; Ma´rtin, H. O. Phys. Rev. E 1994, 50, 600. (20) Hoyuelos, M. Phys. Rev. E 1994, 50, 2597.

Annihilation and Coagulation Reactions

Langmuir, Vol. 12, No. 1, 1996 63

In the model the particles perform a random walk between NN sites of a one-dimensional lattice of length L ) 100 000 with periodic boundary conditions. At t ) 0 each site of the lattice is occupied by a particle with probability ρ0 (ρ0 e 1), which is the initial particle density. After that the diffusion starts. At each Monte Carlo step, one o the n(t) particles present in the lattice at time t, randomly chosen, attempts to jump to any of the NN sites with equal probability 1/2. The following situations may appear: (i) if the chosen site is empty, the particle jumps. (ii) If the chosen site is occupied by another particle, they react with probability p. If successful, the selected particle is removed from the lattice and the number of particles is decreased by one, n(t) f n(t) - 1. (iii) Otherwise the jump is not performed and the selected particle remains at its position. This means that the selected particle jumps to the chosen site, collides with the other particle present at this site, and is reflected back with probability (1 - p) to its original position. In the simulations, a time interval equal to 1 is defined as the time needed for the n(t) particles to have, on average, one chance to jump. Specifically, after each jump attempt, time t is increased by 1/n(t). We define the density ρ(t) as the number of particles per lattice site at time t and is obtained by averaging over many (typically 20-30) samples. B. The Rate Equation. Let us denote by P(n f n 1) the probability for a reaction to take place in a Monte Carlo step. The change of the density in a reaction is δρ ) -1/L and the time increases by δt ) 1/n. Then one has dρ/dt ) (δρ/δt)P(n f n - 1). Now let us consider a given pair of NN occupied sites. The probability for one of these two particles to be selected at random in a Monte Carlo step is 2/n. The probability for the selected particle to react with the other one is 1/2p. Then P(n f n - 1) ) (2/n)(1/2)pn1, where n1 is the number of NN pairs present in the system at time t. Finally one obtains the exact rate equation

dρ/dt ) -pρ2

dρ/dt ) -pΓ1

(4)

where Γ1 ≡ n1/L is the number of pairs of NN occupied sites per lattice site. C. The Scaling Function. For p ) 1, due to the reaction, the probability of finding a particle near another one is less than ρ and Γ1 < ρ2. But for small values of p there is a mixing effect. A particle must collide many times before the reaction takes place. For perfect mixing the probability of finding a particle is ρ independently of if there is, or is not, another particle close to the first one. Then we expect Γ1 ) ρ2 at very short times. At very long times also the particles must collide many times before they react. But, because the density is very small, the distance between particles is so long that the mixing does not occur. We have numerically checked that, for very long times, ρ(t) is independent of p (see Figure 2 below). Then, for all p, dρ/dt ) -πρ3/2, which is the known result21 for p ) 1 and t f ∞. Then, using eq 4, when t f ∞, Γ1 ) πρ3/2p. The intersection of this asymptotical behavior with Γ1 ) ρ2 at short times occurs at the crossover density ρc, given by

ρc )

2 p π

(5)

In summary, it is expected, for small values of t (large ρ) (21) Doering, C. R.; Ben-Avraham, D. Phys. Rev. A 1988, 38, 3035.

(6)

the mean-field behavior (second-order reaction). For large t (small ρ) one has

dρ/dt ) -πρ3/2

(7)

the anomalous behavior (third-order reaction, ρ ∼ t-1/2, see eq 2). The crossover density ρc (eq 5) from mean-field to anomalous behavior decreases as p decreases. In order to analyze the universal (independent of p) behavior of Γ1 as a function of ρ, we propose the following ansatz

{

Γ1 x if x . 1 ) ρf(ρ/ρc), where f(x) ) 2 ρc x if x , 1

(8)

Monte Carlo data, shown in Figure 1, support this scaling ansatz, especially when ρ0 . ρc ) 2p/π. D. An Approximation for ρ(t). A simple function which fulfills the scaling form is

f(x) )

x2 x+1

(9)

(i.e. Γ1 ) ρ3/(ρ + ρc)). We will use f as an approximation to the scaling function (see Figure 1). Integrating eq 4, one obtains the approximate density ρa,

ρa(t) )

ρ0 1 + x(1 + ρc/ρ0)2 + 4p2t/π 2 1 + ρc/2ρ0 + ρ0pt

(10)

In Figure 2 we present log ρa and the Monte Carlo data of log ρ versus log t for different values of p. The agreement between ρa and ρ is reasonably good, even for cases p ) 1 and p ) 0.1 where the condition ρ0 . ρc is not fulfilled. Let us note that the logarithmic scale used is more appropriate than the linear scale for ρ in order to analyze the difference between ρa and ρ for a fixed value of t. Let us now analyze eq 10 as a function of time t. Assuming ρ0 . ρc ) 2p/π we have

ρa )

ρ0 ρ0 , ) 1 + ρ0pt 1 + t/t0

for t , t1 )

π 4p2

(11)

the mean-field approximation, and

ρa )

t-1/2 , xπ

for t . t1

(12)

the anomalous behavior. t1 is the crossover time between these behaviors and t0 ) 1/pρ0 , t1 is a characteristic time in the mean-field regime; log10(ρa) as a function of log10t remains almost constant for t < t0. It is possible to use other approximations to the scaling function f. For example, f(x) ) x3/2/(x + 1/x + 2b)1/2 or f(x) ) x3/2/(xk + x-k)1/2k, with x ) ρ/ρc, and where b and k are adjustable parameters. Nevertheless we prefer to use eq 9 due to its simplicity and because using this equation one can obtain a close form for the approximate density ρa (see eq 10). Moreover, the ρa obtained is a good approximation for the true ρ. Let us finally discuss the case of annihilation reaction A + A f 0. In this case, when two particles meet, both are removed from the lattice and the “exact” rate equation is dρ/dt ) -2pΓ1 (compare with eq 4). For p , 1 we expect Γ1 ) ρ2 at very short times. For p ) 1 it is known that

64

Langmuir, Vol. 12, No. 1, 1996

Hoyuelos and Ma´ rtin

Figure 1. Scaling of log10(Γ1/ρρc) versus log10(ρ/ρc) for different values of p: p ) 1 (circles); p ) 0.1 (squares); p ) 0.01 (full triangles). The respective values of ρ0/ρc are 1.26, 12.6, and 126. The line represents the approximation to the scaling function, see eq 9.

Figure 2. log10 ρ versus log10 t for the coagulation reaction A + A f A in one dimension, for different values of p: p ) 1 (circles); p ) 0.1 (squares); p ) 0.01 (full triangles). The points are numeric results and the curves are the approximation of eq 10. A value of ρ0 ) 0.8 was used for all curves.

ρ(t) ) t-1/2/(2xπ)22 in the asymptotic regime. Then, dρ/dt ) -2πρ3 for all p and t f ∞. Using the same procedure discussed above for A + A f A, we obtain the following approximate density ρ′a for the annihilation reaction

ρ0 1 + x(1 + ρc/ρ0)2 + 4p2t/π ρ′a(t) ) 2 1 + ρc/2ρ0 + 2ρ0pt

(13)

with ρc ) p/π. The agreement between ρ′a(t) and Monte Carlo data of A + A f 0 is similar to the agreement between ρa(t) and numerical results of A + A f A shown in Figure 2. III. The Annihilation Reaction in Percolation Clusters A. The Model and the Monte Carlo Simulation. We have considered the reaction between particles dif(22) See, e.g., Ben-Avraham, Burschaka, and Doering in ref 1.

Figure 3. Density ρ vs time t in log10-log10 scales for annihilation reactions on fractal percolation clusters (q ) qc ) 0.593) on a square lattice of linear size L ) 350, ρ0 ) 0.2 and for different short-range interactions: U ) W ) 0 (O); U ) 4, W ) 0 (0); U ) 4, W ) -3 (+). The straight line has slope -0.66.

fusing on a two-dimensional percolation cluster. In the simulation, the standard percolation model (see, e.g., ref 11) was used on L × L (L ) 350 and 700) square lattices. Each site of the lattice can be occupied with probability q or empty with probability (1 - q). A cluster is defined as a group of occupied sites connected by NN distances. A percolation cluster is a cluster which has its length and width equal to L (a snapshot picture of a percolation cluster is shown in Figure 6 below). In the model the particles perform random walks between NN sites of percolation clusters. Periodic boundary conditions were used. The particle density ρ is defined as the number of particles per cluster site and the reaction takes place between NN occupied sites. At t ) 0 the lattice is filled at random with a particle density ρ0, avoiding pairs of NN particles. After that the diffusion starts. The energy of a given configuration can be written as H/kBT ) 1/2∑ij′Uninj + 1/2∑ij′′Wninj, where ni is the occupation number of site i; ni ) 1 if site i is occupied and ni ) 0 otherwise. Primed sum and double primed sum symbols denote sums over NN and NNN sites, respectively. At each Monte Carlo step one of the n(t) particles, randomly chosen, attempts to jump to any of the four NN sites with equal probability. The following situations may occur: (i) The chosen site does not belong to the cluster and the jump is not performed. (ii) The chosen site belongs to the cluster, the jump is performed with probability w ) min[1, exp(-∆H/kBT)], where ∆H is the energy change in the movement of the selected particle to the chosen site. If the jump is performed and the chosen site is NN of a second particle, both particles react and n(t) is reduced to n(t) - 2. After each jump attempt, time is increased by 1/n(t). This model is known as the blind ant, because particles try to jump to any of the four NN sites regardless of whether they belong to the cluster or not. There is another model, the myopic ant, in which particles try to jump only to sites that belong to the cluster. We also made simulations with the myopic ant model and did not find relevant differences. B. Results. In Figure 3 we plot the time dependence of the density for different values of the interactions for fractal percolation clusters (q ) qc). The case U ) W ) 0 is shown for comparison with other cases.

Annihilation and Coagulation Reactions

Langmuir, Vol. 12, No. 1, 1996 65

In the interval 4 < log10t < 6 our results for U ) W ) 0 can be fitted by a straight line of slope -0.63 ( 0.01. This value (greater than -0.66, see eq 3) is in agreement with other Monte Carlo simulation results.23 The small curvature suggests that the asymptotic behavior will hold for longer times not reached in our simulation. The approach of all curves to the case U ) W ) 0, observed in Figure 3 at t ∼ 106, is an indication that the asymptotic behavior (see eq 3) will be valid for longer times. The asymptotic regime has also been suggested in other related models.9,10 In the next subsection we will analyze the short time regime for each case of interaction. C. With Repulsive NN Interaction. Let us start by considering the case of repulsive NN interaction (U > 0, W ) 0). Due to the repulsive potential, many attempts of reaction occur before a pair of particles finally reacts regardless of the substratum type. Then, at short times (large density and small distances between particles) a mixing effect appears and we expect a random distribution of particles over the substratum with the restriction of no NN occupied sites (remember that the reaction takes place at NN sites). From an analysis similar to the one used in obtaining eq 4 for one dimension, one gets (for more details see ref 19)

-

dρ ) 2e-U ρ2 dt

(14)

Figure 4. log10 R as a function of log10 t in the percolation cluster at q ) qc ) 0.593 for different values of U (and W ) 0): U ) 0 (0); U ) 2 (×); U ) 4 (4); U ) 6 (*). There are also plotted the results for q ) 0.65 for the interactions U ) 0 (O) and U ) 4 (+). For the sake of clarity, the results for q ) 0.65 have been shifted on the ordinate axis 4 higher. The straight lines correspond to the approximation given by eq 15. A value of ρ0 ) 0.2 and a square lattice of linear size L ) 700 were used.

but now eq 14 is not exact and was obtained assuming that eU . 1, a random distribution of particles and using some simple configurations of particles. The structure of the percolation cluster is strongly heterogeneous and the exact expression for -dρ/dt is very difficult to obtain. The alternative that we use is to propose eq 14 as a first approximation to the problem and to check how well it works. We will use eq 14 for all q (qc e q e 1). Integrating eq 14 we have

R≡

(

)

1 1 ρ ρ0

-1

)

1 2e-Ut

(15)

In Figure 4 we plot log10 R versus log10t for different values of interaction U and q ) qc. The straight lines correspond to eq 15; the agreement with the numerical results at short times is apparent. At very short times there is a small deviation of the numerical results with respect to the lines. This is due possibly to the crude approximation used in eq 14. In Figure 4 we also plot the results for q ) 0.65 and different values of U; the results are shifted for clarity. We can see that the approximation of eq 15 for short times also holds in this case. At long times the mixing effect does not appear and we expect that the behavior of eq 3 holds for all U > 0. Let us consider the case q ) qc and compare the mean-field behavior of eq 15 with the anomalous behavior at long times (ρ ∼ t-0.66). We define the crossover time t1 as the time at which the intersection of these behaviors occurs and the crossover density ρc as the density at this time. Assuming 2ρ0 exp(-U)t1 . 1 one has t1 ∼ exp(3U) and ρc ∼ exp(-2U) (for simplicity we have assumed that t-0.66 = t-2/3). As expected, the length of the time interval where the mean-field approximation holds (0 < t , t1) increases as U increases. Defining Rc ) (1/ρc - 1/ρ0)-1 (as ρc , ρ0, Rc ∼ ρc ∼ e-2U), we expect that R/Rc will be a function of t/t1 only, independent of the values of U. In order to (23) Argyrakis, P.; Kopelman, R. J. Chem. Phys. 1984, 81, 1015; J. Phys. A 1988, 21, 2753; Phys. Rev. A 1992, 45, 5814. McCarthy, J. F. J. Phys. A 1988, 21, 3379.

Figure 5. log10(R/e-2U) versus log10(t/e3U), scaling of the data in Figure 4 for q ) qc. The interactions used are U ) 0 (full triangles), U ) 2 (0), U ) 4 (line), and U ) 6 (O).

analyze this universal behavior, in Figure 5 we plot log10 (R/e-2U) against log10(t/e3U). As expected this approximate scaling works reasonably well for large values of U. For small values of U the mixing effect does not appear and eq 15 does not hold, for this reason the data for U ) 0 do not collapse with the results for other cases. Now, comparing the mean-field behavior of eq 15 with the asymptotic behavior ρ ∼ (t/log t)-1 for q > qc, one obtains that the crossover time and the crossover density behave as t1 ∼ exp(eU) and ρc ∼ exp(-eU). Let us note that in this case the crossover is very smooth (from ρ ∼ t-1 to ρ ∼ (t/logt)-1) and time t1 increases very fast with a double exponential in U. We were not able to check this crossover with simulations due to the lack of enough computational capacity. D. With Repulsive NN and Attractive NNN Interactions. Let us consider the case where repulsive

66

Langmuir, Vol. 12, No. 1, 1996

Hoyuelos and Ma´ rtin

Figure 6. NNN structures of particles. A 50 × 50 section of the fractal percolation cluster (q ) qc) on a 350 × 350 lattice is shown. The cluster sites occupied (empty) by particles A are denoted by big (small) circles. The snapshot picture corresponds to t ) 500 with interaction U ) 4, W ) -3, and initial density ρ0 ) 0.2. The particles separated by NNN distance have been joined by lines in order to guide the eyes. One clearly sees that most of the particles belong to NNN structures of particles.

NN (U > 0) and attractive NNN (W < 0) interactions are present. Due to this interaction, the most energetically favorable configuration of two particles is when they form a pair of NNN occupied sites (see Figure 6). A way of characterizing this structure formation is to measure the number of particle pairs at NNN sites. We define θ(t) as the number of pairs of NNN occupied sites per particle present on the percolation cluster at time t. We also define a NNN structure of particles of size s as a group of s particles on the percolation cluster connected by NNN distances. If a particle is isolated, we consider that it is in a NNN structure of size s ) 1. Other way to obtain information about structure formation is measuring the ratio of particles (the relative mass) belonging to NNN structures of size s. Then we define τ(s,t) as the number of particles belonging to all NNN structures of size s, at time t, divided by the total number n(t) of particles present on the percolation cluster. In Figure 7 log10ρ(t) and θ(t) are shown as a function of log10t. The interaction used was U ) 4 and W ) -3. In these figures, θ grows at short times (0 e t < 103 and 0 < t < 103.4, for ρ0 ) 0.05 and 0.2, respectively), indicating structure formation (also see Figure 6). In intermediate times the structures disintegrate. An abrupt decay of θ(t) and log10 ρ(t) takes place (3.4 < log10 t < 4.6 and 4 < log10 t < 4.8, for ρ0 ) 0.05 and 0.2). The decay of θ(t) and log10 ρ(t) is more abrupt for larger values of ρ0. At long times, the structures disappear (θ f 0). In Figure 7b values of θ(t) for the rest of the cases are shown as a function of log10 t. Let us note that only for W < 0, θ(t) increases at short times. For all the other cases θ(t) decreases for all times. In Figure 8 values of τ(s,t) for different times and for cases (a) U ) W ) 0, (b) U ) 4, W ) 0, and (c) U ) 4, W ) -3, are shown. Let us note that for case c, and for s g 6, τ increases when the time increases from t ) 9 to t ) 90 and to t ) 4200. This fact indicates the increase, at short times, of the relative mass of each NNN structure of size s g 6 at the expense of the decrease of the total 5 relative mass of small NNN structures [∑s)1 τ(s,t)]. This is due to diffusion and interaction effects. This fact is only observed for case c, and for all the other cases τ, for s g 2, decreases as t increases for all values of t (see parts

Figure 7. (a) log10 ρ vs log10 t and (b) θ vs log10 t for the interaction U ) 4, W ) -3, q ) qc, (L ) 350), and for the following initial densities: ρ0 ) 0.2 (full triangle), and 0.05 (+). In (b) there are also plotted the following interactions for ρ0 ) 0.2: U ) W ) 0 (0); U ) 4, W ) 0 (O). Note that the case with W ) -3 is the only one in which the value of θ increases in the short time regime.

a and b of Figure 8). As expected for case c at intermediate times (when the abrupt decay of log10 ρ and θ takes place), the NNN structures of s g 2 disappear (see the data of Figure 8c for t ) 48 000). As the reaction takes place at NN sites, it is expected that

-

dρ ∼ Γ2 dt

(16)

where Γ2 is the correlation at NNN sites, i.e., the probability of finding a pair of NNN occupied sites (in eq 4, -dρ/dt ∼ Γ1, and the reaction occurs when two particles occupy the same site). In eq 14, Γ2 was assumed to be equal to ρ2. Let us now consider the value of Γ2 for the present case. From the definitions of θ and Γ2 one has Γ2 ) θρ. The data of Figure 7b show that the increase of θ takes place in a very long interval of time. Then, one can approximate θ by a constant at short times. So Γ2 ∼ ρ and

Annihilation and Coagulation Reactions

Langmuir, Vol. 12, No. 1, 1996 67

Figure 9. log10 ρ versus t in the percolation cluster at qc (L ) 700) for interactions U ) 4, W ) -3 (O) and U ) 4, W ) 0 (full triangles). The upper straight line denotes the linear behavior present when an attractive NNN interaction is used, see eq 17. The lower continuous line is the approximation of eq 15.

is plotted for comparison). We can see at short times (in a range of approximately 3 decades in t) a linear behavior as predicted by eq 17. For long times, the groups of NNN particles disintegrate and the asymptotic behavior without interaction, ρ ∼ t-0.66, is obtained. Similar results are obtained for the percolation cluster at q > qc. IV. Steady-State Regime in Annihilation of Immobile Reactants with Input of Particles

Figure 8. Values of τ, the relative mass, as a function of s, the size of the NNN structure, for different values of time (L ) 350, ρ0 ) 0.2): t ) 0 (b), 90 (×), 4200 (full box), and 48000 (full triangle) (only in part c). The cases plotted are for the interactions (a) U ) W ) 0, (b) U ) 4, W ) 0, and (c) U ) 4, W ) -3. Note, in the last case, the formation of structure of large s as time increases, and the disintegration of these structures for large enough times. There is an interval from τ ) 2.5 to 7.5 not shown in the figures.

we have -dρ/dt ∼ ρ, which means an exponential decay

ρ ∼ e-Rt

(17)

where R is a parameter that depends on the values of U, W, ρ0, and the probability of the percolation cluster q (qc e q e 1). In Figure 9 we plot log10 ρ versus t for U ) 4 and W ) -3 for the percolation cluster at qc (the case with W ) 0

A. The Model and the Equation System. We will describe the calculation for the reaction A + A f A between immobile reactants in a one-dimensional lattice with probability p of reaction and input rate  per lattice site. The method for case A + A f 0 is very similar. In the model the system evolves with a discretized time t. The one-dimensional lattice has L sites with periodic boundary conditions. At each time step, one of the L sites is randomly chosen with equal probability 1/L and time is increased by 1/L. If the site is empty, it is occupied by an A particle with probability . If the chosen site is already occupied by a particle, this particle “looks” at either of the first neighbors with equal probability 1/2, if the neighbor is occupied by another particle, the first one reacts and evaporates with probability p, if not, nothing happens, and the particle remains at its site. Let us denote the occupation number of a generic site i by si. If site i is occupied then si ) 1, otherwise si ) 0. If we have at time t any given configuration {s}, the occupation probability of site i at time t + δt (with δt ) 1/L) is

Pi{s}(t + δt) ) (1 - si-1)(1 - si)(1 - si+1)/L OOO + (1 - si-1)(1 - si)si+1/L OOb + (1 - si-1)si(1 - si+1) ObO + (1 - si-1)sisi+1(1 - p/2L) Obb + si-1(1 - si)(1 - si+1)/L bOO + si-1(1 - si)si+1/L bOb + si-1si(1 - si+1)(1 - p/2L) bbO + si-1sisi+1(1 - p/L) bbb (18)

68

Langmuir, Vol. 12, No. 1, 1996

Hoyuelos and Ma´ rtin

At the right side we show the configuration of sites (i 1, i, i + 1) which corresponds to each term of Pi{s}(t + δt). Symbol b (O) denotes an occupied (empty) site at time t. Simplifying eq 18 and averaging over configurations {s} we obtain

P(t + δt) - P(t) )

  p - P(t) - C2(t) L L L

(19)

Where we have considered that there are not privileged sites in the lattice, so 〈si〉 ) 〈sj〉 ) P for all i, j; and we have denoted 〈sisi+1〉 by C2. The occupation probability per lattice site P(t), now independent of i, is equivalent to the global density of particles, which we will denote by C1(t). Knowing that δt ) 1/L, in the continuous limit, when L f ∞, we have

dC1(t) ) [1 - C1(t)] - pC2(t) dt

(20)

The first term represents the input of particles: the input probability  times the probability of finding an empty site. The second term represents the decrease of C1(t) due to the reaction. In the same manner we can obtain less obvious linear differential equations for the n-particle correlations Cn ) 〈s1s2 ... sn〉

dCn(t) dt

where yn are polynomials in a (except for n ) 0), defined by

n

)

(〈s1 ... sj-1sj+1 ... sn〉 ∑ j)1 Cn) - p[(n - 1)Cn(t) + Cn+1(t)] (21)

for n g 2. If we consider no input of particles, i.e.,  ) 0, then eq 21 is equal to eq 3 in ref 17 for a one-dimensional lattice. For n ) 2 we have dC2(t)/dt ) 2[C1(t) - C2(t)] - p[C2(t) + C3(t)]. For n g 3 we will make an approximation for 〈s1 ... sj-1sj+1 ... sn〉. We can propose as a first step a decorrelation of the product, i.e., 〈s1 ... sj-1sj+1 ... sn〉 = 〈s1 ... sj-1〉〈sj+1 ... sn〉 ) Cj-1Cn-j. This approximation becomes valid when the chain of j - 1 sites is far from the chain of n - j sites. We consider that it could be also a good approximation to take both chains together instead of one far from the other, i.e., 〈s1 ... sj-1sj+1 ... sn〉 = 〈s1 ... sj-1sj ... sn-1〉 ) Cn-1. We will use this last approach because it has the advantage of producing linear differential equations, so the steady state can be easily obtained. Finaly we have

dCn(t) ) n[Cn-1(t) - Cn(t)] - p[(n - 1)Cn(t) + dt Cn+1(t)] for n g 3 (22) In the next section we show numerical simulations that support this approximation. B. The Steady State. To evaluate the steady-state density C1 we have to solve the equation system

aC1 + C2 ) a naCn-1 - (na + n - 1)Cn - Cn+1 ) 0

for n g 2 (23)

where parameter a is equal to /p. After some algebra it can be proved that ∞

Cn ) yn-1

∑

i)n-1

(-1)

i-n+1

Figure 10. Plot of Cn versus a for the reaction A + A f A. From top to bottom n ) 1 (the density), 2, and 4. The lines are evaluated from eq 24. The points were obtained via numerical simulations. In the simulations a time tmax ) 200 was enough to reach the steady state; random initial distributions were used with initial density ρ0 ) 0.5; a lattice with periodic boundary conditions and size L ) 30 000 was used and the results were obtained averaging over 10 samples.

(i + 1)!ai yiyi+1

(24)

yn+1 ) (n + 1)ayn-1 + [(n + 1)a + n]yn

(25)

with y0 ) 1/a and y1 ) 1. For n ) 1 the steady-state density is

C1 )

1

∞

∑

ai)0

(-1)i

(i + 1)!ai yiyi+1

(26)

In Figure 10 we plot Cn as a function of a for different values of n. Although eqs 20 and 21 are valid only for uniform initial particle distributions, the steady state derived from these equations is valid for all kinds of initial conditions, since the input and the annihilation make the particle distribution homogeneous. Therefore, a universal behavior of the Cn’s, in the sense that they are independent of initial conditions, is obtained. Numerical data are also plotted in Figure 10 to confirm eq 24. The agreement between numerical and analytical data supports the approximation made in eq 22 for n g 3. The case with no input of particles, i.e., a ) 0, was solved in refs 15 and 17. If random initial distributions are used, the density behaves as C1(t) ) ρ0 exp[-ρ0(1 - e-pt)], where ρ0 is the initial density. For long times C1 ) ρ0e-ρ0. In our case, when a f 0, C1 f 1/3. For a f 0 the result is completely different than the one for a ) 0. This means that a very small input of particles is sufficient for the system to forget initial conditions and to reach the steadystate density value which only depends on the parameter a. The value C1 ) 1/3 is an exact result independent of the approximation we made for n g 3 because when a f 0, Cn f 0 for n g 2, and the result can be derived from the two first equations of (23) (for n ) 1 and n ) 2) which are exact. Repeating the same process for the annihilation reaction A + A f 0 an expression similar to eq 24 can be obtained for the values of Cn.20 The result for the steady-state density in the limit when a f ∞ is C1 f 1/5. As before, the

Annihilation and Coagulation Reactions

value 1/5 is exact. As expected, this value is lower than the one corresponding to coagulation reaction. V. Conclusions We have analyzed the decay of the density ρ(t) of diffusing particles for two cases: (i) coagulation A + A f A and annihilation A + A f 0 reactions with probability p of reaction in one dimension and (ii) the annihilation reaction with NN and NNN dimensionless interactions, U and W respectively, in two-dimensional percolation clusters. Let us stress that in both cases, the exact results of ρ(t) are not known. For case i we deduced the exact rate equation (4) for a coagulation reaction. From the analysis of the crossover from short time to long time regimes of NN correlation Γ1, the scaling ansatz of eq 8 is proposed. The Monte Carlo data support this scaling, specially when ρ0 . ρc ) 2p/π (see Figure 1). Using a simple function which fulfills the scaling form, we obtain an analytical approximation of the density ρa (see eq 10). This procedure can be easily extended to annihilation reaction and we obtain the approximation ρ′a (see eq 13). The agreement between the approximation and the Monte Carlo data is reasonably good for all times and for different values of p (see Figure 2). For coagulation reaction an alternative analytical diffusion-equation-type approximation24 has been developed to study the crossover between the short and long time regimes. With this method an approximation for ρ(t) is obtained, which improves for large values of p and small initial density ρ0. Our approximate scaling form (eq 8) becomes more accurate in the opposite region, that is, for small values of p and large ρ0. For case ii with repulsive NN potential (U > 0, W ) 0) we obtained the mean-field approximation (14) that yields the behavior given by eq 15. Such behavior is independent of the substratum type and is present at short times (see Figure 4). The lapse during which the effect of the interaction is present increases as the value of the potential U increases. Specifically this time interval extends to t , t1 where t1 behaves as t1 ∼ exp(3U) for q ) qc and t1 ∼ exp(eU) for q > qc. For the case of repulsive NN interaction (U > 0) with attractive NNN interaction (W < 0), it is interesting to remark the formation of groups of particles connected by NNN distances at short times (see Figures 6, 7, and 8) (24) Privman, V.; Doering, C. R.; Frisch, H. L. Phys. Rev. E 1993, 48, 846.

Langmuir, Vol. 12, No. 1, 1996 69

and the abrupt decay of log10 ρ and θ(t) at intermediate times, which means the sudden extinction of NNN structures. Due to the formation of these groups and from eq 16 the approximate exponential decay of eq 17 is obtained at short times for q g qc. This behavior has been checked by Monte Carlo simulations (see Figure 9). Let us comment that, even though case ii is restricted to two-dimensional percolation clusters, we expect that eqs 14, 15, and 17 approximately hold for other lowdimensional substrata. On the other hand, in the shorttime regime, the density is not so low as in the asymptotic regime, where one expects the universal behavior of eq 1. Therefore, from the experimental point of view, the short time regime acquires importance. We hope that our work will encourage new experiments on annihilation reactions between adsorbed particles on disordered media. Experimentally, the random initial deposition can be thought of as the condensation of incident particles onto a cooled substratum. After that, with a sudden increase in temperature, the diffusion starts. Although we do not attempt to predict any specific experiment, it is interesting to comment that, for example, an adimensional interaction U ) 4 at T ) 200 K corresponds to an interaction energy U′ (U′ ) kBTU) of about 6.7 kJ/mol, which is on the same order as that of the adsorbate-adsorbate lateral interaction observed in real systems. For the case of immobile reactants we studied the steady-state regime for the coagulation (A + A f A) and annihilation (A + A f 0) reactions with input of particles (0 f A) in a one-dimensional lattice with probability p of reaction. Approximate analytical results for the n-particle correlations Cn were obtained. In particular, when n ) 1 we have the density of particles. These results are in agreement with Monte Carlo data (see Figure 10). A universal behavior of the Cn’s, in the sense that they are independent of initial conditions, was obtained. In the limit of small input, a f 0, we obtained the exact value of the steady-state density C1 ) 1/3 and C1 ) 1/5 for coagulation and annihilation reactions, respectively. These results are quite different from the cases with no input which was solved in refs 15 and 17. This means that the input changes the behavior of the system, even when the input rate per lattice site tends to zero. Acknowledgment. We wish to acknowledge the Consejo Nacional de Investigaciones Cientı´ficas y Te´cnicas (CONICET) for partial support. LA9408432