5462
J. Phys. Chem. 1989, 93, 5462-5467
Asymptotic Analysis of Dlffuston-Influenced Kinetics with a Potentialt Nicholas J. B. Green and Simon M. Pimblott* Physical Chemistry Laboratory, South Parks Road, Oxford, OX1 3QZ. U.K., and Radiation Laboratory, University of Notre Dame, Notre Dame, Indiana 46556 (Received: August 29, 1988; In Final Form: January 10, 1989)
A transformation is introduced and used to formulate a series of approximations for the time-dependent diffusion-controlled geminate reaction probability. Each approximation has the properties of a probability distribution function. The implications of the approximations are analyzed by use of established relationships to give corresponding approximations for the geminate reaction probability with a radiation boundary condition and for the time-dependent rate coefficient with a Boltzmann initial concentration profile. The first-order approximation, which holds when the interparticle forces are weak, is applied to the recombination of ions in an aqueous solution of a strong electrolyte; its accuracy is assessed by comparison with numerical results.
1. Introduction In a recent series of we have investigated the kinetics of diffusion-controlled reactions between ions in solvents of high permittivity. In ref 1 we presented analytic approximations to the time-dependent diffusion-controlled geminate ion recombination probability; in ref 2 we investigated the consequences of the approximation and derived solutions of comparable accuracy for geminate recombination with a radiation boundary condition and for bulk time-dependent rate constants; in ref 3 we showed how the approximations could be used in a stochastic model to describe reactions of ions in clusters, such as those found at short times in radiation tracks. Although the unscreened Coulomb potential used in ref 1-3 is applicable to some problems, for other situations the interparticle force is modified 'by the presence of other chemically inert ions in solution. The theory of the modified potential energy of interaction, derived by Debye and Hiicke14for low ionic strengths, suggested that the potential should be modified by the inclusion of an exponential screening factor. The screened Coulomb potential is known to be asymptotically correct in the limit of low ionic ~ t r e n g t h . ~Simple extensions for higher concentrations of ions do not alter the functional form. Explicit results presented in this paper for the screened Coulomb potential are, of course, only applicable under conditions where the ionic strength is sufficiently small for the potential to have the correct form. Section 2 is devoted to a description of the geminate recombination problem with a generalized potential. We make a transformation of the geminate reaction probability based on the infinite divisibility of a first-passage time distribution. This well-established mathematical property is discussed in section 2. The transformation leads to an expansion, which is known to be convergent for a screened potential,6 and may be truncated after the first term to give an approximate reaction probability, analogous to that for the Coulomb potential presented in ref 1. In section 3 the basic approximation is used to obtain corresponding approximations for the geminate reaction probability with a radiation boundary condition and for the rate constants of bulk ionic reactions, following the methodology of ref 2. Finally, in section 4 we present explicit results for the screened Coulomb potential, together with comparisons with numerical calculations and simulations which indicate the range over which the approximation will be numerically useful. 2. Geminate Recombination The mathematical analysis of diffusion-controlled reactions is usually performed using the Debye-Smoluchowski equation. It is, by now, well established that the geminate reaction probability, expressed as a function of the initial relative configuration of the pair and of time, obeys a Kolmogorov backward equation, which is mathematically adjoint to the Debye-Smoluchowski equa'Document No. NDRL-3126 from the Notre Dame Radiation Laboratory.
0022-365418912093-5462$0 1.50/0
tion.'-1° In this paper we will assume that both the diffusion and the potential energy of interaction of the pair are spherically symmetrical, so that the reaction probability W(x,a,t) is a function only of the initial separation of the pair x, the time t, and the reaction distance a. The backward equation for Wcan be written
where U is measured in units of kT, D is the relative diffusion coefficient of the pair, and the following boundary conditions are employed: W(x,a,O) = 0
(x
> a)
W(a,a,t) = 1
(I
> 0)
lim W(x,a,t) = 0 x--
where the variable a represents the reaction distance of the pair. The second of these conditions is the Smoluchowski boundary condition, in which reaction is assumed to occur instantaneously on encounter." Following the procedure described in ref 7 the backward equation is reduced to its canonical form:
aw - = - - -1 a aw at
2 aM
as
(3)
where S is a transformed distance scale and M is the speed measure.12 s(x) = -x-zeWx)
m(x) = 1/[2Ds(x)] S(x) = l x s ( x ) dx M ( X ) = l x m ( x ) dx
(4)
( 1 ) Clifford, P.; Green, N. J. B.; Pilling, M. J. J. Phys. Chem. 1984, 88, 4171. (2) Green, N. J. B. Chem. Phys. Lett. 1984, 107, 485. (3) Clifford, P.; Green, N. J. B.; Pilling, M. J.; Pimblott, S.M. J . Phys. Chem. 1987, 91, 4417. (4) Debye, P.; Huckel, E. Phys. 2.1923, 24, 185. (5) Stell, G.; Patey, G. N.; Hmye, J. S.A h . Chem. Phys. 1981, 48, 183. (6) Sibani, P.; Pedersen, J. B. Phys. Reo. Lett. 1983, 51, 148. (7) Karlin, S.;Taylor, H. M. A Second Course in Stochastic Processes; Academic: London, 198 I . (8) Kolmogorov, A. N . Math. Ann. 1931, 104, 415. (9) Sano, H.; Tachiya, M. J . Chem. Phys. 1979, 71, 1276. (IO) Clifford, P.; Green, N. J. B.; Pilling, M. J. Chem. Phys. Lett. 1982, 91, 101. ( 1 1 ) Smoluchowski, M . yon Z. Phys. Chem. 1917, 92, 129.
0 1989 American Chemical Society
The Journal of Physical Chemistry, Vol. 93, No. 14, 1989
Asymptotic Analysis of Diffusion-Influenced Kinetics The limits on these integrals are not usually specified, but the definitions of s and S adopted here are particularly convenient m and U - 0 then S ( x ) l/x because they ensure that as x asymptotically. The backward equation can be solved exactly for the probability that reaction will ultimately occur, W=(x,a),and the solution is given in terms of the scale measure S defined above:
-
+
W A a )= S(x)/S(a)
(5)
Since S ( x ) asymptotically approaches l / x at large x , we can treat 1 / S ( x ) as an effective distance xOffand write instead
wco(x,a)=
5463
s*(x) = -S(X)/S(X)Z
m*(x) = - S ( X ) ~ / ~ D S ( X ) S * ( x ) = s x s * ( y ) dy = l / S ( x ) = x,ff M * ( x ) = Jxm*Cy) dy
- -
-
We note that as x m , S* x and M* x/2D. Taking the Laplace transform of (9) and making the transformation of (8), we find
aeff/xeff
which points out the relationship with the well-known solution for radical recombination. The time-dependent reaction probability, however, is much more difficult to analyze. Very few exact solutions are available at present; indeed in the time domain the only exact solution is for the case U = 0. The Laplace transform of W is available for the Coulomb potential,I3 but its numerical computation is a considerable undertaking, and the Laplace transform must be inverted numerically. In order to obtain an approximation for Wthat has the properties of a probability distribution function, we use the infinite divisibility of the distribution of a first-passage time for a one-dimensional diffusion process. This property can be understood in a simple physical manner. Diffusion processes have sample paths that are continuous (by definition); hence, in order to arrive at the state a a diffusion must pass through all intermediate states. The first-passage time for x to a can be decomposed into the sum of the first-passage time to some intermediate state b and the subsequent first passage from b to a: if w represents the first-passage time density function, this can be written w(x,a,t) = J‘w(x,b.t’) w(b,a,t-t’) dt’
(a
< b < x)
(6)
Formally this decomposition arises because the diffusion obeys the strong Markov property and the intermediate first-passage time is a Markov time.’ Each first passage time can be decomposed into the sum of arbitrarily many component first-passage times, hence the distribution is infinitely di~isib1e.l~ Feller (ref 14, section XIII.7, Theorem 1 ) has proved that the Laplace transform of any infinitely divisible distribution can be written in the form (7)
which is a Riccati equation for 4. The equivalent equation for the special case of the Coulomb potential is derived in ref 1 . Equation 10 is strongly related to the Riccati equation for the Laplace transform of the bulk rate kernel derived by Sibani and Pedersen.6 The origin of this relationship will be discussed in section 3. We now attempt to find an expansion for C#J in powers of p’12 following the procedure of ref 1 and 6: 4(X,P) = 41P1/z + 42P
This expansion leads to the following system of differential equations for the coefficients &: &’ - 4‘ d In s*/dx = 0
4; - d2 d In s*/dx =
&’ - & d In s*/dx =
at In consequence the scale and speed parameters of the conditioned process are given by (12) ItS, K.; McKean, H. P. Diffusion Processes and their Sample Paths, 2nd ed.; Springer: Berlin, 1974. (13) Hong, K. M.; Noolandi, J. J . Chem. Phys. 1978,68, 5163. (14) Feller, W. An Introduction to Probability Theory and its Applications; Wiley: New York, 1971; Vol. 2.
n
X~~J~C$~~
The system of equations for 4ncan be solved recursively: the first two solutions are &(x) = s*(x)/D112
Some care is necessary in evaluating 42 since both terms in the integrand tend to 1 at large x. For reasonable potentials, however, the integral is defined. The first-order approximation to W is therefore
erfc Xdf
aw
-p / D
i- 1
where p is the Laplace transform variable conjugate to t. Furthermore, cP(x,a,O) = 0, and for the problem of interest the inner boundary condition prescribes that +(a,a,p) = 0; hence we write
in the same way as in ref 1, where 4 = a@/ax. Equation 8 can be interpreted in the following way. The reaction probability is equal to the ultimate reaction probability W , multiplied by a function containing the time-dependence, w * ( x , a , t ) . W represents the reaction probability conditioned on ultimate reaction and obeys its own backward equation whose parameters are related to those of (1).Io It is shown in ref 10 that the backward equation for W is
+ ...
(-)
Xeff
- acff
(12)
( 4 D t ) I/z
This approximation is intuitively very satisfactory as it is obviously related to the solution in the absence of any potential by a simple change of scale. The second order approximation corresponds to shifting the time origin by a constant QZ.The approximations will be valid for a given a if t is sufficiently large, or for all interesting t if a (and therefore x ) is sufficiently large. Numerical evaluation of the approximation for the Coulomb potential has been presented in ref 1 , and for the screened potential results are shown in section 4. The advantage of this method is that successive approximations all represent proper probability distributions with the required infinite divisibility. Most previous analyses have neither property; however, for applications such as the IRT simulation technique of ref 3 it is essential to have an accurate approximation to the whole distribution function. Pedersen and Sibani16 have suggested the same approximation, among others, by using different methods, and have recommended the form derived here as the “most convenient”. We have shown in this (15) Pedersen, J. B. J . Chem. Phys. 1980, 72, 3904. (16) Pedersen, J. B.; Sibani, P. J . Chem. Phys. 1981, 75, 5368.
5464
The Journal of Physical Chemistry, Vol. 93, No. 14, 1989
section, however, that this approximation arises naturally from the fundamental mathematical properties of the diffusion approximation and is therefore preferable to the alternatives.
3. Consequences of the Approximation The mathematics relating the diffusion-controlledrecombination probability to the corresponding probability for partially diffusion-controlled reactions and to various time-dependent rate coefficients is well-knonw, and its application to the transformed reaction probability presented in section 2 is straightforward. Pedersen and Sibani16 have already presented the consequences of the first-order approximation in most of the cases considered below. However, we also present the same procedures applied to the exact transformed solution in order to trace the functions Q and through the various relationships introduced. This enables us to consider the effects of various alternative approximations in section 4. 3.1. Radiation Boundary Condition. The mathematical relationship between the reaction probability for the absorbing boundary and the corresponding probability for a radiation boundary condition (elastic boundary) has been derived by Pedersen.15 If the radiation boundary condition on the adjoint equation is formulated9
+
U
= - [ 1 - W(x,a,t)] then the Laplace transform of the reaction probability subject to this reactive boundary condition can be written in terms of the functions Q and defined in section 2 by
+
mrS/(x,a,P;u) =
(14) The term in square brackets is the probability of ultimate reaction, so that the remainder of the right-hand side represents the Laplace transform of the conditioned reaction probability. Equation 14 is exact, but we may not invert the Laplace transform unless we have an explicit expression for 4. Substituting in the first-order approximation to from section 2, we find for the corresponding radiation boundary condition in terms of xeff = l/S(x):
+
Wr&(x,a,t;u) =
Green and Pimblott approximation is necessary before a rate constant can be derived. Although the geminate survival probability may be an exact solution within the framework of the diffusion approximation, the formulation of bulk rate constants depends on the ‘independent pairs” approximation. In the “concentration gradients” approach to the rate constant the concentration is assumed to obey a diffusion equation in a coordinate system in which the central particle is stationary. This assumption, identified in the review of N ~ y e s , ~ is equivalent to assuming that all the reactive distances in the ensemble evolve independently of each other. The same approximation is made, either explicitly or implicitly, in the derivations of this relationship, but is only strictly correct if the central particle is truly stationary and there are no interparticle forces. In particular, in the case where the traps are static and no forces are present, the approximation cannot be accurate a t long times because in three dimensions the survival probability of the diffusing particle behaves asymptotically as exp(-ht0.6).21 The only test of the approximation when ionic forces are present has been published by Clifford et al.,3 who report that it is accurate for Coulombic forces between small clusters of ions in aqueous solution. Although the time-dependent rate coefficient can be derived for any initial extended reactant concentration profile, in this paper we will only derive explicit solutions for an initial Boltzmann distribution, where the concentration of reactant about the central ion is equal to coe-u(r),and co is the bulk concentration. We make this restriction here because the primary application considered is to a shielded Coulomb potential and, in order for this potential function to be applicable, the appropriate initial distribution of ions must be the equilibrium distribution. In a case where an ion is suddenly formed from a neutral particle and subsequently reacts by diffusion control, the transient part of the decay occurs on the same time scale as the formation of the ionic atmosphere and so the potential of mean force does not take on its equilibrium form until the transient is virtually complete. Consequently we will only apply the approximations to the initial Boltzmann equilibrium radial distribution. For an initial Boltzmann distribution the rate constant can be expressed in terms of the reaction probability:2+6 k ( t ) = - 4 ~ D a ~ e - ~dW/dxlxla (“)
(16)
Thus we can express the Laplace transform of the rate constant for a Smoluchowski boundary condition using the function introduced in section 2.
+
where the dimensionless groups a, /I, and 6 are given by
p = (1
+ l/6)(Df)’12/Ueff
6 = -Ds(a)/uS(a) = Da,rpu(a)/va2 If (1 2) is an adequate approximation to the first-passage time distribution, then (15) will be correspondingly accurate for a radiation boundary condition. (See, for example, the comparisons in ref 2.) 3.2. Boltzmann Rate Constant. The relationship between geminate survival probability and the rate constant for homogeneous reaction is well established. The earliest derication appears to be due to Steinberg and Katchalski,” who considered only the case of free diffusion. Subsequent generalizations and applications of the method have been reported by Tachiya18 and Green.2 A similar relationship, between geminate survival probability and concentration gradients, has been derived by Gosele.19 The relationships are often presented as but an additional (17) Steinberg, I. 2.;Katchalski, E. J . Chem. Phys. 1968, 48, 2404. (18) Tachiya, M. Radial. Phys. Chem. 1983, 21, 167. (19) GBsele, U.Chem. Phys. Lefr. 1980, 69, 332.
The first term is the Laplace transform of the steady-state rate constant, and the second term represents the transient part of the rate constant. We can now understand the origin of the similarity between the Riccati equation for derived in section 2 and the Riccati equation for the Laplace transform of the rate kernel derived by Sibani and Pedersem6 The two approaches are seen to be cquivalent, in that there is a simple relationship between @ and k(p) so that a_ny expansion for in powers of p1I2implies an expansion for k and vice versa. Consequently we can apply the first-order approximation to and obtain, after inversion, the same result as Sibani and Pedersen:6
+
+
+
F)
k ( t ) = 47rDuefl( 1 + *Dt)’I2
We note, however, that this expression does not extend to the case of the initially uniform concentration profile. For example, for (20) Noyes, R. M. Prog. React. Kine?. 1961, I , 129. (21) Donsker, M. D.; Varadhan, S. R. S. Comm. Pure Appl. Math. 1975, 28, 5 2 5 .
.
.
*
I . , . , .
.
- ,.
4
..
-7.
Asymptotic Analysis or uirrusion-~nr~uencea ninetics 1
the case of the Coulomb potential the time-dependent part of the uniform rate constant must be independent of the sign of the interparticle force,2 which is clearly not the case with (18). The relationship between the geminate survival probability and the time-dependent rate coefficient can also be applied to the case of the radiation boundary condition, (13). The independent pairs rate constant can be expressed in terms of 4,
-
kr&)
:;;['
= -- - - u
The Journal of Physical Chemistry, Vol. 93, No. 14, 1989 5465 I
2.5
u
+ D (4(a) - s(a)/S(a))
Employing the first-order approximation for 4, we obtain 0.0 0.0
where /3 and 6 are the same as in (15). The advantage of employing this form for the rate constant in the presence of the radiation boundary condition is that although it is the result of a first-order approximation for 4 it is actually accurate over a much wider range of times than a simple asymptotic expansion of k. The asymptotic expansion falls off toward the steady-state rate constant as t-II2, but the use of such an expansion at short times leads to an unphysical divergence.z2 In contrast, at small times (20) reduces to the correct limit, 4aa2e-u(a)u.Comparisons of (20) with numerical solutions for the Coulomb potential were reported in ref 2. The expression was shown to be a great improvement on that of Rice et
I .O
2 .o
x/nm Figure 1. Comparison of effective distance scales: (- - -) pure Coulomb potential, ro = -0.7 nm; (-) screened Coulomb potential, r, = -0.7 nm, rD = 3.0 nm, a l = a2 = 0.3 nm; (---) zero potential.
In this work we use an alternative analytic method of evaluating S ( x ) . The required integral is c
which can be expressed as a sum of exponential integrals
4. Screened Coulomb Potential The use of a screened potential in the diffusion equation describing concentration gradients for ionic reactions in solution was first suggested by D e b ~ e .The ~ ~ screened potential he employed arises from the Debye-Huckel theory4 and can be written
where r, is the Onsager length and rD is the Debye screening length: rc = z , z 2 e 2 / 4 a ~ k T rD-2= (e2/ck7')xciz?
l
ea,/% + = 2l l 1 +ea2/z'D a l / 2 r D 1 + a2/2rD
zl, z2, a,, and a2 are the charges and radii of the two ions. For the limiting law of low ionic strengths the constant y equals 1.
Logan has suggested weighting the two terms in the expression for y in the ratio of the diffusion coefficients of the ions.24 The choice of y does not affect any of the analysis that follows. In our numerical calculations we have used the Debye-Hiickel prescription and in addition have assumed that al = a2 = a / 2 for convenience. 4.1. Steady-State Solutions: The Scale Function. The steady-state solutions, both for the ion pair escape probability and for the rate constants, can be expressed in terms of the scale function, S(x), which is related to the potential by ( 4 ) . Unfortunately, if the screened Coulomb potential is used in the expression for S(x) in (4), the resulting integral cannot be found in closed form. This led Debye to consider only the asymptotic form of the steady-state rate coefficient for large values of rD. Montr01l~~ also derived an approximation, but this method was equivalent to ignoring the exponential damping of the potential and decreasing the dielectric constant of the solvent. It therefore has the effect of strengthening the potential rather than weakening it. In later work Loganz4 has used numerical integration to evaluate the steady-state rate coefficient, and this method has been used in more recent ~ o r k . ~ ~ , ~ ' (22) Rice, S. A,; Butler, P. R.; Pilling, M. J.; Baird, J. K. J . Chem. Phys. 1979, 70,4001. (23) Debye, P. Trans. Electrochem. SOC.1942, 82, 265. (24) Logan, S . R. Trans. Faraday SOC.1966, 62, 3416. (25) Montroll, E . W. J . Chem. Phys. 1946, 14, 202.
where E&) = J;x-me-vx dx.28 The exponential integrals can be computed efficiently according to the algorithm of G a ~ t c h i ~ ~ and the series converges rapidly. The corresponding value of xeff can also be calculated and is shown in Figure 1. This figure makes a comparison with the effective distance scales for the pure Coulomb potential (Le., neglecting screening) and for r, = 0 (Le., neglecting the electrostatic interactions completely). The screened scale function is intermediate between the Coulomb scale and the radical scale and decays toward the latter, as expected. The solution can now be used in the limiting values of (12) and (1 5 ) for the geminate reaction probability and in the limiting values of (18) and (20) for the rate coefficient. This is the first time an exact analytic solution has been presented for these quantities. 4.2. Geminate Recombination Probability. Only one analysis of the time-dependent problem is known to the authors. This approach, due to Montroll,zs uses an inappropriate perturbation technique,30 with the consequences described in section 4.1. However, Montroll also introduced an approximation for the function 4, which he assumed to take the same value as it does for radical recombination, i.e. 4 = l/D1/z. If Montroll's approximate scale function is denoted SM,then his approximation is equivalent to
This is compared with the approximation derived in this work, (1 2), and with a numerical solution of the backward equation for W(x,a,t)in Figure 2. Although for the screening length employed here (3.0 nm) the approximation of Montroll is satisfactory, as the screening becomes more substantial, the Montroll result becomes progressively worse. This arises mainly from his inappropriate approximation for the escape probability, rather than from the approximation for the time scale, as can be seen from (26) Byakov, V. M.; Petuchov, R. Radiochem. Radioanal. Lett. 1982,51, 143. (27) Logan, S. R. Radiochem. Radioanal. Lett. 1982, 54, 259. (28) Abramowitz, M.; Stegun, I. A. Handbook of Mathematical Functions; Dover: New York, 1972. (29) Gautschi, W. C.A.C.M., 1973, 16, 471. (30) Rice, S. A. Diffusion Limited Reactions; Comprehensive Chemical Kinetics 25; Bamford, Tipper, Compton, Eds.; Elsevier: Amsterdam, 1985.
5466 The Journal of Physical Chemistry, Vol. 93, No. 14, 1989
Green and Pimblott
0.41
'""I
400
T
*E
t
O.I 0.oL
0.0
,-
A
I
.o
kA
/
1
I
2.0
3.0
I
4.0
Figure 2. Geminate recombination probability W(x,u,t)of an ion pair in a solution of a strong electrolyte: rc = -0.7 nm, r~ = 3.0 nm, a1 = a2 = a / 2 = 0.3 nm, x = 2.0 nm, D = 1.0 X IO-* m2 s-l. (+++) numerical solution of backward equation; (-) approximation of this work, from (12); (---) approximation of M ~ n t r o l lfrom , ~ ~ (24).
I
0
I
\
\ -9.0
-8.0
-7.0
log ( x / m ) Figure 4. Infinitesimal mean of the conditioned diffusion equation: r, = -0.7 nm, rD = 3.0 nm, a, = u2 = 0.3 nm, D = 1.0 X IO-* m2 s-'. (-) pure Coulomb potential; (+++) screened Coulomb potential. forms for @ and 4 into (14),we obtain for the radiation boundary condition aeff [erfc ( a ) + exp(2aP) erfc ( a + p)] xed1 + 6 ) (26) where a is given above: u(1 + 6) /3 = -(Dt)'i2 Ds+(a)
0.4
1
0.0 0.0
I
W+(x,a,t,u)=
3 0.2
loot
-10.0
log ( t i m e / ps 1
0.8
I\
30011
I .o
6 = Da,fpu(4)/ua2 2 .O
3.0
4.0
l o g ( t i m e / ps 1
Figure 3. Geminate conditioned reaction probability, w*(x,a,t): rc = -0.7nm,rD=3.0nm,y= 1.011,a=0.6nm,x=2.0nm,D=l.OX m2 s-I. (+++) numerical solution of conditioned backward equation, screened potential; (- - ) approximation for pure Coulomb potential;' (-) approximation for screened potential, from (12); (---) approximation of M ~ n t r o l l from , ~ ~ (24).
Figure 3, where the approximations for the conditioned reaction probability are shown together with the corresponding function for the pure Coulomb potential. With most sets of parameters we have analyzed we find that there is little difference between the three time scales. The screened potential has the best asymptotics and is mathematically most elegant, but the zero potential solution (Montroll) is as good for most of the curye and the pure Coulomb scale is, if anything, even better. The similarity of the screened and pure Coulomb scales can be rationalized by comparing the infinitesimal means in the conditioned backward equation, (9). This is shown in Figure 4 for a value of rD = 3.0 nm, which is appropriate for an ionic strength of 0.01 mol dm-3. The two curves are almost indistinguishable and the drifts are small in magnitude, with the result that the scale measures for the conditioned processes are also very similar. The conclusion is that although we must use the correct scale for the ultimate reaction probability, the conditioned process, which determines the time scale of the reaction, is much less sensitive to the approximation used. The best agreement is found for the conditioned pure Coulomb scale, except for very short distances, where the zero potential scale is better, but any of the three approximations is good enough to be applied under most realistic conditions. 4.3. Radiation Boundary Condition. If a different conditioned scale is used for the time dependence in the geminate reaction probability, rather than xeffderived in section 2, then the modification carries through to the corresponding solution with a radiation boundary condition. Thus, if we use the scale S+(x), with density,'s instead of the scale S*(x), the geminate recombination probability with a Smoluchowski boundary condition becomes v ( x , a , t )= (aeff/xeff) erfc ( a )
(25)
where a = (S+(x) - S ' ( ~ ) ) / ( 4 D t ) l /Substituting ~. the appropriate
Comparisons for the pure Coulomb potential showed that this approximation was as good as that for the Smoluchowski boundary conditiom2 The similarity of the infinitesimal parameters for the two conditioned processes shown in Figure 4 guarantees that this will also be the case here. 4.4. Boltzmann Rate Constants. The same substitutions for @ and 4 can be made in (17) and (19), and the Laplace transforms can be inverted to give approximations for the rate constants. These approximations are
for the Smoluchowski boundary condition and
for the radiation boundary condition, where p and 6 are the same as in (26). We note that these approximations will have the same regions of validity as the basic geminate recombination approximation. We also note that if the conditioned scale S* is used for S, then these equations reduce to (18) and (20). The limiting values at long and short times are both independent of the choice of s+. 5. Conclusions In this paper we have used well-established mathematical properties of the first-passage time distribution function to develop a transformation of the geminate recombination probability. The transformation reduces the backward equation for the reaction probability to a Riccati equation. Recursive solution of the Riccati equation gives rise to a series of approximate solutions for the geminate reaction probability as a function of time. Each solution is a probability distribution function with the property of infinite divisibility and the correct asymptotics at long times. The first two orders of approximation can be obtained, and we have examined the implications of the first-order approximation for the case of geminate recombination with a radiation boundary condition and for the time-dependent rate coefficients with an initial Boltzmann concentration profile. The approximations will only be valid if the forces between the particles are weak, and we have therefore concentrated our explicit work on the screened Coulomb
J . Phys. Chem. 1989, 93, 5461-5414 potential in aqueous solution. We have presented an exact analytic solution for the ultimate escape probability and the steady-state rate constants by finding the scale function of the process. We have also examined the geminate reaction probability as a function of time, and to assess the validity of the approximate treatment we have compared our results with numerical solutions of the backward equation. We find that it is necessary to employ the correct scale function in the asymptotic reaction probability but that the time dependence of the reaction probability conditioned on ultimate reaction is rather insensitive to the conditioned scale chosen. We have shown how a modification of this conditioned scale carries through to the corresponding geminate reaction
5467
probability with a radiation boundary conditon and to the Boltzmann rate constants. All these approximations hold under conditions where the extended Debye-Huckel limiting law holds in aqueous solution, Le., where the screened Coulomb potential is applicable. Furthermore, the approximate rate constant for the radiation boundary condition holds exactly not only in the limit of long times but also in the limit of short times. Acknowledgment. We acknowledge helpful discussions with Drs. M. J. Pilling, P. Clifford, and A. Mozumder. This work was supported by the Office of Basic Energy Sciences, of the U.S. Department of Energy, and A.E.R.E., Harwell, U.K.
Amplitudes and Phases of Small-Amplitude Belousov-Zhabotinskii Oscillations Derived from Quenching Experiments P. Graae Sarensen* and F. Hynne* Department of Chemistry, H . C. 0rsted Institute, University of Copenhagen, Universitetsparken 5, DK-2100 Copenhagen, Denmark (Received: October 6, 1988; In Final Form: February 2, 1989)
We report a series of accurate experimental quenchings of harmonic limit cycle oscillations near a supercritical Hopf bifurcation of the Belousov-Zhabotinskii reaction system. In a quenching experiment, the concentrations of some chemical species are suddenly changed by addition of substances or by dilution so that oscillations temporarily stop. Quenching by addition was successful with Ce4+,HBr02, HBrO, Brz, and Br- but not with any of several organic compounds tested. The results can be used for chemically relevant quantitative tests of models of the reaction. They suggest that at least four species are necessary for a description of the main features of the oscillatory system. The amplitudes and phases of the oscillatory parts of [Ce4+], [HBr02], and [Br-] have been calculated from three quenching experiments on the assumption that all other species oscillate with sufficiently small amplitudes but without any assumption about the dimension of the phase space. In a strictly three-dimensional Oregonator-type phase space, we have calculated the average concentrations from one additional quenching by dilution, but the result is physically unacceptable in the present case. Indeed, the phase of that quenching shows that a species not tested in the present study must be significantly involved in the quenching by dilution.
Introduction In a previous paper in this journal,' we reported in brief the possibility of obtaining quantitative information about harmonic chemical oscillations near a supercritical Hopf bifurcation through the use of quenching experiments. At that time, the experiments were not sufficiently precise to allow reliable conclusions to be drawn about the oscillatory behavior. In the meantime, the perturbation experiments have been completely automated as described in the Experimental Section. As a result, we can now determine, with some confidence, the small-amplitude oscillations of the Belousov-Zhabotinskii (BZ) system. The purpose of this paper is to report the more accurate quenching results now available for a number of substances, to select on the basis of these experimental results a set of species that appear to be fundamental to the oscillatory behavior, to calculate the amplitudes and phases of these crucial species, and to provide all the theoretical and experimental details that were omitted in the preliminary paper for lack of space. Quenching experiments are made with limit cycle oscillations in an open system at a flow rate chosen sufficiently close to a supercritical Hopf bifurcation to make the oscillations almost harmonic. A typical experiment is carried out by addition of a small amount of a substance (which by assumption is an important participant of the oscillatory reaction) within an interval of time that is small compared to the period of oscillation. By suitable choice of substance, of the amount added, and of the phase of the (1) Hynne, F.; Graae Sarensen, P. J . Phys. Chem. 1987, 92, 6573.
0022-3654/89/2093-5467$01.50/0
oscillation at the instant of the addition, it is possible to quench the oscillation, Le., to temporarily stop it. Figure 6 shows typical successful results of quenchings where the concentration of Ce4+ is monitored. Recordings like those shown in Figure 6 give an incomplete picture of the state of the reaction system during a quenching experiment. They represent the time evolution of the projection of the phase point onto the [Ce4+]-axis. What is going on in the entire phase space during an experiment is best explained by Figure 1, which applies to some hypothetical three-dimensional phase space. During the limit cycle oscillation before the quenching addition, the phase point circulates the closed curve A, which is essentially an ellipse lying on the unstable manifold r of the saddle focus F; in the small part of the phase space of concern, I? is essentially a plane. From some point on the ellipse, specified by the phase of the oscillation, the phase point is shifted tb a point T o n the stable manifold Z of the saddle focus (which is essentially a straight line) by changing the concentrations through addition of one or more of the chemical species participating in the reactions. A state placed precisely on Z by the perturbation will decay exponentially to F along Z if the system is left completely undisturbed. However, in practice, one cannot hit I: exactly, only get very close; so the phase point will spiral out while it decays toward the plane containing the limit cycle. (Figure 1 is meant for illustration of the principle only: The actual rate of decay along Z in relation to the rate of spiraling out is important to the appearance of experimental results like those of Figure 6.) Note that the proximity of the Hopf bifurcation is essential to the 0 1989 American Chemical Society
