7680
J. Phys. Chem. 1996, 100, 7680-7690
Structure and Dynamics of Solvent Landscapes in Charge-Transfer Reactions Vitor B. P. Leite† and Jose´ N. Onuchic* Department of Physics, UniVersity of California, San Diego, La Jolla, California 92093 ReceiVed: May 11, 1995; In Final Form: December 7, 1995X
The dynamics of solvent polarization plays a major role in the control of charge-transfer reactions. Although in principle solvent dynamics looks extremely complicated, the success of Marcus theory describing the solvent influence via a single collective quadratic polarization coordinate has been remarkable. Onuchic and Wolynes have recently proposed (J. Chem. Phys. 1993, 98 (3), 2218) a simple solvent model demostrating how a many-dimensional complex system composed of several dipole moments (representing solvent molecules or polar groups in proteins) can be reduced under the appropriate limits into the Marcus model. This work presents a dynamical study of the same model. It is shown that an effective potential, obtained by a thermodynamic approach, provides an appropriate dynamical description. At high temperatures, the system exibits effective diffusive one-dimensional dynamics in this effective potential, where the Born-Marcus limit is recovered. At low temperatures, a glassy phase emerges with a slow non-self-averaging dynamics. At intermediate temperatures, we will discuss the concept of equiValent diffusion paths and polarization-dependent effects. The equivalent paths are necessary to reduce the problem into the Marcus picture. A discussion of how these different regimes affect the rate of charge transfer is presented.
I. Introduction Electron transfer (ET) reactions are among the most common, simplest, and best-studied chemical reactions.1 Particularly, the effect of the solvent environment in ET reactions has been extensively investigated over the past three decades.2 Marcus theory successfully explained the role of solvent in ET as a single collective reaction coordinate.3 It divides the environment degrees of freedom coupled to the reaction in two groups. “Inner sphere” coordinates are intramolecular vibrational modes (mostly high-frequency modes). The other group is represented by a single “outer sphere” continuous classical collective coordinate, which describes the rearrangement of the remaining degrees of freedom and accounts mostly for the solvent polarization. The total solvent polarization around the donor and/or acceptor cavity is used as the electron-transfer reaction coordinate. A typical (unidimensional) representation of polarization reaction coordinate x is shown in Figure 1. The position where the donor and acceptor diabatic surfaces cross, called x*, defines the transition state for the reaction. The classical rate in the nonadiabatic limit is proportional to
[
2 k ∝ TDA exp -
]
E(x*) - E(xeq) kBT
where TDA is the matrix element between the donor and acceptor electronic states. E(xeq) and E(x*) are the energies at equilibrium and reactive state polarizations, respectively. Usually polarization is treated as an effective overdamped harmonic mode. Similar approaches have also been used to compute other charge-transfer reactions such as proton transfer.4 Although widely used, the question of why this description has been so successful is still open. Solvent dynamics happen on a very high-dimensional surface, and this complex landscape is populated by a large number of minima.5 How can this landscape be accommodated into a single-coordinate picture? This question has been recently addressed by Onuchic and Wolynes.6 † Current address: Departamento de Quimica, Universidade Federal de Sao Carlos, SP 13965 905, Brazil. X Abstract published in AdVance ACS Abstracts, April 1, 1996.
S0022-3654(95)01316-5 CCC: $12.00
Figure 1. Typical diabatic surfaces for donor and acceptor electronic states in electron-transfer reactions. x is the solvent polarization reaction coordinate, and x* and xeq are the reactive state and equilibrium polarizations, respectively.
They considered a polar solvent interacting with a charged cavity, representing the donor or accepor site for ET. A simple analytical discrete model was used to treat each solvent molecule independently. For every temperature, this model has two regions in polarization space: a glasslike region with a slow dynamics and another one exhibiting normal diffusive dynamics where the linear response picture of Born-Marcus is recovered. Onuchic-Wolynes (OW) studied the polarization-dependent thermodynamic phase transition for this model. Since the reactive region has a slow dynamics below the transition temperatures, they argued that this phase transition probably influences the dynamics of the system. Here a procedure to obtain a dynamical phase transition is presented, and a comparison between the dynamical and thermodynamic phase transitions is discussed. The organization of this paper as follows. Section II gives an overview of the OW model and the random energy model (REM). In section III, the kinetics rules are defined and the polarization space is studied assuming no interaction energy (entropic effects only, i.e., at infinite temperature). In this limit the reduction to an unidimensional coordinate representation is exact, and the kinetic effective potential is similar to the thermodynamic one. In section IV, energy interactions between dipoles are considered (finite temperature). The transition from unidimensional to glassy regime is discussed. The concepts of local and global phase transitions and equiValent diffusion paths © 1996 American Chemical Society
Solvent Landscapes in CT Reactions
J. Phys. Chem., Vol. 100, No. 18, 1996 7681 where n+ (n-) is the number of dipoles oriented inward (outward) of the charge cavity. The average solvent energy is
E(p) ) pξ(q) +
zJ 2 p N
(5)
or Figure 2. Single-shell OW model for solvent dipoles around a cavity with charge q.
are presented. Three kinetic critical temperatures are discussed. Section V summarizes the results and presents comments on the implications and relevance to real solvents. For clarity, the first paragraph of each section gives a summary of its results. The reader not interested in the technical details might read just the first paragraph of each section and skip to the conclusion. II. Onuchic-Wolynes Thermodynamic Model for Solvent Section summary: The OW model considers a single shell of solvent molecules with simple rotational dynamics. They are represented by dipoles pointing only in two directions, inward and outward, i.e., as Ising spins (Figure 2). The interaction energies between dipoles and between each dipole and the charge cavity are assumed to be rather simple, allowing for an analytical solution. The total polarization is used as the reaction coordinate for the solvent, and the free energy depending parametrically on the total polarization is used as the effective potential. Although this model is far from representing details of real solvents, it includes the basic features of a roughenergy landscape: multidimensional degrees of freedom, with each solvent molecule being treated independently, a disordered energy landscape with multiple minima, and a polarizable medium around a charge cavity. Above the thermodynamic glass transition, it recovers the continuum dielectric limit. At high temperatures, the system exibits effective diffusive onedimensional dynamics, where the Born-Marcus limit is recovered. At low temperatures, a polarization-dependent glassy phase (polarization-dependent phase transition) appears and a slow non-self-averaging dynamics is expected. For the single shell of N Ising spinlike dipoles, the total number of states is 2N. The solvation energy is
i)1
(1)
〈ij〉
where the first term is the dipole-charge interaction and the second one is the dipole-dipole interaction. q is an index associated with the cavity charge. The charge-dipole interaction (ξi(q)) is assumed to obey a Gaussian distribution, with average ξ(q) and standard deviation ∆ξ(q). The dipoledipole interaction is assumed to be of the form
Jij[σz(i), σz(j), ri,rj] ∼ Jij(ri,rj) σz(i) σz(j)
where x ) p/N is the average polarization per dipole. Following the REM7,8 approximation to evaluate the solvent energies, the energies are given by Gaussian-distributed random variables. The standard deviation of the solvation energy is assumed to be independent of x
∆E2 ) N∆2 ) N[∆ξ2 + z∆J2]
(6)
The total-polarization dynamics is considered slow compared to the other relevant time scales of the solvent. Therefore, the system is studied parametrically on x. The energy probability distribution g(x,E) at polarization x is
g(x,E) )
1
x2π∆E
[
exp -
(E - E(x))2 2∆E2
]
(7)
leading to an average number of states with polarization x and energy between E and E + dE of
〈n(x,E)〉 ) Ω(Nx) g(x,E) dE
(8)
N! N! ) n+!n-! [N(1 - x)/2]![N(1 + x)/2]!
(9)
where
Ω(Nx) )
is the total number of states with polarization x. For 〈n(x,E)〉 . 1, one can approximate 〈log n(x,E)〉 by log〈n(x,E)〉, and the entropy becomes
S(x,E) = log〈n(x,E)〉 = log Ω(Nx) -
(E - E(x))2
(10)
2∆E2
At the critical energy
N
Esolv ) -∑ ξi(q) σz(i) + ∑ Jij[σz(i), σz(j), ri,rj]
E(x) ) N[xξ(q) + zJx2]
(2)
where each dipole interacts with only z neighboring dipoles via a random interaction, which also obeys a Gaussian distribution with average J and standard deviation ∆J. The total polarization p is defined by
p ) n+ - n-
(3)
n+ + n- ) N
(4)
and
Ec(x) ) E(x) - ∆E(2 log Ω(Nx))1/2
(11)
the entropy vanishes, and the approximation above is not valid. The system becomes frozen into a small number of states with low energies. In the thermodynamic limit for E < Ec(x), instead of eq 10, S(x,E) ) 0.7 A similar approach has been used by Bryngelson and Wolynes for the protein-folding problem.9 The model predicts that for every polarization x there is a critical temperature, the polarization-dependent phase transition (Tc(x)), given by
βc(x) )
|
1/2
[2S*(Nx)] 1 ∂S ) ) ∆E Tc(x) ∂E E)Ec(x)
(12)
where S*((Nx) ) log Ω(Nx) is the configuration entropy. In the large N limit, using the Stirling’s approximation
[ (1 +2 x) log(1 +2 x) - (1 -2 x) log(1 -2 x)] (13)
S*(Nx) = N -
At a particular T, such that T ) Tc(x0), for |x| < |x0| the system has a behavior like the standard Born-Marcus model.
7682 J. Phys. Chem., Vol. 100, No. 18, 1996
Leite and Onuchic
Figure 3. Free energy of the OW model as a function of total polarization for three different temperatures, eq 14. The dotted line corresponds to the temperature in which the system is frozen for x > 0.5. The dashed is the high-temperature limit. The continuous line is at the temperature in which the entire system is frozen.
As it hits x0, the dynamics becomes glassy. The average free energy is
{
F(x) ) ∆E2(q) for T > Tc(x) (|x| < |x0|) - TS*(Nx) 2T Ec(x) ) E(x) - ∆E(2S*(Nx))1/2 for T e Tc(x) (|x| > |x0|)
E(x) -
(14)
The free energy depends parametrically on the total polarization, and it is used as a smooth effective potential above Tc(x). There is a competition between the entropic and energetic terms (roughness of the potential). Since the entropic term is a function of the total polarization, there is a polarizationdependent phase transition. At T e Tc(x), the system becomes frozen in a few low-energy states, and flutuations become important. The average free energy of eq 14 is the average result over several local minima. Figure 3 shows several plots for the average F(x), and no drastic behavior is observed at x0. Above Tc(x), the average effective potential F(x) is well described by a quadratic curve, as in the Marcus theory. The question remains regarding the fluctuations of this free energy below the critical temperatures when fluctuations are dominant. The procedure to estimate the entropy, following eq 10, can be used for any physically relevant order parameter or reaction coordinate. In eq 10, the total polarization was used. Thermodynamically, however, there is only one phase transition, here called global phase transition, which is associated with the total number of states (2N). It is characterized by
Sg(E) ) log[2Ng(E)] ) 0
(15)
where Sg(E) is the total entropy and g(E) is the energy probability distribution. Using eqs 6 and 12 1/2
βgc )
1 [2N log 2] ) ∆E Tgc
)
(2 log 2)1/2 ∆
(16)
which is similar to the result for βc(x)0). The critical energy Egc where the entropy of the entire system becomes zero is given by24
Egc
) -N∆(2 log 2)
1/2
(17)
III. Kinetic Rules and Polarization Coordinate Section summary: The kinetic rules allow only a single dipole flip per elementary move, with move acceptance based on a
Monte Carlo procedure. The polarization space is studied assuming no interaction energy between dipoles. The multidimensional master equation can be reduced into a unidimensional diffusion equation for the total polarization. The resulting dynamics corresponds to a particle diffusing in an effective almost quadratic potential. In this limit (noninteracting dipoles), all diffusion paths are equiValent and the unidimensional equation is exact. The mean first passage time (MFPT) to evolve between two polarizations is used as a parameter of comparison between the kinetic and thermodynamic (section II) effective potentials. The thermodynamic effective potential reasonably describes the exact one. The idea that entropic potentials can successfully describe kinetic effective potentials emerges from this model. To address dynamics, the connectivity between states has to be defined, i.e., which configurations are connected by elementary moves. Only a single dipole flip per elementary move is allowed. Dipole motion is therefore asynchronous. This simple scheme is attractive because it is easily realized by numerical simulations. Therefore for a system with N dipoles, each state is connected to N other solvent states. Each state with polarization p is only connected to n- states with polarization p - 2 and n+ with p + 2 (single flips changes polarization by (2). To complete this dynamical scheme there is a need for interconversion rates between states. A reasonable dynamics, which satisfies detailed balance, is the Monte Carlo kinetics, which is equivalent to a rate from a state A to a connected state B of
{
R0 exp[-(EB - EA)/T] for EB > EA R ) RN 0 for EB < EA N
(18)
where R0 is associated with the flipping rate of any dipole. This kinetics is related to our numerical simulations in the following way. Allowed moves are obtained by picking up one dipole at random, and the moves are accepted with a Monte Carlo probability. The time unit (1/R0) is defined as one attempt to change a dipole state. In this section, the limit of noninteracting dipoles is studied (this is equivalent to infinite temperature). All interconversion rates are R0/N. States with the same polarization p are equivalent since they are connected to the same number of neighboring states with polarization p ( 2. As an illustration, a connectivity map for all states is shown for a system of four dipoles (Figure 4). The master equation is easily written, and the transition matrix is (2N) × (2N). We project all the 2N states in a total polarization 1-d coordinate and get an effective unidimensional master equation for the total polarization (see Appendix A). The 1D equation has 2N + 1 polarization states, and all transitions are local. Since the effective potential is a slowly varying function of x, a Smoluchowski equation might be a good approximation for sufficiently large systems. Two tests for the reasonableness of this kinetic approach can be made. First, the equilibrium distributions Ve(x) obtained from thermodynamics and kinetic approaches have to be the same. The equilibrium distribution of the master equation, obtained by recursion (Appendix A), is given by the binomial distribution
Ve(x) )
N! ) exp[-F(x)/T] n+!n-!
(19)
where F(x) ) -TS*(x) is the entropic free energy in the noninteracting dipole limit. The second relevant test is the comparison between the MFPT given by the Smoluchowski equation and the one given by the
Solvent Landscapes in CT Reactions
J. Phys. Chem., Vol. 100, No. 18, 1996 7683 Notice that no separation of time scales is needed to validate the free energy as an effective potential; the equilibrium distribution was obtained from a local master equation. If that was not the case, the system would have to probe the polarization space for each polarization x, keeping x nearly fixed (x always in the interval x ( δx, where δx , x). Using simple arguments, one can show that it is not reasonable. As an extreme example: given any state with zero polarization (n+ ) n-), there is its symmetric state which consists of flipping every single dipole. The agreement between the kinetic and thermodynamic effective potentials is due to the local geometry (or connectivity) that reflects the global entropy of the system. IV. Kinetic ModelsSystem with Interactions
Figure 4. Connectivity map for a system with four dipoles. Each dot corresponds to a state. The states that are converted into each other by a single elementary move (change of a single dipole state) are connected.
full master equation. Starting from the Smoluchowski equation10
[ ( ) ]
∂ 1 ∂ ∂F(x) ∂2 + 2 V(x,t) V(x,t) ) D ∂t kBT ∂x ∂x ∂x
(20)
where V(x,t) is the polarization probability distribution. The MFPT, τ(xf,xi), from an initial state of polarization xi to a final polarization xf is given by11,12
τ(xf,xi) ) ∫x dx∫x dx′ xf
x
i
i
exp[β(F(x) - F(x′))] D
(21)
The discrete version of the equation above is13 xf
τ(xf,xi) )
x
∑∑
x)xi x′)xi
δx2 D
exp[β(F(x) - F(x′))]
(22)
where D ) 2/N2 and δx ) 2/N. τ(xf,xi) can be compared with the exact MFPT using the full master equation (Appendix A). As expected, as the system size increases, the two results become identical. In the small-x limit,6 the Smoluchowski equation is
[
]
2R0 ∂ 2R0 ∂2 ∂ Vt(x,t) ) x + 2 2 Vt(x,t) ∂t N ∂x N ∂x
(23)
where Vt(x,t) is the total polarization time probability distribution, with drift coefficient -2R0x/N and diffusion coefficient 2R0/ N2. Its effective potential corresponds to a particle diffusing in a quadratic effective potential10
(
F(x) ) TN -log 2 +
)
x2 2
(24)
This potential has to be linear with T, since it is fully entropic (note, eq 23 is T-independent). In terms of “routes” or “paths” for the diffusion, every state with the same x feels the same environment. Therefore all paths between two polarizations are equivalent. This justifies the representation of the process in terms of a single coordinate. The concept of equiValent diffusion paths (or simply equivalent paths) is used in subsection IV.c.
Section summary: In this section interactions between dipoles and between dipoles and the charge cavity are included. The energy distribution is described by two parameters, the average energy as a function of x, E(x), and its random fluctuation with width ∆E (associated with the roughness of the landscape). We investigate their contributions separately. The first term gives rise to a system with a simple Arrhenius dependence, subsection IV.a. The latter one produces several temperature-dependent kinetic regimes, including a glassy one. The following kinetic phase transitions are associated with it. (a) The local phase transition corresponds to a temperature where the system starts to get trapped in local minima. It is responsible for a spread in the distribution of escape times from local minima, and it gives rise to the breakdown of the equivalent diffusion path picture, subsection IV.b. Below this temperature, escaping from local minima starts to be dominated by the smallest barrier, disfavoring other routes. (b) The global phase transition is the temperature such that below this temperature a small set of lowenergy states are probed most of the time, and only few are visited in each event. This is associated with the breakdown of ergodicity. Also, different systems with the same average energy distribution have different behaviors; this is called nonself-averaging. Thermodynamically, this is the only phase transition (subsection IV.b). (c) The polarization-dependent phase transition, obtained from the OW model, is associated with the distribution of energies at each polarization, subsection IV.c. The relative strength between the two parameters, (E(x) and ∆E), determines the kind of kinetics that is involved at low temperatures, either following an Arrhenius behavior or presenting various glassy-like regimes. Dipole-dipole or charge-dipole interactions are now included. What happens when these interactions are “turned on”, with energies following the distribution given by eq 7? As the interaction increases (or equivalently, as temperature decreases), trapping in local minima starts to occur. To probe the different dynamical regimes, we vary the two parameters that characterize the energy distribution, the average energy as a function of x, E(x), eq 5, and its random fluctuation with width ∆E, eq 6. There are two limit cases to be studied: (a) ∆E ) 0 and E(x) * 0. There is no disorder, and a constant energy value is associated with each x. (b) ∆E * 0 and E(x) ) 0. There is no correlation between polarization x and energy. States for all x’s obey the same energy distribution. We first discuss case (a). Sections IV.b and IV.c investigate case (b), following two characteristic times of the problem, the mean escape time to move out of any state (independent of being a minimum or not) and the MFPT defined in the previous section. In section IV.b, we estimate the escape times and describe the concepts of global and local phase transitions. Section IV.c discusses how the local phase transition affects
7684 J. Phys. Chem., Vol. 100, No. 18, 1996
Leite and Onuchic
Figure 5. Mean first passage time τ(0.4,0.0) from the equilibrium polarization (xi ) 0) to the final polarization xf ) 0.4, as a function of the inverse temperature (in units of average energy parameter J) for a system with 20 dipoles. The continuous line is the result using the free energy as an effective potential in eq 22. The dashed line is the result using the exact matrix kinetic formalism.
Figure 6. Average fraction f computed as the ratio of the fastest neighboring escape rate and the total escape rate (eqs 30 and31), as a function of inverse of temperature β. Calculations were performed for a system with ∆E ) 201/2 and 20, 190, and 1140 neighboring states. The vertical lines indicate the position where f changes curvature, defined as βcloc.
kinetics, using the equivalence of paths picture. We also elaborate on the consequences of the polarization-dependent phase transition. Section IV.d addresses the general case where ∆E * 0 and E(x) * 0. IV.a. Non-zero Average Energy Interaction with No Disorder (E(x) * 0 and ∆E ) 0). In the limit that E(x) * 0 and ∆E ) 0 (no disorder), all states with the same polarization are equivalent, i.e., all neighboring states with the same polarization have the same energy, as in section III. The multidimensional master equation can be exactly projected into a simple one-dimensional (total polarization) equation (discrete form of a diffusion equation). The matrix elements for the unidimensional master equation are given by eqs A.12 and A.13. Again the MFPT, τ(xf,xi), is calculated using the exact unidimensional matrix formalism and using the free energy F(x) as an effective potential, eq 22, but with a drift coefficient E(x). Figure 5 shows a comparison between the two approaches for τ(0.4,0.0) for a system with 20 dipoles, using E(x) as25
approximated by
E(x) ) NJx2
(25)
The agreement of the free energy approach with the exact solution is remarkably good. Both present a simple exponential behavior. As the interaction Jβ increases, the deviation becomes smaller. At very low temperatures, the system is almost frozen, but nevertheless it is not initial condition dependent. No phase transition exists. IV.b. Local and Global Phase Transitions and Escape Times. Local Phase Transition. This transition corresponds to a temperature where the system starts to get trapped in local loc minima, Tloc cl (1/βcl ). Below this temperature, escapes from local minima go preferentially through the neighbor with the smallest barrier instead of overcoming a typical barrier. We assume that each state has M neighbors (in the simple case considered in section III, for a single dipole flip, M ) N ) number of dipoles). The average number of neighboring states 〈nl(E)〉 for a given state with a energy between E and E + dE is
〈nl(E)〉 )
M
x2π∆E
[
exp -
]
E2 dE 2∆E2
(26)
As in the section II, if 〈nl(E)〉 . 1, a local entropy can be
Sl(E) = log〈nl(E)〉 = log M -
E2 2∆E2
(27)
There is a critical energy Eloc c , where the entropy becomes of ) -∆E(2 log M)1/2)) and order zero (Eloc c
βloc c )
|
∂Sl (2 log M)1/2 2 log M 1/2 1 1 ) ) ) loc ∂E E)Elocc ∆E N ∆ Tc (28)
(
)
the system spends most of the time at local For β g βloc c minima. The local phase transition is determined only by the local distribution of energies.26 Tloc c is what was defined by Bryngelson and Wolynes as Thigh (see eq 29 of ref 12). They high) is have pointed out that one of the consequences of Tloc c (T that below this temperature, the rate distribution is flattened. The rate Ri for the system to leave a state i is given by the sum of the individual rates to all its neighbors M
Ri ) ∑rij
(29)
j)1
where rij is the individual escape rate from i to j, eq 18. To investigate the concept of “preferential escape” below the local phase transition, we compute the ratio fi between the fastest () max(rij)) and the total rate Ri individual rate rmax i
fi ) rimax/Ri
(30)
The ratio fi shows the importance of the fastest rate (smallest barrier) on the overall escape rate. When the fastest rate starts to dominates, a single “kinetic” route becomes preferential. This is similar to the ergodicity breaking observed for spin glasses14 and random heteropolymers.15 To probe this transition, we measure the average value of fi
f ) ∑ fi exp[-βEi]/∑ exp[-βEi] i
(31)
i
Figure 6 shows f as a function of the inverse of temperature β. f was calculated for three different systems, with ∆E ) 201/2,
Solvent Landscapes in CT Reactions
J. Phys. Chem., Vol. 100, No. 18, 1996 7685 and the average escape time 〈tesc〉 is simply its inverse
TABLE 1: Local Phase Transition Results M neighbor states
theory βcloc
simulation βcloc
20 190 1140
0.55 0.72 0.84
0.46 0.68 0.80
having M ) 20, 190, and 1140 (it corresponds to different connectivity rules of flip of one, two, or three dipoles simultaneously, respectively; note that only for f calculations, it is used as a connectivity different than the standard single-flip one). We average f over 20 system’s configurations. For β ) 0, f ) 1/M, since all individual rates are the same. As β increases, f increases until it saturates. Since there is no sharp transition for a finite system, there is not a unique way to define the transition temperature. We define βloc c as the inflection point for f. Table 1 shows a comparison between the phase transition temperatures obtained through simulation and the thermodynamic limit given by eq 28. Even for such a small system, our kinetic definition is in qualitative agreement with eq 28. At βloc c the fastest individual rate starts to dominate the escape rate. Escape Times and Global Phase Transition. An analytical estimate of the escape time tesc is obtained using an approach similar to that of Bryngelson and Wolynes for the protein-folding problem.12 At a temperature T, the occupancy of states with energy E is
〈tesc〉 )
1 〈R〉
One might note that these calculations are done on the inverse time coordinate (“rate domain”). The same calculations can be done equivalently on the direct time coordinate (“time domain”). The average escape time tesc(E) from states with energy E is the inverse of its escape rate R(E), tesc(E) ) 1/R(E). Since low (high) energy states take a longer (shorter) time to escape, the probability distribution of visiting a state with energy E, h(E) (which is proportional to the number of escapes from states with energy E) is
h(E) ) Khg(E)/tesc(E)
g(E) ∝ e
g(E)
(32)
where g(E) is the probability of having a state with energy E and is given by eq 7 with E(x) ) 0. Since E and x are not correlated, there is no polarization dependence: g(E,x) ) g(E). The normalized g(E) is
{
}
(33)
K ) [∫-Ecg g(E) exp(-βE) dE]-1
(34)
g(E) ) K exp -
[E + β∆E2]2 2∆E2
where K is Eg
c
Kh )
[∫
Egc
-Egc
g(E) dE tesc(E)
]
-1
) [∫-Ecg g(E) R(E) dE]-1 ) 〈R〉-1 Eg
[
∞
]
(E′ - E) + T
R0N∫-∞ dE′ g(E′) (35) E
Substituting g(E), we obtain
(
R(E) ) R0Nζ
) [
]
E + (1/2)β∆E2 E + β∆E2 + exp T 21/2∆E
(
R0Nζ -
)
E (36) 2 ∆E 1/2
where
ζ(x) ≡ 1/2[1 - sgn(x) erf(|x|)]
(37)
The average escape rate is given by
〈R〉 ) ∫-Ecg g(E) R(E) dE Eg
c
(41)
c
The average escape time 〈tesc〉 is given by
〈tesc〉 ) ∫-Ecg dE h(E) tesc(E)
(42)
c
Substituting eqs 40 and 41 into eq 42, we obtain eq 39 as well. Both formalisms are valid, but since our simulations are done on time domain, we use the later approach for the numerical work. As β approaches βgc , the global phase transition, the escape time distribution is flat. To probe the dynamics around βgc , we compute the mean and the standard deviation of the log of escape times, 〈log tesc〉 and ∆log tesc ) [〈log2 tesc〉 - 〈log tesc〉2]1/2. The average of the log of escape times 〈log tesc〉 can be estimated by
〈log tesc〉 ) ∫-Ecg dE h(E) log tesc(E) Eg
(43)
c
The average escape rate for a state with energy E is (using eqs 18 and 29)
R(E) ) R0N∫E dE′ g(E′) exp -
(40)
where Kh is a normalization constan
Eg
-βE
(39)
(38)
A more rigorous calculation of 〈log tesc〉 is shown in Appendix B. 〈tesc〉 and 〈log tesc〉 are compared with simulation results as a function of β. Figure 7 shows average escape times for 10 different runs, each with 104 escapes for a system with 20 dipoles. Even though the number of escapes (states) probed in this simulation is much smaller than the total number of states (220-106), it is much larger than any number of escapes involved in the MFPT (∼102). Runs differ only by the initial conditions. For β < βloc c , 〈tesc〉 is in good agreement with this simple theory. For β > βloc c , fluctuations are large and the results vary depending on initial conditions, providing evidence for the local phase transition. Nevertheless, the theory still gives the right asymptotic behavior at large times. Because βgc > βloc c , fluctuations due to the local phase transition occur before the system reaches βgc , and one cannot use fluctuations on 〈tesc〉 to probe the global phase transition. Figure 7b,c shows 〈log tesc〉 and ∆log tesc. There is a strong dependence on the initial conditions for both 〈log tesc〉 and ∆log tesc if β∆ g βgc ∆ ) 1.1. This is in agreement with the REM expected phase transition temperature βgc ∆ ) (2 log 2)1/2 ) 1.18. This is a characteristic of a glass phase transition. (Differently than for 〈tesc〉, 〈log tesc〉 does not give the correct
7686 J. Phys. Chem., Vol. 100, No. 18, 1996
Leite and Onuchic
Figure 8. System with all states connected, M ) 2N. Average escape time 〈tesc〉 as a function of inverse of temperature (in units of disorder energy parameter ∆) β∆, for N ) 20. This result is scaled with N, but the transition temperature βcg is N independent. The saturation value is around the total number of states 2N ()106), as expected. The continuous line is the theoretical result, eq 42. Using the change of curvature point as parameter for the inverse of critical temperature, we obtain βcg ) 1.0. The symbols are the simulation results, average escape times over 104 escapes for six different systems. For β > βcg the results are system dependent.
Figure 7. Average escape times for a system with 20 dipoles over 104 escapes as a function of the inverse of temperature (in units of disorder energy parameter ∆). Ten runs were performed with different annelead initial conditions using a kinetic rule where a single dipole flip is allowed (M ) 20). (a) Average escape time 〈tesc〉. (b) Average of the log of escape time 〈log tesc〉. (c) Standard deviation of the log of the escape time ∆log tesc ) [〈log2 tesc〉 - 〈log tesc〉2]1/2. O in a and b is the theoretical estimate using eqs 42 and 43, respectively. ] in b is the rigorous calculation using eq B.2. Since a large number of escapes are averaged over in each run, the numerical error bars of the means are irrelevant (this is true for the remaining Monte Carlo simulations). Notice that in the nonglass regime the results are essentially the same.
asymptoptic behavior at large times and low temperature; see Appendix B). Of the two transitions, the local phase transition is the weak one, and its effects can be minimized by long simulation due to averaging over multiple escapes, particularly when we measure 〈log tesc〉. During the time evolution of the simulation,
multiple minima are visited, and this averaging takes place. For the global phase transition (the actual phase transition), there is a break of ergodicity. It depends on global properties of the system, and its effects are more severe. For temperatures below Tgc , fluctuations due to initial conditions or different samples are large compared to the average values. From the analytical results it is observed that 〈tesc〉 converges to a saturation value equal to the total number of states of the system. An interesting case to illustrate the saturation behavior of 〈tesc〉 is when the number of neighboring states M is equal to the total number of states 2N (all states are connected to each other, a totally connected system). The dynamics resembles the globally connected system that has been studied by Saven et al.16 In this limit βgc ) βloc c . 〈tesc〉 as a function of β has a sigmoid-like behavior. Around the critical temperature, 〈tesc〉 changes its curvature similarly to f (eq 31 and Figure 6). At βgc , however, it is still much smaller than 2N (Figure 8). IV.c. Equivalence of Paths and Polarization-Dependent Effects. In this subsection we focus on time evolution of the polarization coordinate. All results are from numerical simulations. Since we are interested on the evolution between two different polarizations, the first passage time (FPT) is the characteristic time to describe the events in the total polarization space. All dynamical information can be obtained from the distribution of FPT. We define the stochastic variable ω(xf,xi) as the number of escapes during the evolution between xi to xf (the equilibrium polarization is taken as xi ) 0). Notice that the number of escapes is different than the number of steps (time), since only successful attempts are considered. For a free random walk, the FPT is given by ωrw(xf,xi), and every step is “successful”. The random walk MFPT is simply the average 〈ωrw(xf,xi)〉. As the temperature is lowered, since it takes longer for escapes to occur, the MFPT increases. For β < βloc c , the paths are equivalent, and the MFPT is simply the product between the average number of escapes and the average escape time
τ(xf,xi) ) 〈ωrw(xf,xi)〉〈tesc〉
(44)
This result is equivalent to eq 22, if D ) 2/N2〈tesc〉. As β g βloc c , the system escapes preferentially through some particular
Solvent Landscapes in CT Reactions
J. Phys. Chem., Vol. 100, No. 18, 1996 7687 For a system evolving from polarization xi to xf, the average change in energy is
δE(xf,xi) ) E(xf) - E(xi) ) NJ[xf2 - xi2]
(46)
This is the potential energy associated with the reorganization energy of the polarization coordinate. Until the glassy regime is reached, the local barriers due to the roughness are of the order of ∆E2/T.12,18 There are three different regimes, according to the relative values of roughness and change in potential energy. We define γ that measures the ratio between them
γ) Figure 9. ∆Ω(xf,β), (eq 45), is plotted as a function of the inverse temperature (in units of disorder energy parameter ∆). Averages of the log of number of states visited 〈log ω(xf,β)〉 are subtracted (normalized) by the random walk result 〈log ωrw(xf)〉. Therefore, the different final polarization xf results are easily compared.
neighbors. This leads to large fluctuations, the paths are not equivalent, and τ(xf,xi) deviates from eq 44. To probe the equivalent path idea, ω(xf,xi) is investigated as a function of temperature, ω(xf,xi,β). If the diffusion had all paths equivalent, ω(xf,xi,β) would be temperature independent. Since the distribution of ω(xf,xi,β) is flat, 〈log ω(xf,xi,β)〉 is calculated instead of 〈ω(xf,xi,β)〉. To quantify the deviation from the random walk result, we used
∆Ω(xf,xi,β) ) 〈log ω(xf,xi,β)〉 - 〈log ωrw(xf,xi)〉 (45) The results obtained are shown in Figure 9. If all paths were equivalent, ∆Ω ) 0. Again, the transition is smooth, and there is no unique way to define βloc c . The qualitative transition from equivalent to nonequivalent diffusion paths can be inferred. Notice that due to trapping ∆Ω is always a positive number. Also, a small dependence on the value of xf is observed. According to the thermodynamic result of the OW model, the polarization-dependent phase transition should occur for larger values of xf. This is in qualitative agreement with our data. Nevertheless, we do not have sufficient data to fully explore the polarization-dependent phase transition. Our approach using FPT to probe the phase transition is computing intensive, particularly close to βc(xf) (βc(xf)∆ ∼ 1.). βc(x) is particularly < βc(x). The polarizationhard to observe because βloc c dependent phase transition becomes blurred by the local phase transition. In other systems, with different local energy distribu> βc(x), and βc(x) can be easily tions, one can have βloc c observed.17 IV.d. General Case (E(x) * 0 and ∆E * 0). The general case has both a nonzero average energy (E(x) * 0) and disorder (∆E * 0). Any realistic system has both of these nonzero. Each of them contributes differently to the MFPT. E(x) given by eq 25 is considered. As discussed, ∆E modulates the escape rates (diffusion coefficient), and E(x) adds a quadratic potential that the system has to overcome when diffusing from xi to xf. In the previous subsections, the two limit cases were studied. Case (a) presents a well-behaved system with no phase transitions, section IV.a. Case (b) presents phase transitions with different dynamical regimes, sections IV.b and IV.c. In this subsection we propose a parameter (γ(xf)) to measure the competition between cases a and b.
δE(xf,xi) ∆E2/T
(47)
Three distinct dynamical behaviors are observed depending on the value of γ. (i) For γ , 1, the distribution of energies g(x,E), eq 7, is slightly changed in a single dipole move. In other words, g(x,E) = g(x + δx,E) and is polarization independent. g(E) is a good approximation for the energy distribution of states connected to a particular state Nx (see eq 32). Since ∆E does not depend on x, the diffusion coefficient is also independent of x. The MFPT is given by xf
τ(xf,xi) )
x
∑∑
x)xi x′)xi
δx2 D
exp[β(F(x) - F(x′))]
(22)
where D ) 2/N2〈tesc〉. This is the limit discussed in section IV.b. The only difference is that the free energy has also a drift coefficient E(x), but this is a minor effect. One expects to g observe the phase transitions βloc c , βc(x), and βc . (ii) For γ . 1, the other extreme, the disorder is very small compared to the reorganization potential energy. τ(xf,xi) is approximated by eq 22, with diffusion coefficient D ) 2/N2, as in section IV.a. In this limit, phase transitions are not expected. Cases (i) and (ii) provide upper and lower bounds for τ(xf,xi), respectively. (iii) For γ ∼ 1, there is no simple analytical form for τ(xf,xi), and no simple expression for the diffusion constant is available. The phase transitions βloc c and βc(x) are not precisely defined, nevertheless βgc is the same. Using β∆E ) 1.79, xi ) 0.0, and xf ) 0.4, simulations were performed with varying βJ (Figure 10). For βJ < 1.0 (γ < 1.0), the system falls in case (i). As βJ increases, the MFPT starts to move toward case (ii). For βJ ) 3.4 (γ ) 3.4), the system already has a behavior similar to that of regime (ii). In the intermediate regime (iii), γ ∼ 1 (βJ ∼ 2), τ(xf,xi) is bounded by the two results. For all limits, as β g βgc the system is frozen into local minima states, and path equivalence is not valid. The single reaction coordinate description breaks down. V. Conclusion Solvent polarization plays a major role controlling chargetransfer reactions. Exploring a simplified model, we have shown how kinetically and thermodynamically a complex solvent landscape can be reduced into a global one-dimensional reaction coordinate. Such a reduction, which validates the Marcus theory for electron transfer, is possible as long as equivalent diffusion paths for the solvent dynamics exist during the solvent dynamical evolution in the polarization coordinate. As the system evolves between two polarization values, it may do so by different paths that, at high enough temperatures, are kinetically
7688 J. Phys. Chem., Vol. 100, No. 18, 1996
Figure 10. Mean first passage time τ(0.4,0.0) from equilibrium polarization (xi ) 0) to a final polarization xf ) 0.4, as a function of the inverse temperature (in units of average energy parameter J), for a system with 20 dipoles. The continuous and dot-dashed lines are the results using the free energy as an effective potential for β∆E ) 0.0 and 1.8, respectively (cases (ii) and (i), respectively). ] is the simulation result for β∆E ) 1.8.
equivalent. In this regime, this “complex solvent” shows a Debye--like behavior, and Marcus theory is recovered. As the temperature is lowered the equivalent path description breaks down. The kinetic behavior depends on (i) the average energy and the density of states for each polarization and (ii) the energetic roughness of the landscape. The first term determines the average polarization-dependent free energy, which is the effective potential in Marcus theory. The roughness controls the diffusion constant. Two characteristic times are used to describe the dynamics. One is the escape time, which is the average time to leave a state (to make a successful move). The other one is the mean first passage time, τ(xf,xi), between polarizations xi and xf. As the roughness is increased (or temperature reduced), three different critical temperatures are reached. Below the local phase transition, trapping in local minima occurs, and when escaping from them, the system starts to favor particular routes. The concept of equivalent paths starts to break down. Below this temperature, fewer paths dominate the τ(xf,xi), and a deviation from Debye-like solvent behavior is observed. The second critical temperature is the polarization-dependent phase transition, and it is associated with a different density of states as a function of the total polarization. The last critical temperature is the global phase transition, which occurs when the system is frozen in a few very energy-deep states. This is associated with breakdown of ergodicity and is the only real thermodynamic phase transition. At this point the escape times are strongly dependent on the initial condition, and a strong failure of the single reaction coordinate picture is observed. In this limit, different samples show large deviations in their kinetic properties. The equivalent path idea is associated with a regime where the dynamics can be described as simple diffusion along a single effective reaction coordinate. All roughness effects are included in the diffusion constant. In the glassy regime, paths are not equivalent, and this simple diffusion picture breaks down. Few paths dominate the diffusion process, and nonexponential behavior as a function of the temperature is expected. This work gives insight into how the free energy can be used as an effective potential for solvent dynamics. Even though real molecules have more complicated dynamics than the one proposed here, we expect the basic features to be preserved. The random energy term associated with the landscape roughness may be composed of several contributions when extended
Leite and Onuchic to real solvents. Geometrical disorder, fluctuations in dipole strength, and hydrogen bonding are some of the possible examples. Since the dynamics is strongly dependent on the roughness of the system, one should investigate these effects in real solvents. One expects strong hydrogen bonding to give rise to large values of roughness. Such solvents (like alcohols) should present the glassy behavior. A possible experimental evidence might be the low-temperature results for alcohols, which deviate dramatically from the equilibrium predictions.19 The biggest pitfall of this theoretical model is that it does not include energy correlation between states, and that is clearly not realistic. Nevertheless, the general features are expected to hold for real solvents. One does expect to observe the limiting cases: (1) the high-temperature limit with Marcus-like behavior, and (2) the low-temperature glassy behavior. The techniques applied for this model can be extended for systems where there are energy correlations between states.17 It is interesting to notice that some of the features observed for this solvent model are expected for other problems that also involve complex landscapes. For example, in section III, it was shown that entropic potentials can successfully describe kinetic effective potentials. This general result explains why some entropic potentials are successfully used to describe complex systems (rubber-band-like potentials).20 The transition from multiple to single paths and the different kinetic behaviors may be generalized to other complex problems. The protein folding problem is one of the cases where such connection may be possible.21 The connectivity map for protein folding is very complicated and not trivially obtained as here, but the equivalent paths description can be related to the funnel concept.22 The solvent problem is highly simplified because of its connectivity, and its description in terms of a simple reaction coordinate is a consequence of this property. Acknowledgment. We thank Peter Wolynes, Bruno Zimm, and Phil Kiefer for the helpful discussions and suggestions. We also thank the referee for the careful comments and suggestions that substantially improved this paper. J.N.O. is a Beckman Young Investigator. V.B.P.L. was partially supported by the Brazilian Agency Conselho Nacional de Desenvolvimento Cientı´fico e Tecnolo´gico, CNPq. This work was funded by the National Science Foundation (Grant No. MCB-9316186). Appendix A: Derivation of the Effective Unidimensional Master Equation Given a system with N dipoles (2N states), the master equation that describes the time evolution of this system is given by
∂ b V (t) ) WV b(t) ∂t
(A.1)
where b V(t) is the time probability distribution. b V(t) is a 2N component vector, and W is the transition matrix. We can rewrite eq A.1 as
∂ ∂t
Vn(t) ) ∑ [Wn,n′Vn′(t) - Wn′,nVn(t)]
(A.2)
n′
We first consider the dipoles without interaction (this corresponds to study of the system at infinite temperature). In this case all states are equivalent, and the transition matrix can be written as
Wn,n′ )
{
R0/N if n and n′ are connected 0 otherwise
(A.3)
where R0/N is the flipping rate of a dipole. To relate the results
Solvent Landscapes in CT Reactions obtained here with the numerical simulations, the time unit (1/ R0) is defined as one attempt to change a dipole state. Given a state n with polarization x0, one can rewrite eq A.2, separating the contribution of states of polarization x0 - (2/N) and x0 + (2/N)
∂ ∂t
R0
∑ Vn′,x′(t) δ(dn,n′ - 1) δ(x0 - (x′ - δx)) + N n′ R0 N
Vn′,x′(t) δ(dn,n′ - 1) δ(x0 - (x′ + δx)) ∑ n′
(A.4)
where δ(dn,n′ - 1) restricts the sum to states that are connected to state n (“state distance” dn,n′ ) 1). δ(x0 - (x′ ( δx)) restricts the sum to states with polarization x0 ( δx (δx ) 2/N). The first and second sums correspond to contributions with polarizations x - δx and x + δx, respectively. Finally, we can sum eq A.4 over all states with same polarization x0
δ(x - x0)[eq A.4] ∑ n′′
(A.5)
which gives
∂ ∂t
∑n Vn,x (t) δ(x - x0) ) -R0 ∑n Vn,x (t) δ (x - x0) + 0
R0 N R0 N
0
Vn′,x′(t) δ(dn,n′ - 1) δ(x0 - (x′ - δx)) δ(x - x0) + ∑ n′,n
Vn′,x′(t) δ(dn,n′ - 1) δ(x0 - (x′ + δx)) δ(x - x0) ∑ n′,n
(A.6)
If we define the total polarization time probability distribution of the system as
Vt(x0,t) ) ∑ Vn,x(t) δ(x - x0)
(A.7)
(
Wt is the transition matrix of the form
Wt )
R0 N
Vn,x0(t) ) -R0Vn,x0(t) +
J. Phys. Chem., Vol. 100, No. 18, 1996 7689
-N N
1 -N N-1
2 -N N-2
0 3 -N
· · ·
0
·
·
·
· · -N · 1
N -N
)
(A.10)
The equilibrium polarization probability distribution b ve(n+) can be obtained as
Wtb V e(n+) ) 0
(A.11)
and the solution is
Ve(n+) )
N N ) n !n n+!(N - n+)! + -!
as in eq 9. The Smoluchowski or master equation (eq 20 or A.8) can be used to calculate τ(xi,xf). In this work we use the master equation in its matrix form, so we can calculate τ(xi,xf) without any approximation. τ(xi,xf) can be calculated for any xi and xf using the method described in ref 23. When dipole interactions are included, one substitutes R0 by 〈R〉x,x′, that it is polarization dependent. As discussed in section IV, there are two limiting cases. For ∆E ) 0 and E(x) * 0, each term on the right-hand side of eq A.8 changes to include drift and R0 ) 1. Equation A.8 becomes
1 + x + δx ∂ Vt(x0 + V (x ,t) ) -s0(x0) Vt(x0,t) + s+(x0) ∂t t 0 2 1 - x + δx Vt(x0 - δx,t) (A.12) δx,t) + s-(x0) 2
(
(
)
)
n
where after some manipulation of eq A.6 we obtain
[
(n+ + 1) ∂ Vt(x0,t) ) R0 - Vt(x0,t) + Vt(x0 + δx,t) + ∂t N (n- + 1) Vt(x0 - δx,t) N
]
or
∂ 1 + x + δx Vt(x0 + δx,t) + V (x ,t) ) R0 - Vt(x0,t) + ∂t t 0 2 1 - x + δx Vt(x0 - δx,t) (A.8) 2
[
)
(
(
)
]
where n+ and n- are related to x0 by Nx0 ) n+ - n-. Equation A.8 can be written in a matrix formalism like eq A.1
∂ V t(t) b V (t) ) Wtb ∂t t
(A.9)
where b Vt(t) is the time probability distribution of the total polarization of the system (a 2N + 1 component vector) and
s((x0) )
{
exp[-β(E(x0) - E(x0 ( δx))] for E(x0) < E(x0 ( δx) for E(x0) > E(x0 ( δx) (A.13)
1
and s0(x0) ) s+(x0) + s-(x0) in order to satisfy detailed balance. For ∆E * 0 and E(x) ) 0, R0 is substituted by 1/〈tesc〉 (calculated in subsection IV.b). Appendix B Technical remarks on the averages of the log of escape times are considered in this appendix. A more rigorous calculation of 〈log tesc〉 has to take into account that
log tesc(E) * log tesc(E)
(B.1)
The correct expression for 〈log tesc〉 is given by
〈log tesc〉 ) ∫-Ecg dE h(E)log tesc(E) Eg
c
(B.2)
7690 J. Phys. Chem., Vol. 100, No. 18, 1996
Leite and Onuchic
where ∞
log tesc(E) ) ∑ pr(t,tesc(E)) log t
(B.3)
t)1
where pr(t,tesc(E)) is the distribution of escape times with average escape time tesc(E) (recall t is a discrete variable). The escape rate of a state with energy E is 1/tesc(E). The probability of escaping at time t is equal the probability that it has not escaped before ((1 - 1/tesc(E))t-1) times the probability it escapes at time t, so
pr(t,tesc(E)) )
1
(
1-
tesc(E)
1
)
t-1
(B.4)
tesc(E)
Equation B.3 can be approximated by a continuous integral, so
log tesc(E) ) ∫1 pr(t,tesc(E)) log t dt ∞
(B.5)
which can be evaluated analytically
log tesc(E) )
[
(
- (tesc(E) - 1) log 1 -
1
)] [ ( -1
tesc(E)
Ei log 1 -
1
)]
tesc(E)
(B.6) ∞ where Ei is the exponential integral function (Ei(z) ) - ∫-z exp(-t) dt/t). As can be seen in Figure 7b, eq B.2 represents a minor correction compared to 〈log tesc〉 given by eq 43. For β < βloc c , most of escapes happens through a typical barrier, tesc(E) is a good estimate for the escape time at energy E, and eq B.6 gives the correct value for 〈log tesc〉. For β > βloc c , the energy distribution around a given state becomes important for the escape rates. The distribution of escape times at energy E is broad and system dependent (glassy). Even though its average value tesc(E) may be correct, log tesc(E) will be overstimated because
log tesc(E) >
1
M
∑ log ti(E) M i)1
(B.7)
with M being the number of of escape events simulated. That explains why theoretical 〈log tesc〉 overestimates numerical simulation results, and it gives another indication for the existence of βloc c . References and Notes (1) (a) Ulstrup, J. Charge Transfer Process in Condensed Media; Springer: Berlin, New York, 1979. (b) Onuchic, J. N.; Beratan, D. N.;
Hopfield, J. J. J. Phys. Chem. 1986, 90, 3707. (c) Photoinduced Electron Tranfer; Fox, M. A., Cannon, M., Eds.; Elsevier: Amsterdam, 1988; Parts A-D. (d) Mikkelsen, K. V.; Ratner, M. A. Chem. ReV. 1987, 87, 113. (e) De Vault, D. Quantum Mechanical Tunneling in Biological Systems; Cambridge University: New York, 1984. (f) Newton, M. D. Chem. ReV. 1991, 91, 767. (2) (a) Weaver, M. Chem. ReV. 1992, 92, 463-480. (b) Heitele, H. Angew. Chem., Int. Ed. Engl. 1993, 32, 359-377. (c) Maroncelli, M.; McInnis, J.; Fleming, G. R. Science 1989, 243, 1674. (d) Simon, J. D. Acc. Chem. Res. 1988, 21, 128. (e) Barbara, P. F.; Jarzeba, W. AdV. Photochem. 1990, 15, 1. (f) Bacgchi, B. Annu. ReV. Phys. Chem. 1989, 40, 115. (3) (a) Marcus, R. A.; Sutin, N. Biochim. Biophys. Acta 1985, 811, 265. (b) Marcus, R. A. J. Chem. Phys. 1956, 24, 966. (4) Borgis, D.; Hynes, J. T. Chem. Phys. 1993, 170, 315. (5) (a) Stillinger, F. H.; Weber, T. Science 1984, 225, 983. (b) Daggett, V.; Levitt, M. Annu. ReV. Biophys. Biomol. Struct. 1993, 22, 353. (c) Elber, R.; Karplus, M. Science 1987, 235, 318. (6) Onuchic, J. N.; Wolynes, P. G. J. Chem. Phys. 1993, 98 (3), 2218. (7) Derrida, B. Phys. ReV B 1981, 24, 2613. (8) Derrida, B. Phys. ReV. Lett. 1980, 45, 79. (9) Bryngelson, J. D.; Wolynes, P. G. Proc. Natl. Acad. Sci. U.S.A. 1987, 84, 7524. (10) van Kampen, N. G. Stochastic Processes in Chemistry and Physics; North Holland Publishing: Amsterdam, 1981. (11) Ha¨nggi, P.; Talkner, P.; Borkovec, M. ReV. Mod. Phys. 1990, 62, 251. (12) Bryngelson, J. D.; Wolynes, P. G. J. Phys. Chem. 1989, 93, 6902. (13) Weiss, G. H. J. Stat. Phys. 1981, 24, 587-594. (14) Mezard, M.; Parisi, G.; Virasoro, M. A. Spin glass theory and beyond; World Scientific: Singapore, 1987. (15) Shakhnovich, E.; Gutin, A. J. Chem. Phys. 1989, 93, 5967. (16) Saven, J. G.; Wang, J.; Wolynes, P. G. J. Chem. Phys. 1994, 101, 11037. (17) Leite, V. B. P.; Onuchic, J. N. To be published. (18) Bryngelson, J. D.; Wolynes, P. G. Biopolymers 1990, 30, 177188. (19) (a) Knoshita, S.; Nishi, N. J. Chem. Phys. 1988, 89, 6612. (b) Maroncelli, M. J. Mol. Liq. 1993, 57, 1-37. (20) Kubo, R. Statistical mechanics; North Holland Publishing: Amsterdam, 1965. (21) (a) Onuchic, J. O.; Wolynes, P. G.; Luthey-Schulten, Z.; Socci, N. D. Proc. Natl. Acad. Sci. U.S.A. 1995, 92, 3626-3630. (b) Bryngelson, J. D.; Onuchic, J. N.; Socci, N. D.; Wolynes, P. G. Proteins: Struct., Funct., Genet. 1995, 21, 167-195. (22) Leopold, P. E.; Montal, M.; Onuchic, J. N. Proc. Natl. Acad. Sci. U.S.A. 1992, 89, 8721. (23) Leopold, P. E.; Shakhnovich, E. I. Proceedings of the 26th Annual Hawaii International Conference on System Sciences: Architecture and Biotechnology; Hunter, L., Ed.; IEEE Computer Society Press: Los Alamitos, CA, pp 726-735; Vol. 1. (24) Finite corrections on Ecg take into account the nonexponential terms that make the entropy vanish on eq 17; the first-order correction for the REM is given by (see eq 33 of ref 7) Ecg ) -N∆[(2 log 2)1/2 - (log N/(N(2 log 2)1/2))]. In our rate calculations, this value is used as the critical energy. (25) If one studies the general energy E(x) case, given by eq 5. The only change is in the equilibrium polarization xeq. For ξ(q) * 0, xeq is nonzero and temperature dependent. (26) The local distribution might be different than ∆E of eq 6. For example, if locally the energies are correlated and do not follow g(E) given by eq 7, then βcloc is not given by eq 28, and it might be larger than βcg. We return to this point when we discuss real solvents in the conclusion.
JP951316+
