Solvent effect on activated rate processes: On the validity of the GLE approach

Solvent effect on activated rate processes: On the validity of the GLE approach

Chemical Physics 152 ( 199 1) 153- I67 North-Holland Solvent effect on activated rate processes: on the validity of the GLE approach * Gilles Tarjus ...

1MB Sizes 0 Downloads 13 Views

Chemical Physics 152 ( 199 1) 153- I67 North-Holland

Solvent effect on activated rate processes: on the validity of the GLE approach * Gilles Tarjus ’ and Daniel Kivelson of Chemistry and Biochemistry University of CaliJbmia. Los Angeles, CA 90024, USA

Department

Received 10 September 1990

The solvent effect on the prefactor of activated rate processes is usually described in terms of a Langevin equation or a generalized Langevin equation (GLE). The validity of such GLEs in the presence of a strong “external” potential (activation energy barrier) has been verified only for unphysical limiting models. Starting from a microscopic description of the system, we obtain a set of conditions under which a GLE is valid. For effective particles moving along a reaction coordinate, A’, this GLE is applicable to those particles which are initially at fixed position, Xc, and a different GLE with a different friction, C”(Xo, t), is required for particles in different initial positions. Thus, though a GLE may be suitable for describing the reactive motion over a restricted region (piecewise validity), it is of doubtful validity when applied to the whole reaction path.

to the solvent. The friction and random force are related by a fluctuation-dissipation relation

1. Introduction

The dynamic effect of the solvent on activated rate processes in solution was first addressed and consistently treated by Kramers [ I]. In Kramers’ approach, the reaction system is idealized as an effective Brownian particle moving across a one-dimensional potential barrier, V(X), where X represents a reaction coordinate. The influence of the solvent on the reaction dynamics is described stochastically through a Fokker-Planck equation for the distribution of position X and velocity k of the Brownian particle or, equivalently, by a Langevin equation: Mji(t)=-~~(X(t))-M~~(t)+R(t))

(1)

where A4is the mass of the effective Brownian particle, the “external force” ( - V,) is given by -V,(X)=-av(x)/ax,

(2)

and [and R ( t ) are the friction constant and the random force, respectively, associated with the coupling * ’

Supported in part by the National Science Foundation and by the Centre National de Recherches Scientitiques. Permanent address: Laboratoire de Physique Thkorique des Liquids, Universitt Pierre et Marie Curie, 4 Place Jussieu, 7552 Paris Cedex 05, France

030 I-O IO4/9 1 /$03.50

=2W-’

d(t, -tz) f,

(3)

where j?= (k,T) -’ and ( ) indicates an equilibrium ensemble average. The random force satisfies the relations =O ,

(R(t)& (R(f))=O,

(4) (5)

and is, moreover, assumed to be a Gaussian process. Using this description, Kramers was able to find the steady-state escape rate out of the potential well. Kramers’ theory is Markovian since it assumes that the characteristic time scale of the solvent fluctuations is very short compared to the time scale of the Brownian particle motion. Arguing that this latter assumption breaks down for the rapid motion of the effective particle near the top of the barrier, Grote and Hynes [ 2 ] proposed a generalization of Kramers’ theory. Their description starts from a generalized Langevin equation ( GLE ) :

0 I99 1 - Elsevier Science Publishers B.V. (North-Holland

)

G. Tarjus. D. Kivelson / Validity of the GLE approach

154

M/Y(t)= - v,(x(t)) -A4

I

dt, [(t-t,)k(t,)+R(t),

(6)

0

ict,-b)


(7)

3

in which the random force satisfies eqs. (4,5 ) and is assumed to be a Gaussian process. In the GroteHynes treatment only the motion in the vicinity of the barrier top, X=0, is relevant, and the potential V(X) is approximated by V(X) = Vo- jMf.$P

)

(8)

where W, is a frequency related to the curvature of the top of the potential barrier. Combining eqs. (48) with the Gaussian property of R(t), Grote and Hynes derived an expression for the barrier crossing rate constant; the transmission coefficient characterizing deviations from transition state theory is predicted to be AX/o,, where I, is the divergent solution of the GLE for the region near the barrier peak and is given by ~X=&(~X+e(~X))>o

>

(9)

where [(A,) is the Laplace transform of c( ). See also refs. [ 3-8 1. Other authors [ 9- 15 ] have used the GLE approach to calculate the rate constants without restricting its validity to the vicinity of the barrier top, a procedure which, as we shall discuss, we find questionable. Though the description of the motion of the effective particle along the reaction path in terms of a GLE seems reasonable, it is not rigorous, and one encounters difficulties in its use. First of all, it has not been generally proven that in the presence of an arbitrary potential V(X), which we call the external potential, eqs. (6,7) can hold simultaneously with eqs. (4,5). Secondly, the random force R ( t ) is usually not wellspecified and may itself be dependent upon the dynamics of the effective particle in the external field; therefore, in using a GLE one must make additional assumptions of questionable validity, such as the Gaussian nature of R ( t ) and a model description of the friction c(t) in terms of the properties of the liquid bath. So before considering an extension of the functional form of [(X, t), we should focus on the basic applicability of a GLE. The derivation of a Langevin or generalized Lant

gevin equation from exact microscopic equations is a problem that has received much attention [ 16-27 1. Here we reexamine the question, focusing on activated rate processes in solution. For GLE to be of practical value, ( 1) the dynamics entering in the time-dependent friction c(t), should, at the very least, be independent of the externaljeld, V( A’). Additionally, (2 ) the analysis is simplified if the dynamics of the friction i(t) are independent of the motion of the effective particle. Though not a requirement ( 3 ) the GLE can be used most effectively if the dynamics of i(t), in contrast to its value at t =O, depend only upon the dynamics of the pure bath and are independent of the very presence of the stationary tagged particle. In the generalized Stokes-Einstein-Debye treatment of translational and rotational diffusion this is the case; there is no external potential V(X), the friction [depends upon the properties of the bath through the coefficient of shear viscosity q of the pure solvent or bath, and the coupling of the particle motion to the bath is described by a coupling parameter which depends upon the static boundary conditions between particle and bath, but is independent of the particle dynamics. If these conditions do not hold, i.e. if the friction [ depends upon the external field and upon the particle motion, the usefulness of the Langevin approach and the very concept of friction may be of questionable value. It should be noted that the introduction of a time-dependent friction i(t), as in the Grote-Hynes theory, need not necessarily violate any of these requirements; such a time-independent friction can still depend upon the “pure bath” properties through the viscosity, but through a frequency-dependent and k-dependent viscosity, and the static coupling parameter can also be k-dependent. Zwanzig [ 161 and later Lindenberg and co-workers [ 25-271 have shown that, indeed, for a very special model the GLE holds exactly, irrespective of the magnitude of V(X) and of the random (coupling) force R (t ). See also ref. [ 23 1. This model is one in which the bath consists of harmonic oscillators, and in which the potential coupling the effective particle to the bath is bilinear in X and the coordinates (B) of the bath. This appears to be the only case in which the GLE satisfying the above conditions ( 1) and (2 ) holds exactly [ 28 1. On the other hand, Kim and Oppenheim [ 2 1] have shown that for sufficiently small

G. Tarjus. D. Kivelson / Validity of the GLE approach

external field (I’,), for high bath relaxation frequency (rk)-‘, and for large tagged particle mass (M), the GLE also holds. We have sought general conditions under which a GLE is valid and useful. Our studies recapture the exact conditions of Lindenberg and co-workers, as well as the Brownian diffusion limit of Kim and Oppenheim. We also obtain criteria indicating how far a system can deviate from these limits and still be usefully described by a GLE. And, under conditions for which a GLE is valid, we examine the effect on the friction C(t) of the X dependence ofR(X, B).

2. Our approach In this section we outline our approach to the problem, introducing a number of parameters of special interest. In the next four sections we present our results in considerable detail. After that, in section 8 we carry out the step-by-step derivations which lead us to the aforementioned results; by deferring the rather involved procedural details to a late stage in the article we hope to motivate some readers to follow these details, and to save those readers, who either have confidence in our algebra or who lose interest after reading our results, the time and effort needed to read through it all. We then end up with a short discussion of “barrier crossing” in section 9 and a summary of our “conclusions” in section 10. In studying reaction dynamics in liquids in terms of friction, we follow the standard approach and attempt to isolate a reaction coordinate X which is largely decoupled from both the other modes involving the intramolecular motions of the reacting species and those associated with the solvent. As in the Kramers theory, we treat the reaction as the motion of an effective particle moving across an activation barrier along the reaction coordinate, all the other motions being collectively taken as the “bath”, and one regards friction as the result of the coupling of the reactive motion to the bath. It may be that no single, well-separated reaction coordinate can be identified, in which case we can extend the concept of reacting modes to include a set Xof reaction coordinates which are largely decoupled from all other motions. If this cannot be done with a small set of such reaction coordinates, alternative approaches must be

155

sought. We shall focus on the case where one reaction coordinate X suffices. The “bath modes” are not true normal modes of the solvent (modes taken in the generalized hydrodynamic sense or perhaps normal modes such as those introduced by Seeley and Keyes [ 291 and Xu and Stratt [ 301) because they are evaluated in the presence of the effective particle; although motion along the reaction coordinate is explicitly excluded from the “bath modes”, it follows, as we shall discuss, that the bath modes depend upon the presence, but not the motions, of the effective particle. In the gas phase the bath consists entirely of internal modes; in liquids the bath consists of both internal and solvent modes. At least at low frequencies, the behavior of the solvent modes can be treated as that of a hydrodynamic bath and related to wellknown macroscopic liquid properties such as viscosity q(w) and dielectric permeability E( 0); thus the frictional effects can be associated with these properties. At high frequencies, such as those associated with passage over the peak of the activation barrier, the bath may lose its collective, long wavelength behavior, and the relevant frictional properties may not be clearly associated with the macroscopic solvent properties. And, of course, the friction arising from coupling to internal modes probably has little to do with properties such as solvent viscosity and dielectric permeability. Although our formal treatment encompasses friction from both solvent and intramolecular sources, our discussion here will focus largely upon solvent friction. It is useful to describe the bath in terms of two relaxation times, a short time 7; related to the Einstein frequency, and a long time 7% related to the diffusive processes; rf, is a collective property associated with the bulk behavior of the solvent whereas 7: is strongly dependent upon the short-range interaction between the effective particle and the nearby solvent. This isF an oversimplified picture but does seem to incorporate the principal features observed in MD simulations. It should be noted that the bath dynamics, along with these characteristic times, are those for the bath in the presence of a stationary effective particle. Other times which are of interest are il;’ and o;‘, the latter given by eq. (8 ) or more generally by

(10)

G. Tarjus, D. Kivelson / Validity of the GLE approach

156

and the former being approximated by eq. (9); o, is thus related to the curvature of the activation barrier at its peak (X=0), and A, is the actual frequency of barrier-peak passage in the presence of friction. Still another significant time is T, the smallest relevant time after which one is no longer interested in the details of the dynamics affecting i(t) and the response of the bath; this time is given by r= minimum

of r%or I;’

.

(11)

To understand this, note that if, as required, c(t) is independent of the external potential and of the motion of the effective particle, then in a time t> rg the friction is totally relaxed and we need no longer follow its temporal behavior; on the other hand, we are not at all interested in the dynamics after barrier-peak passage has been achieved, and the time for this is 1; ’ . Two interesting limits in which we can examine the GLE are the Kramers-Langevin limit in which the relevant times are very long so that ?=T% )

L-(MP)-l/‘exP(>r)-l

x

1+ F [

-x r3,.

by Hynes

(13)

We find it useful to examine the GLE for a tagged (effective) particle initially placed at a given position X0 along the reaction coordinate. This position might be at the peak of the barrier (XC 0) or any place else, even at the minimum of a potential well. Restating the GLE in this way does not restrict its applicability because once a GLE of this kind is derived, it can be solved for any initial position, and one can then integrate the solutions over any desired distribution of initial positions. To fix the initial position in this way, we introduce an average (A) Oin a lim-

x

,X,,(PM)l/2

1.

(15)

Similar considerations apply to the case where X0 is close to the minimum of a potential well. It is convenient to start with a Hamiltonian of the form

H=[fA4x2+W(X0,6X)]+

C fm,b:+U(X,B), Y

(12)

and the frozen solvent limit, first mentioned and co-workers [ 8 1, in which 74;’

reactive motion in this way because the interaction of the particle with the surrounding bath may be strongly dependent upon the position along X, and an average over all X might mask the interesting dynamics near the barrier peak. The fact that the resulting GLE is localized in this way leads to the concept of piecewise validity, i.e. segmental validity of the GLE. Thus the dynamic length L, the distance traveled by the effective particle in time r, must be short for a GLE to hold; as we shall show, for X0 close to the peak of the barrier,

(16) where W is the gas phase potential, 6X=X-X0, the bath coordinates B= {b,}, and U( X, B) represents all the potential energy contributions over and above W. Next we wish to transform the bath coordinates to get the best normal modes with the tagged particle frozen at some position Xc,. (The actual process for doing this may be a bit vague since in a liquid there is no unique configuration about which the system will make small oscillations, and, consequently, no unique set of modes [ 291. ) The transformed coordinates 6, are linear combinations of the b,‘s, and the new masses M, and the frequencies Q, are functions of the frozen-particle position X0. We can now rewrite U(X, B) as

ited initial ensemble;

0

=

(A4X-&I) <&X-X0) > ’

(14)

where ( ) indicates an average over the complete equilibrium ensemble. In addition to localizing the initial position, we localize the motion along the reaction coordinate to a very small range IX-X0Ix L about X0, where L, the dynamic length, is determined by the maximum time T over which we need follow the dynamics in detail. It is important to localize the

C r,(xo)6,+c,,~bbx(xo, Y +Q~~F~“*(X~, B, 6X)

, 1

B)

(17)

where ebbs bbbis the anharmonicity of the bath when the tagged particle is frozen at X0, and all terms multiplied by A represent couplings between the bath

G. Tarjus, D. Kivelson / Validity of the GLE approach

variables and the particle displacement 6X from X0: that with Fbbx is the coupling term linear in 6X and nonlinear in by’s, and that with Fbxx is all the coupling contributions which are nonlinear in 6X. The terms with e’s represent corrections to the Zwanzig model studied by Lindenberg et al. (harmonic oscillator bath and coupling term bilinear in bath variables and reaction coordinate). The functions F are “normalized” such that ( (Fbbb)*),,= ,, and ((Fb”)*)=((Fbbx)*)= $ ( (Cm,a;E)*) ( (,W,6,)*); the dimensionless strength parameters c are estimated as follows:

(18a) (18b)

157

potential that give rise to deviations from the Zwanzig model. Thus ebbb,which measures the anharmonicity of the bath modes, is described in terms of the cubic term in the potential; the corresponding length @,b-+m if the anharmonic terms vanish, but in the more likely situation that the bath potential depends upon 1b,- b,, I -“, the length o&b/n should be typiCa of interparticle separations, ob, in the bath. If the Xdependence of the potential U(X, B) enters exclusively through variables of the form Ib,- XI “‘, as is likely in diffusive problems, then we would expect the lengths ebb, and cb, to be quite comparable to the characteristic bath distance o&b, but in some chemical reactions, particularly those involving charge transfer, cb,, could be considerably shorter than a&b and frbb,.As a consequence, in the discussion below we will set f&,x c&b, but will retain the distinction in ub,,.

3. Results

(19) The dynamic length L is defined in eq. ( 15 ) and the lengths obbxT obrn nbbb and lb(r) are defined as follows:

GOa)

We find that under appropriate conditions a GLE similar to that in eq. (6) is valid for effective particles initially at X0: MX(t)=

- K(X(l)) I

dt, l’(t--t,)k(t,)+R’(t),

-M

(23)

0

(2Oc)

together with the fluctuation-dissipation (R’(t, )RO(tz) >o =MB-‘CO(t,

-t2)

relation (24)

,

&,(~)=~(l/~tib)“*,

(21)

and the necessary auxilliary relations

@‘ai&=a’

(22)

(R’(t)k),=O,

(25a)

(R”(t))o=O.

(25b)

2

where the subscripts on U indicate derivatives with respect to the reaction coordinate (X) and bath variables B, and ( )Bx is an equilibrium ensemble average taken over all bath coordinates B with X fixed. ab is taken as a typical intermolecular separation in the bath, while &.(7) represents the distance traveled in the time r by a typical effective bath quasiparticle of characteristic effective mass fib. While the coupling parameter li depends upon the actual motion of the system, the coupling parameter @does not. In specifying the magnitudes of the e parameters we have focused on the leading contributions in the

The random force R”(t) is given by the relation

R”(t)=[(exp[(~=,t)l}R(X,B)l~=Xo,

(26)

where R(X, B) = - u,(X, B) + ( u,(X, B)

)Bx

,

(27)

and LB, is the Liouville operator describing the motion of the bath with the effective particle stationary at a point X. V(X) is the potential of mean force. An equation similar to eq. (23) was derived by

G. Tarjus, D. Kivelson / Validity of the GLE approach

158

Adelman [ 3 1 ] by means of a “partial clamping approximation” which consists of two main assumptions: that the coupling between the tagged particle and the solvent is linear in SX, and that the response of the solvent to this perturbation is itself linear. However, Adelman only touched upon the conditions under which these approximations are appropriate; see e.g. ref. [ 321. In our picture R O( t ) and Co( t ) are independent, by construction, of the external potential V(X). Below we show that eqs. (23-27) hold provided the following conditions are satisfied

the friction C”(t ) in the GLE is independent of ebxx,but linearly dependent upon tbbb and $,bx: Under these conditions

~“(t)-(Mj?L2)-‘~2{1+Ebbb+$jbx}.

It is interesting

to examine

;i;‘>rpb,

(30)

with the result that ~=r%. The dynamic eq. ( 15 ) becomes

L- (Mj3-“2r;K(X,)

)

length L in

(31)

where K(X,)=

[

1+

F

1

,xo,(/3A4)L12 .

x

(32)

In order to simplify the discussion we assume that the characteristic lengths ob, o&b, a&, and ab,, are all comparable as, for example, in simple diffusive problems (see above). We assume, too, that K( X0) is not much larger than 1; this is true provided X0 is not too far from an extremum of the potential of mean force. The two inequalities in eqs. (28a,b) become

(29)

Eq. (29), together with the conditions in (28a-c), constitute our principal result. Recalling that the e’s measure deviations from the Zwanzig model and that n measures the strength of the coupling of the effective particle to the bath, our results imply that we can incorporate significant corrections to the Zwanzig model within a GLE, these corrections corresponding t0 Signifkant VdUeS Of ebb6 and tbbx, subject of course, to the limitations above. More specifically we conclude that, provided the coupling n is sufficiently small, the bath anharmonic terms e&b and terms t&x which describe the contributions to the coupling interaction which are nonlinear in bath variables B, can be appreciable without invalidating the GLE. On the other hand, tbxx, the nonlinearity in X in the coupling interaction, must be small for a GLE to be applicable. Similar conclusions were also obtained by Lindenberg and co-workers [ 25-27 ] in their study of model systems in which the bath is described as a collection of harmonic oscillators.

4. Discussion of results: Kramers-Langevin

number of different limits. One extreme case is the Kramers-Langevin limit in which

In the preceding formulas the factor within the square brackets is a measure of e&b (and of eb6x) and thus gives the size of the anharmonicity constituting the coupling between the “normal modes” ofthe bath. In simple liquids where all potentials depend only on the relative distance between particles, there is no reason for the factor in square brackets to be small; it is likely that this factor is greater than one. In the Kramers-Langevin limit $(A,) xC”(0)7pb, which, in conjunction with eq. (9) and the expressions above, leads to

(34) (1,/w,) is a measure of the “frictional strength”, the friction being weak if this ratio is of order unity, i.e. if

limit

the GLE problem

whereas that in (28~) becomes

in a

159

G. Tarjus,D. Kivelson/ Validityof the GLE approach

and strong ifl,lo,c<

and eq. (28~) becomes

1, i.e. if

(35b) Two cases can be distinguished: ( 1) (oxr%) R 1. Since we assumed that (n,r%) K 1, a GLE is valid only if the “frictional strength” is large, i.e. if 1,/w,<< 1. Comparing eqs. (33a,b) with eq. (35b), we conclude that under these conditions a GLE is valid only when the Zwanzig model holds, i.e. when anharmonicity and coupling nonlinearities in 6, and 6X are negligible and when the coupling strength @is enormous; these conditions are unphysical, and so we cannot rely on the GLE if W,?Ekl. (2) (w,rpb) << 1. Eqs. (33a, b) simplify to

(3)‘”@[ (mb;;,,20Ja 1? ($q2[ (mb;;l,20J~ 1.

(36a) (36b)

5. Discussion of results: frozen-bath limit

I,’

with the consequence L-

limit,

<<~f
(Mp)-“*1;‘K(X,)

(37) that .r=I;‘. .

where [ (rf&#) ‘I2 @.,] - ’ is an estimate of ebxxand of the anharmonicity e&b of the coupling between the “normal modes” of the bath. In the frozen-bath limit, which is an extreme limit ofthe r=l;’ case, one finds (‘(2,) x[O(O)/n,, and the various expressions above, combined with eq. (9), lead to

($5 l-(3 42[ (fibp)f”obwJ (40) )

provided that the second term on the right-hand side is less than 1. Once again we take the ratio (AJo,) as a measure of the “frictional strength”, the friction being weak if l/2

0 [ -fib

Here the “frictional strength” may be either weak or strong. Eqs. (36a,b) suggest that the Kramers-Langevin GLE is valid if ( Mb/M)1’2 and @(#r&4) “’ are sufficiently small, and provided the anharmonicity indicator, [T”,/ (tit& “*)ob] N cbbbN ebbx,is not too large. The friction c’(t) in this case corresponds to the usual Brownian friction which is calculated for a frozen tagged particle, and one essentially recovers the conclusions drawn by Kim and Oppenheim [ 2 11.

In the frozen-bath

(39b)

It follows that (38)

we nOW assume that ob, r&,bb,and abb, are comparable, but we let ab,, be different. As before, we assume, too, that K( X0) is not much greater than 1. The inequalities in eqs. (28a,b) become

M

1

@ (&$)“*~bW,

1

a ’ ’

(4la)

and strong if

(4lb) By comparing eq. (41 b ) with eq. (39a), and noting that the exponential factor in eqs. (39a,b) becomes small, one concludes that in the strong friction limit the conditions for a GLE to be valid are unlikely to be satisfied; they may be satisfied only in the limit of the Zwanzig model for which the anharmonicities of the bath and the nonlinearities in 6, and 6X of the coupling potential are negligible. The coupling strength @must, therefore, be enormous unless (@,/ M) is itself enormous (compare eq. (41 b) with eqs. (39a,b)), and with such an enormous solute-bath coupling, it is doubtful that the bath could be frozen on the relevant time scale 2; ’ . In the weak-friction limit, ;1,-w,, eqs. (39a,b) yield

(42a)

G. Tarjus, D. Kivelson / Validity of the GLE approach

160

with which, together with eq. (41a), are the conditions for a GLE to be valid. Therefore, a GLE in the frozenbath, weak-friction limit is best satisfied by large o,, and the restrictions on (r~$lM)“*, @,and (r~~~~/cr~b) become less severe the larger the w,. In physically meaningful situations, [ 7$/ (H&,/3) 'I* ab] 2 1, and so the frozen-bath-limit requirement that A,- o,., together with eq. ( 37 ), yields

and eqs. (42a,b)

7s,
Often it also appears liquids,

I

dtr”(t)=

qdtiQ(r).

0

viscous

(48)

0

Under the conditions of eq. (44) we can get a rough picture of the friction by setting

E) x

(43b)

that for sufficiently

cc

fO(n,)r

become

(47)

b

(49)

which enables us to rewrite eq. (9) as

(43c) It should be noted that the inequality incorporated in that in eq. (43a).

in eq. (43b) is (50) where n,r% > 1 and [ (t?i&) “2 $,&] - ’ - cbbb- cbbxis the anharmonicity parameter. Under the conditions of eq. (45 ),

6. Discussion of results: intermediate cases Here we examine our results in situations intermediate between the Kramers-Langevin and the frozen-bath limits. We shall assume that r=l;‘, with a subdivision into two cases, one for which rs,c
(51) which enables us to rewrite eq. (9) as

(44)

and one for which A;'
(45)

The first of these represents an extension of the Kramers-Langevin model, the latter an extension of the frozen-bath model. In order to get a simple enough picture, we again assume that the scale lengths ab, r&b, and ebb, are comparable, but we let obXXbe different, and we assume too that K(&,) is not much greater than 1. The inequalities in eqs. (39a,b) hold in both these cases since they were obtained under the conditions r=A; ’. It appears that in some cases, at least, one can rewrite Co( t ) as i"(~)=rS(~17Sb)+~4cQ(~17%),

(46)

={

(~bp)!,2ubaJ++ l+(g@‘[ 3

(52) where A,7$ > 1. Generally, the inequalities (39a and b) are most easily satisfied if - (fib/M) is small; - $ and (cbxr/$) are SInall; - co, is large; and - ,%,/w, is not too small, i.e. the friction is not too great. And, indeed, by comparing eqs. (39a,b) with eqs. (50) and (52), we see that for the GLE to be valid

161

G. Tarjus, D. Kivelson / Validity of the GLE approach

in the strong friction limit, # must be enormous and the anharmonicity parameter [ ( P@,/?)‘I2a*&] - ’ must be negligibly small. This brings us back to the limit of the Zwanzig model. In the discussion above, we made some simplifying assumptions. We assumed that there is no marked difference between the scale lengths 06, f&b, and cbbx, and that the factor K(&,) entering in the dynamic length L is not much larger than 1; but, if X,, is taken far from an extremum of the potential of mean force, K( X,,) may be large and the conditions (28a-c) more difficult to satisfy. We assumed too that there is no appreciable difference between the “back” and the “forward” coupling, i.e. that the magnitude of the effect exerted by the effective particle on any bath particle is comparable to that of the effect of all bath particles on the effective particle. In so doing we neglected the fact that the order of magnitude of the two effects may actually differ by a factor of NC,NC being the effective number of bath particles interacting significantly with the effective particle.

7. Relation to viscosity A striking observation is that the prefactor (i.e. the transmission coefficient) of the rate constant in a number of reactions in solution, especially conformational changes, seems to vary monotonically and smoothly with the solvent viscosity 4( 0) over a reasonably wide range of viscosities. In some cases, the transmission coefficient is well described by a Kramers formula in which the zero-frequency reactive friction [( 0 ) is taken as proportional to the zero-frequency solvent viscosity Q(O) [ 33,341. This association can be understood if one assumes that both the viscosity r,r(t) and the friction, c’(t) have large characteristic low-frequency components which are proportional to each other and small, well-separated highfrequency components, as indicated in eqs. (46)(48 ). Thus if the frequency of interest il, is very low, then f”( 0) N 4(O). In other cases, however, the transmission coefficient can be related to the zero-frequency solvent viscosity Q(0)) but with a weaker, nonlinear dependence [ 33,35-401 (e.g., it may vary as Q(0) --a with 0
from the fact that the frequency of interest A, may be sufficiently low that the contributions of the high frequency components of both co (0) and rj( w ) are unaltered but the low frequency contributions of both are somewhat diminished from their zero-frequency values; in this case $(~,)=P(O)+f(4(0)&

7

(53)

where f is a nonlinear function of q( 0). An expression for f( 6) has been obtained by introducing frequency-dependent transport coefficients into a hydrodynamic description [ 411, but it has not been successful in quantitatively reproducing experimental data [ 34-36 1, see also ref. [ 42 1.

8. Procedure

In this section we reformulate the equation of motion for effective particles initially at X=X0 in such a way that this equation takes a form analogous to a GLE. The system under study consists of an effective particle, characterized by a one-dimensional coordinate X, coupled to a bath, and subject to an external potential V(X) which is independent of the bath. Actually, since V(X) is the potential of mean force acting on the effective particle, it contains the equilibrium or averaged effect of the bath, and so we also refer to it as the “renormalized” external potential. We start by writing the classical Hamiltonian as in eq. ( 16), and we introduce the mean force V,(X) acting on the effective particle (see eq. (2) ) and the residualforce R (X, B) acting on the effective particle (see eq. (27 ) ). The equation of motion for the effective particle then takes the form: Mx(t)=-V,(X(t))+R(X(t),B(t)).

(54)

By construction, the residual force satisfies the following equation: (R(X B) )~x=0 3

(55)

where, as indicated above ( ) BXis an average taken only over the bath coordinates Band velocities b. As already mentioned, we seek a limited GLE in which only the trajectories of the effective particle starting initially at a given position X0 are considered. Accordingly, we introduce the average ( ). specified in eq. (14).

G. TarJus, D. Kivelson / Validity of the GLE approach

162

Following the Zwanzig-Mori approach [ 16,17 1, we can, with the aid of a projection operator Qo, divide theforceR(X(t),B(t))appearingineq. (54),into a “random” force and a “dissipative” term. We detine Qo as an operator which when acting on A yields QoA=A-

s

k=~-hfj3(~k),/i:.

In order to obtain the approximate (23-27), we rewrite the exact equation eq. (59) as am=-V,(X(I))-Mldt,

(56)

(‘(t-tr)k(t,)

f

0

A “random”

+R’(+iUJdt,

force R + ( t ) can then be defined as

GLE in eqs. of motion in

6i+(t-t,)%(t,)+6R+(t), 0

B) ,

R+(t)=exp(iQ&)R(X,

(61)

(57)

where iL is the full Liouville operator associated with the full Hamiltonian H given in eq. ( 16). We have set ( k2) = 1/A4p. With this formalism, the equation of motion in eq. (54) can be rewritten as eq. (6), but with R(t) replaced by R + (t) and i(t) replaced by

wherelO and (27),

c+(t)=-P(

and

(iLR+(t))k)o.

(58)

Note that in the preceding formula we cannot simply permute the iL operator so that it acts on Malone be). imcause the averaging process indicated by ( plies that _%is multiplied by 6(X- X0). We thus have the GLE KY(t)=-V,(X(l)) f -M

s

dt, l+(t-t,)k(t,)+R+(t)

0

(59) along with the auxiliary (R+(t)k)o=O,

condition

isdefinedineq.

ineqs.

(26)

SR+(t)=R+(t)-exp[(iL,,t)]R(B,X) =R+ (t)-R(X,B’“‘(t))

s[+(t)=--P(

,

(62)

[iL 6R+(t)]k),=i+(t)-r’(t), (63)

I?@) (t) is defined by eq. (62). Recall that the conditions in eqs. (25a,b) still hold and the effective particle has been taken with initial position X0. Our task is to show that the correction terms SC’ (t ) and 6R + (t) in eq. (6 1) are negligible under the conditions in eqs. (28a-c), and that CO(t) is independent Of tbXXbut Contains terms linear in e&b and t&x. By using the fact that R and iL& are independent of jr, as well as the identity exp[ (C+D)t]

=exp(Ct)

I

(60)

which is analogous to eq. (4). Although eqs. (59) and (60) are exact and fulfill some GLE requirements, they have severe shortcomings: the “friction coefficient” [’ (2) defined by eq. (58) is not related to the “random” force R + (t) by a fluctuation-dissipation relation such as eq. (7), and, furthermore, c+(t) may involve the external force VX(X). Therefore, an approximate GLE is needed. We note, however, that if V,= 0, we can then introduce the projection operator Q for which the equilibrium ensemble averages are not restricted to X=X,, and in this case eq. ( 58 ) does indeed have the form of a fluctuation-dissipation relation; however, c’ (t) may still depend strongly on the motion of the effective particle.

(24),R”(t)

+

I

drexp[(C+D)(t-r)]Dexp(Cr),

(64)

0

it is a simple matter to reexpress the “additional dom force” 6R + ( t ) as

ran-

I

6R+(t)=

jdt,

{exp[iQ&(t-tr)l}*

0

(65) where the symbol A,{A} denotes the fluctuating of A with respect to the average ( ) o:

part

&#=A-

(66)


163

G. Tarjus, D. Kivelson / Validity of the GLE approach

After some algebraic manipulations and the use of eqs. (63 ) and (65 ), we obtain the following expression for the “additional friction coefficient” SC’ (t):

and L sin(ti,t) cos(f&z) +!I, _ W”

&‘(t)=6, +Ar +

I’

8 Y f

+

dtr ( [iL,exp[iQ&=(t-t,)]

sin[Wt-h)l

dt

0

1-cos(U) rn,dj:

1

m,@

1 {-EbbbFb,bb(B’“‘(tl))

(MB)“*k

0

+/it,,F;xx(@“‘(t,),

x~{(R(X,~‘~‘(~,))R),)I(~P)“*~:)O, (67) where the Liouville operator iLti is defined as

a iL&ax+M-‘[-vx(X)]

a z.

Eqs. (65) and (67) are convenient starting expressions for order of magnitude estimates of the effects on a GLE of deviations from the Zwanzig model. We first rewrite the fluctuating term appearing in eq. (65) as

aw.(xm, =-Ao{UXX(x,~'"'(t))}+Ao ax - ;4{u,(xJ'")(t))ypj.

(69)

U,, U, and U,, are all obtained by differentiating eq. ( 17)) and they are all proportional to A. The following two formulas are convenient for studying the effects of the deviations from the Zwanzig model:

8x)8x)

,

(71)

where r,, Fbbb, Fbbx, and Fbn depend upon the initial position X0, and the subscript v on F, indicates a derivative with respect to 6,. We treat all contributions corresponding to the Zwanzig model (i.e. all terms that are independent of the t’s) exactly, and we seek order of magnitude estimates for the additional contributions. The analysis is simplified by the fact that the time over which one has to study the equation of motion is bounded (t> r), as are the positions of the effective particle. As mentioned earlier, we have only to follow the motion of an effective particle over a distance L, which we called the dynamic length. For X0 close to the peak of the barrier (X=0), an estimate of the dynamic length traveled during the time 7 can be obtained by assuming a parabolic form for the potential of mean force and neglecting all details of the motion of the effective particle except the locally diverging contribution characterized by a reactive frequency A,. This gives rise to eq. ( 15 ). By making use of the dynamic length and all the dimensionless strength parameters introduced in eqs. ( 18 ) and ( 19 ), an order-of-magnitude estimate is readily obtained: (72)

BLR’(t)-/l{l+tbbb+~bbx},

fiL6R + (t) -A{$,,

+/it&b +&,b,}

Eq. (29) is obtained Furthermore,

.

in a similar

(73) fashion.

G. Tarjus, D. Kivelson / Validity of the GLE approach

164

(74)

x~2{~b,,+~~bbb+A~bbx}

(75)

.

The fact that R’(t) and [O(t) do not depend on e& is an exact result. Finally, by inserting eqs. (29 ), (74) and (75 ) into eq. (67), and by making the further simplifying statements that

(76a)

iLoxexp[i~,tlQo-w,exp(o,t),

=

and (Mj?)“*Xw

[ 21) also contains a finite-width transition region, it differs quite markedly, in principle, from our measurement function description; the “stable states” picture corresponds to a discontinuous measurement function (ax(X) being essentially equal to 1 within the transition region and to 0 outside) whereas no such restrictions are made in our approach. We thus must confront the problem of applying the piecewise valid GLE to the evaluation of ( ax( t),k( t)&aX). We can make use of the relation

I dXo(&X-&d

>aAXo) (ax(t)%t)&09 (78)

1)

(76b)

the magnitude of the “additional be approximated as (M8L*)~i+(t)-~*{~bxx

friction” SC+ (t) can

+Acbbb+Acbbx)

exp(&r)

.

(77) Comparison of eqs. (72) and (73) and eqs. (29) and (77) leads to the conclusion that for the equation of motion to reduce to the GLE in eq. (23), the inequalities given in eqs. (28a-c) must be satisfied.

9. Barrier crossing As we discuss elsewhere [28] the first order rate constant for a process in which a barrier is crossed is proportional to Jr dt (a,( t)x( t)h,), where a(X) is the measurementfunction, and trx=da/dx. (This a is to be distinguished from the lengths introduced above). The measurement function is related to the species specification function; it is an indicator of the probability that at a certain position along the reaction coordinate, X, the effective particle should be taken as the “left-hand” species. The function a(X) is often taken as a step function and the associated a,(X) as a delta function, but, as discussed elsewhere, we believe it should actually be a narrow function, but one with finite width AX about the barrier peak located at X= 0. (More rigorously AX should be function the width of the taken as a,(X) exp ( - /3V( X) ), where V(X) is the potential of mean force. ) Although the “stable states” picture introduced by Hynes and co-workers (see e.g. ref.

and we assume that for all relevant positions X0, the dynamic length L=L(X,) satisfies eqs. (28a-c), so that a GLE is valid for each X0. The piecewise GLE then yields (a,. t)k( t)k),, and two distinct limits need be considered: ( 1) If AX< L,the friction coefficient co ( X0, t ) does not vary over AX, and a unique GLE with a friction coefficient calculated at the top of the barrier can be used in evaluating the expression in eq. (78); it is then easy to show in the parabolic-barrier approximation that the Grote-Hynes formula for the transmission coefficient is recovered (see e.g. ref. [ 7 ] ) . (2) If AX> L, the friction coefficient co ( X0, t ) changes along the relevant portion of the reaction path, and one must integrate the solution of the equation of motion as indicated by eq. (78). This provides a simple extension to the Grote-Hynes treatment. According to the preceding analysis, the rate constant obtained through the integration indicated in eq. (78), may be dependent upon AXand, therefore, upon the choice of the measurement function. For barrier crossing processes, however, physical considerations suggest that AX should not be much larger than the barrier width (PM) --L’2w; I. The case where a GLE is valid with AX>> L is thus more easily attained in the Kramers-Langevin limit when the barrier is broad. From previous studies [ 13- 15,43,44], we expect a piecewise treatment of barrier crossing processes (with or without a GLE) to be valid except in two extreme cases: extremely low friction and extremely high friction in a very slow relaxing bath. In both cases, a description of the reactive dynamics over

G. Tarjus, D. Kivelson / Validity of the GLE approach

the whole reaction path is required. Our analysis of barrier crossing can break down in the extremely low friction limit because we have not considered the possibility of long-range coherence and oscillation in the motion along X, a possibility that can become a reality if the frictional damping becomes negligible. This problem was addressed first by Kramers. In any case, our treatment indicates that if the friction is extremely high and the relaxation of the bath extremely slow, there may be difficulty in applying a piecewise GLE. Finally, if for any relevant position X0 the dynamic length L does not satisfy the conditions (28a-c), no simple GLE approach is possible. By introducing a measurement function and its derivative o,, we can guarantee a piecewise analysis of the motions, but there is no advantage in choosing a very narrow a,. If the AX associated with a, is much less than the dynamic length L, then the very short time dynamics of ( a,( t )k( t )kuJ may be correctly described by a GLE, but it is strongly affected by the presence of a,(t) so that knowing this very short time dynamics is not sufficient to determine the rate constant. This latter factor requires us to follow reactive trajectories during times of the order of A; ’ , i.e. over distances of at least L.

10. Conclusions

We have examined conditions under which a generalized Langevin equation (GLE ), which incorporates the concept of friction, can serve as the effective equation of motion for the study in liquids of activated rate processes along a reaction coordinate. Since a GLE holds exactly for the Zwanzig model in which the liquid is described by harmonic modes and the coupling of the effective particle and the liquid bath is bilinear, we examined the effect of deviations from this model on the GLE and on its validity. In liquid solutions there is no a priori reason to expect deviations from the Zwanzig model to be small. In our approach we studied GLE’s governing trajectories of effective particles all starting from some given point X0 along the reaction path; as discussed, this is not a limitation since one can ultimately average over a distribution of starting points X0. Because the friction C”(X0, t) varies over the reac-

165

tion path, we question whether a single GLE can be used to describe reactions in real physico-chemical systems. Appreciable variations in (‘(X0, t) are expected if the activated rate process involves an intrinsically nonuniform environment, as in processes taking place in biological systems (such as ligand migration in biomolecules [45,46] and in membrane transport processes [ 47 ] ), or in the desorption of an adatom from a surface [ 48-501. But in other chemical reactions, as well, such effects as hydrodynamic interactions, screening, or charge migration along the reaction path make the friction position dependent. Possible exceptions are conformational changes in which the mobile segment is much smaller than the fixed one but much bigger than the solvent particles. Qualitative information on the position dependence of the friction are available for a few simple models [20,49,51-541, but a quantitative test has been directly provided only by molecular dynamics simulations [ 55,561; these calculations indicate that for molecular conformational changes in solution, the friction varies along the reaction path. As shown by Straub et al. [ 551, a friction coeficient calculatedfor the bottom of one well would be different from that calculated for the top of the barrier or the bottom of the other well (this was first pointed out by Grote and Hynes [ 2 ] ). As a result of the dependence of C”(X0, t) upon position, we believe that a GLE must be applied piece-wise to describe the reactive motion in activated processes. Besides the positional dependence of the friction Co( X0, t ), there are other reasons to eschew the use of a single GLE along the entire reaction path. In developing conditions under which a GLE is valid, we set the factor K( X0), defined in eq. (32), to be of order of one; this is reasonable when X0 is taken close to an extremum in the potential of mean force, but not if X0 is taken far from these extrema where V,(X,) is large. Large values of K( X0) increase the difficulty of satisfying the inequalities (28a-c), inequalities needed for a GLE to be valid. The approximate expression for the friction Co(t) which is associated with a GLE appropriate to a system not described by the Zwanzig model is given in eq. (24) and the conditions under which such a GLE is valid are given in (28a-c). These conditions are dependent upon a number of relevant physical parameters, and they can most easily be understood in

166

G. Tarjus, D. Kivelson / Validity of the GLE approach

various limiting cases. We find, not unexpectedly, that if

and the anharmonicity indicator rpb/(@$) ‘/‘a, is not too great, then the simple Kramers-Langevin theory with frequency-independent friction is valid (see eqs. ( 36a,b) ). The frictional strength in this case may be either weak or strong. Under the conditions

we have a system in the frozen bath limit which can ordinarily also be described by a GLE, but only for weak enough friction.

Note added in proof

Since this paper was submitted, work by Beme and co-workers has been published [ 571 in which a study is presented of the dynamic friction on a harmonic oscillator dissolved in a Lennard-Jones fluid, they show that for this system the frozen-bond friction is the asymptotic high-frequency limit of the true dynamic friction. Because of our focus on the validity of the GLE we have omitted references to many articles in which the emphasis is placed on the application of the GLE to reaction problems.

References [ 1 ] H.A. Kramers, Physica 7 ( 1940) 284. [ 2 ] R.F. Grote and J.T. Hynes, J. Chem. Phys. 73 ( 1980) 27 15. [ 31 J.T. Hynes, Ann. Rev. Phys. Chem. 36 ( 1985) 573; ibid., in: Theory of Chemical Reaction Dynamics, ed. M. Baer (Chemical Rubber, Boca Raton, FL, I985 ), Vol. IV, p. 17 1; ibid., J. Stat. Phys. 42 ( 1986) 149, and references therein. [ 4 ] P.G. Wolynes, Phys. Rev. Lett. 47 ( 198 1) 968. [ 5] P. Hanggi and R. Motjabai, Phys. Rev. A 26 ( 1982) 1168. [6] E. Pollak, J. Chem. Phys. 85 (1986) 865; E. Pollak, S.C. Tucker and B.J. Beme, Phys. Rev. Letters 65 (1990) 1399. [ 71 D. Borgis and M. Moreau, Mol. Phys. 57 (1986) 33. [8] B.J. Gertner, K.R. Wilson and J.T. Hynes, J. Chem. Phys. 90 ( 1989) 3537, and references therein.

[9] B. Carmeli and A. Nitzan, Phys. Rev. Lett. 49 ( 1982) 423; ibid., J. Chem. Phys. 79 ( 1983) 393; ibid., Phys. Rev. A 29 (1984) 1481; ibid., J. Chem. Phys. 80 (1984) 3596. [IO] F. Marchesoni and P. Grigolini, J. Chem. Phys. 78 ( 1983) 6287. [ 111 E. Guardia, F. Marchesoni and M. San Miguel, Phys. Letters A 100 (1984) 15. [ 121 M.M. Dygas, B.J. Matkowsky and 2. Schuss, J. Chem. Phys. 83 (1985) 597; ibid., SIAM J. Appl. Math. 46 (1986) 265. [ 13 ] J.E. Straub, M. Borkovec and B.J. Beme, J. Chem. Phys. 84 ( 1986) 1788; addendum and erratum, ibid. 86 ( 1986) 1079. [ 141 R. Zwanzig, J. Chem. Phys. 86 (1987) 5801. [ 15 ] P. Talkner and H.-B. Braun, J. Chem. Phys. 88 ( 1988) 7537. [ 161 R. Zwanzig, Phys. Rev. 124 (1961) 983. [ 171 H. Mori, Progr. Theoret. Phys. (Kyoto) 33 (1965) 423. [ 181 R. Kubo, Rept. Progr. Phys. 29 ( 1966) 255. [ 191 P. Mazur and I. Oppenheim, Physica 50 ( 1970) 241. [20] J.M. Deutch and I. Oppenheim, J. Chem. Phys. 54 ( 1971) 3547. [21] S. Kim and I. Oppenheim, Physica 57 (1972) 469. [ 22 ] U. Mohanty, K.E. Shuler and I. Oppenheim, Physica A 115 (1982) 1. [23] G.W. Ford, M. Kac and P. Mazur, J. Math. Phys. 6 ( 1965) 504. [24] R. Zwanzig, J. Stat. Phys. 9 (1973) 215. [ 25 ] K. Lindenberg and V. Seshadri, Physica A 109 ( 198 1) 483. [ 26 ] K. Lindenberg and B.J. West, Physica A 126 ( 1984) 489. [27] E. Cortes, B.J. West and K. Lindenberg, J. Chem. Phys. 82 (1985) 2708. [28 ] G. Tarjus and D. Kivelson, to be submitted. [29] G. Seeley and T. Keyes, J. Chem. Phys. 91 (1989) 5581. [ 301 B.-C. Xu and R.M. Stratt, J. Chem. Phys. 92 ( 1990) 1923. [ 3 11 S.A. Adelman, J. Chem. Phys. 7 1 ( 1979) 4471; ibid., Advan. Chem. Phys. 44 (1980) 143; ibid., Advan. Chem. Phys. 53 (1985) 61. [ 321 C.L. Brooks and S.A. Adelman, J. Chem. Phys. 76 ( 1982) 1007. [ 331 G.R. Fleming, S.H. Courtney and M.W. Balk, J. Stat. Phys. 42 ( 1986) 83. This paper reviews most of the literature on isomerization reactions in solution published before 1985. [ 341 D.P. Millar and K.B. Eisenthal, J. Chem. Phys. 83 ( 1985) 5076. [35] G. Rothenberger, D.K. Negus and R.M. Hochstrasser, J. Chem. Phys. 79 (1983) 5360. [ 361 S.R. Flom, A.M. Brearley, M.A. Kahlow, V. Nagarajan and P.F. Barbara, J. Chem. Phys. 83 (1985) 1993. [ 371 J. Lee, S.-B. Zhu and G.W. Robinson, J. Phys. Chem. 91 (1987) 4273. [38] D.M. Zeglinski and D.H. Waldeck, J. Phys. Chem. 92 (1988) 692. [39] S.K. Kim and G.R. Fleming, J. Phys. Chem. 92 ( 1988) 2168. [40] C.-L. Xie, D. Campbell and J. Jonas, J. Chem. Phys. 92 (1990) 3736. [41] B. Bagchi andD.W. Oxtoby, J. Chem. Phys. 78 (1983) 2735. [42] R.F. Grote, G. van der Zwan and J.T. Hynes, J. Phys. Chem. 88 (1984) 4876.

G. Tarjus, D. Kivelson / Validity of the GLE approach [43] S. Okuyama and D.W. Oxtoby, J. Chem. Phys. 84 (1986)

5830. [44] J.L. Barrat, Chem. Phys. Letters 165 (1990) 551. [45] B. Gavish, Phys. Rev. Letters 44 (1980) 1160. [46] P. Han&, J. Stat. Phys. 30 (1983) 401. [47] M. Vertenstein and D. Ronis, J. Chem. Phys. 85 (1986) 1628. [48] E.G. d’Agliano, P. Kumar, W. Schaich and H. Suhl, Phys. Rev.Bll (1975)2122. [ 491 C. Caroli, B. Roulet and D. Saint-James, Phys. Rev. B 18 (1978) 545. [ 501 B. Carmeli and A. Nitzan, Chem. Phys. Letters 102 ( 1983) 517. [ 5 1] DC. Knauss and G.T. Evans, J. Chem. Phys. 72 ( 1980) 1499; G.T. Evans, J. Chem. Phys. 72 ( 1980) 3849.

167

[52]G.Moro,Chem.Phys. llS(1987) 167and181. [53] R.I. Cukier, R. Kapral and J.R. Mehaffey, J. Chem. Phys. 73 (1980) 5254; M. Schell, R. Kapral and R.I. Cukier, J. Chem. Phys. 75 (1981) 5879. [54] H. Posch, F. Vesely and W. Steele, Mol. Phys. 44 ( 1981) 241. [ 55 ] J.E. Straub, M. Borkovec and B.J. Beme, J. Phys. Chem. 9 I (1987) 4995. [ 561 S.-B. Zhu, J. Lee and G.W. Robinson, J. Phys. Chem. 92 ( 1988) 2401; ibid., J. Chem. Phys. 88 (1988) 7088; ibid., J. Phys. Chem. 93 (1989) 164. [57] B.J. Beme, M.E. Tuckerman, J.E. Straub and A.L.R. Bug, J. Chem. Phys. 93 (1990) 5084.