J . Phys. Chem. 1985,89, 4758-4761
4758
Diffusion-Controlled Reaction Kinetics on Fractal and Euclidean Lattices: Transient and Steady-State Annihilation L. W. Anacker, R. P. Parson,+and R. Kopelman* Department of Chemistry, The University of Michigan, Ann Arbor, Michigan 48109 (Received: June 3, 1985)
+
-
We derive a simple kinetic rate law for A A P diffusion-controlledbinary reactions on both fractal and Euclidean lattices (one-, two-, and three-dimensional). The same rate law describes the relaxation into the steady state and some of the short-time behavior in the presence of a source. In the low density limit the reaction order X is given in terms of the spectral dimension d,, X = 1 + 2 / d , for d, C 2 but X = 2 for d, > 2 . The relaxation time is proportional to XI.Very large scale simulations of reacting random walkers (A + A A) verify these results on a cubic lattice (d, = 3), a linear lattice ( d , = l), and a Sierpinski gasket (d, = 1.36). High density deviations, obtained by simulations, reveal an opposite trend between fractal and Euclidean reactions. This is discussed in terms of topology and many-body terms. -+
1. Introduction “Heterogeneous reactions present a challenge to several branches of science-chemical kinetics, surface and solid-state physics, and surface chemistry”.’ The standard textbook equations1s2are only applicable to 3-dimensional homogeneous media. Reactions in low-dimensional, disordered, or finite-size media have recently been of much intere~t.~”The solutions to such problems have been either ad hoc or limited to idealized simple models, such as 1-dimensional binary reaction^.^*^-^ These very interesting special cases (whose solutions are far from trivial) have only hinted at the difficulties involved in treating catalytic reactions,10reactions in porous media,” and reactions on or in disordered solids or viscous fluids, etc.12 A more general approach to chemical reactions in disordered media may be given by fractal models, involving fractal spatial and/or fractal time structure^.'^^^^ While this has been realized intuitively?s1’ any serious quantitative approach has had to await the differentiation between the static fractal (Hausdorff) dimension,I3 d f , and the dynamic fractal (spectral) d i m e n ~ i o n , l ~ - ’ ~ d,. While the static fractal dimension yields extrinsic properties (mass, volume, density), the dynamic dimension characterizes the spectrum of the Laplacian operator, and is thus associated with such intrinsic properties as the density of vibrational modes or the solutions of the diffusion equation. It follows that diffusion-controlled processes are described by both d f and d, where generally d, < d f < d ( d is the Euclidean dimension), and only for homogeneous media are they described by a single dimension (d, = d f = 6). Moreover, for natural fractal media, Le., random fractals such as percolation c l ~ s t e r s ’or ~ .diffusion-limited ~~ agg r e g a t e ~one , ~ ~always has d , C 2, i.e., a truly low-dimensional behavior.*O The spectral dimension is relevant to the description of diffusion-limited chemical kinetics in disordered media and has been applied to reactions involving two particles of the same A + A, or of different type^,'^,^^ A + B. One way to measure the spectral dimension is to monitor the random walk of a single particle. Originally, a conjecture was given15-16for the number of distinct sites visited, S , for a random walk in the asymptotic limit: N t ) ) tf f = ds/2 d, < 2
f = l ds>2 (1.1) This has now been widely tested by using Monte Carlo simulations of random walks on deterministic fractals16 and on random fractalsI8 and it is found that eq 1.1 holds within =tl% or better. Equation 1.1 has also worked well for systems with random site energy distributions and their accompanying Boltzmann weighted hopping times where the exponent f is some fractional power used Department of Chemistry, University of Washington, Seattle, WA 98195.
to describe a fractal-like distribution of effective hopping times.’* Equation 1.1 has higher order correction terms which are very important for some systems, e.g., ds = 2. Nonetheless, in practice the asymptotic regime is reached fairly quickly so that (1.1) works remarkably well even at fairly short times.I8
Berry, S.; Rice, S. A,; Ross, J. “Physical Chemistry”; Wiley: New 1980. Noyes, R. M. Prog. React. Kinet. 1961, I, 128. de Gennes, P. G. J . Chem. Phys. 1982, 76, 3316. Torney, D. C.; McConnell, H. M. J . Phys. Chem. 1983, 87, 1441. (5) Hatlee, M. D.; Kozak, J. J. Proc. N u f l . Acad. Sci. U . S . A . 1981, 78, 972. (6) Calef, D. F.; Deutch, J. M. Annu. Rev. Phys. Chem. 1983, 34, 493. (7) Torney, D. C.; McConnell, H. M. Proc. R . Soc. London, Ser. A . 1983, 387, 147. (8) Elyutin, P. V. J . Phys. C 1984, 17, 3606. (9) Agranovich, V. M.; Gallanin, M. D. “Electronic Excitation Energy Transfer in Condensed Matter”; North-Holland: Amsterdam, 1982. (10) Hughes, R. “Deactivation of Catalysts”; Academic Press: London, 1984. (1 1) ‘Physics and Chemistry in Porous Media”; Johnson, D. L., Sen, P. N., Eds.; American Institute of Physics: New York, 1984; AIP Conf. Proc. No. 107. (12) Klymko, P. W.; Kopelman, R. J . Phys. Chem. 1982, 86, 3686. (13) Mandelbrot, B. B. ‘The Fractal Geometry of Nature”; Freeman: San Francisco, 1983. (14) Schlesinger, M. F., Mandelbrot, B. B., Rubin, R. J., Eds. “Proceedings of a Symposium on Fractals in the Physical Sciences”; J . Stat. Phys. 1984, 36, Nos. 516. (15) Alexander, S.; Orbach, R. J . Phys. (Paris) 1982, 43, L625. (16) Rammal, R.; Toulouse, G. J . Phys. (Paris) 1983, 44, L13. (17) Kopelman, R. In “Radiationless Processes in Molecules and Condensed Phases”; Fong, F. K., Ed.; Springer-Verlag: Berlin, 1976; p 297; Top. Appl. Phys. Vol. 15. (18) Keramiotis, A,; Argyrakis, P.; Kopelman, R. Phys. Reu. E 1985, 31, 4617. (19) Witten, T.; Sander, L. M. Phys. Reu. E 1983, 27, 5686. (20) Kopelman, R.; In “Fractal Aspects of Materials: Metal and Catalyst Surfaces, Powders and Aggregates”; Mandelbrot, B. B., Passoja, D. E., Eds.; Materials and Research Society: Pittsburgh, PA, 1984; p 21. (21) Kopelman, R.; Klymko, P. W.; Newhouse, J. S.; Anacker, L. W. Phys. Reu. B 1984, 29, 3747. (22) Klymko, P. W.; Kopelman, R. J . Phys. Chem. 1983, 87, 4565. (23) Anacker, L. W.; Kopelman, R.; Newhouse, J. S . J . Stat. Phys. 1984, 36, 591. (24) Anacker, L. W.; Kopelman, R. J . Chem. Phys. 1984, 81, 6402. (25) Newhouse, J. S.; Argyrakis, P.; Kopelman, R. Chem. Phys. Left. 1984, 107, 48. (26) Kang, K.; Redner, S . Phys. Reu. L e f t . 1984, 52, 955. (27) Meakin, P.; Stanley, H. E. J . Phys. A 1984, 17, L173. (28) Argyrakis, P.; Anacker, L. W.; Kopelman, R. J . Stat. Phys. 1984,36, 579. (1) York, (2) (3) (4)
0022-365418512089-4758$01.50/0 0 1985 American Chemical Society
Kinetics on Fractal and Euclidean Lattices
The Journal of Physical Chemistry, Vol. 89, No. 22, 1985 4159
-
For the elementary binary reaction, A + A products, a scaling argumentZ2for compact exploration replaced the classical kinetic rate laws for transient systems (-dp/dt = Kp2) with
and to relate the parameters in them to microscopic models. In the usual phenomenological approach to chemical reactions one considers a rate equation for the spatially averaged density p(t) of reactants dp/dt = F(p,t)
(2.1)
Equation 1.2 has also been tested extensively, and both simulat i o n and ~ ~e x~p ~e r ~i m~e n t ~are ~ ~well ~ ~described ~ by it once initial conditions are no longer significant. Here the integrated form gives us:22
and interprets the function F i n terms of the motion and collisions of individual molecules. If one knows beforehand that the reaction occurs through binary collisions and is diffusion controlled, one expects to find
tf, f
dp/dt = -K(t)p2
PA-’
For the transient A
+B
- 0:
(1.3)
0 reaction a parallel expression26
is obtained in the asymptotic limit, given that PA(t) = pB(t) for in the exponent of eq 1.4 is due to local all t. The factor of density fluctuations which cause the formation of clusters of A and clusters of B particles.26 Trapping is a special case of A + B reactions where A B B; the trapping problem is now well understood even out to very long timesZ9(if one is allowed to adjust the individual diffusion coefficients DA, DBso that D B= 0; this makes a difference in the extreme long time limit). The above led us to an investigation of the steady state, if one exists, on fractals and in Euclidean spaces. To the A A reaction described in eq 1.2, we add a constant source of particles at rate R. The rate law for this system was conjecturedZ4to be
+
-
+
dp/dt = R - kpX X = 1 +f‘
(1.5)
where X i s the effective reaction order. This has now been investigated by using large scale Monte Carlo simulations on the linear lattice,24the simple cubic lattice,24the planar Sierpinski gasket,24 and the two-dimensional percolation cluster.30 All simulation results have demonstrated the existence of a steady state, pss, which follows the relation
R
0:
pssx
X = 1 +f-’
(1.6)
for small R and consequently small p,,. To be sure, when R is set equal to zero, eq 1.5 and 1.6 are inconsistent with eq 1.2 and 1.4. This indicates a singular perturbation problem; we discuss this further in the next section. Here we present the results for the approach to the steady-state population and for deviations from eq 1.6 for large steady-state concentrations. We show that the approach to steady state can be described by
and that the high density steady state can be described by R
Psf
(1.8)
where X‘> X in Euclidean dimensions while X’< X and X’2 for the fractal Sierpinski gasket. Section 2 gives the development of the kinetic analysis while the simulation methodology and results are described in sections 3 and 4. Section 5 examines deviations from the theory and ascribes them to high density effects.
2. Kinetic Analysis of Steady State and Relaxation Reduced dynamical equations play an important role in the theory of interacting many-body systems. One hopes that such a contracted description, involving a few macroscopic variables, may be applied to a wide class of systems, so that by specifying the values of a few parameters and solving the reduced equations one may obtain the dynamics of the system of interest. The form of the reduced equations is often suggested by intuition; theoretical studies attempt to establish the range of validity of the equations (29) Klafter, J.; Blumen, A.; Zunofen, G. J . Star. Phys. 1984, 36, 561. (30) Newhouse, J. S.; Kopelman, R. Phys. Reu. E 1985, 31, 1677.
in which K ( t ) is related to the diffusive motion of the reactant molecules. Accepting the form of eq 2.2, one can calculate K ( t ) by considering the mutual diffusion of a pair of molecules; a mean-field argument relates the rate at which the pairs react to the overall reaction rate. This program was carried out by S m o l u c h o ~ s k ifor ~ ~a ~model ~ ~ ’system ~ ~ ~ ~of~diffusing ~ spheres, in a three-dimensional medium, which react on contact. Smoluchowski’s result is K(t) = 87R0D[l + R , / ( ~ T D ~ ) ’ / ~ ]
(2.3)
in which Ro is the radius of the sphere and D the diffusion coefficient. For t >> Ro2/D, K ( t ) becomes independent of time. If the same analysis is applied to reactions in one- or two-dimensional systems, one finds that K ( t ) goes to zero at long times so that a rate constant does not exist; for example, in one dimension one finds K(t) (D/At)ll2. This sort of behavior33has been called “nonclassical”. Other model systems, such as reacting random walkers on ordered lattices, show similar behavior in appropriate limit~.~~J~ Recently, “nonclassical” rate coefficients have attracted attention through the study of diffusion-controlled reactions on fra~tals.~O-~O These objects are characterized by effective dimensions which need not be integers; the value of this dimension depends to some extent on the problem being studied. For example, the solution to the diffusion equation on a fractal resembles the expression one obtains for a uniform medium, by replacing, in the fundamental solution (autocorrelation function), the dimension of the medium by the “spectral dimen~ion”,’~,’~ d,, of the fractal. Considerations of this sort lead to the postulated expression
-
K(t)
- K0tf-l,
t
-
(2.4)
for the rate coefficient, in which KOis a constant (involving a dimensionally appropriate replacement of DRo) and
2
(2.5)
For d, = 2 there are logarithmic corrections to eq 2.4; in what follows we restrict ourselves to d, # 2. An elegant intuitive argument for eq 2.4 involving the ”exploration volume” swept out by a reactant molecule has been given by de G e n n e ~as , ~well as a more explicit treatment for the trapping p r ~ b l e m . ~ ~ , ~ ~ Equations 2.2 and 2.4 have been shown to be quite accurate by computer simulation^^^^^ as well as by exact calculations in one d i m e n ~ i o n . ~These ~ ’ ~ ~verifications have all involved transient reactions, in which eq 2.2 is the only mechanism for changing the density. If eq 2.2 is correct, then one should be able to obtain p(t) in the presence of an external source by adding a source term: dp/dt = R
- K0tf-I~’
(2.6)
The source term R is taken to be uniform in space and independent (31) Smoluchowski, M. V. 2.Phys. Chem. 1917, 92, 129. (32) Chandrasekhar, S . Reu. Mod. Phys. 1943, 15, 1 (reprinted in: “Noise and Stochastic Processes”; Wax, N., Ed.; Dover: New York, 1954). (33) Suna, A. Phys. Reu B 1970, 1 , 1716. (34) Kenkre, V . Phys. Rev. B 1980, 22, 2089. (35) de Gennes, P. G. C. R. Acad. Sci. (Paris) 1983, 296, 881. (36) Evesque, P.; Duran, J. J . Chem. Phys. 1984, 80, 3016.
4760
The Journal of Physical Chemistry, Vol. 89, No. 22, 1985
of time. Forf # 1 one quickly sees that this implies an unbounded increase in the density p(t); explicitly, one finds p(t) a t(l-n/* as t m . (This result may be obtained by transforming eq 2.6 to a linear second-order equation and applying standard asymptotic method^.^' In the special casef= (d, = l), an exact solution involving Airy functions can be obtained.) Since eq 2.2 and eq 2.4 are valid in the limit of low densities, where multiple particle correlations may be neglected, eq 2.6 must eventually break down.33 One may say that the solution to eq 2.6 is unstable with respect to many-body effects. In this paper we show by direct simulation that eq 2.6 does not describe the reaction dynamics under steady-source conditions, and propose an alternative. Our replacement for eq 2.6 is indicated as follows: we desire pss a differential equation having a steady-state solution, p(t) as t m, since simulations do show such behavior. The equation should then be autonomous dp/dt = R - E ( p ) (2.7)
-
-
-+
where E has no explicit time dependence. Now, eq 2.7 cannot be correct for all R since for R = 0 it contradicts eq 2.2. However, we can choose E(p) so that when R = 0 the solutions to eq 2.2 and 2.7 agrees as t OD. In this way, a smooth crossover between eq 2.2 and 2.7 can be constructed, although we do not pursue that here. The solution to eq 2.2 and 2.4 is
-
If we then set dp/dt = R - Klpl+l/f
(2.9)
so that when R = 0 (2.10) we have the desired agreement in the limit when p ( t )
