Interfacial bond-breaking electron transfer in mixed water-ethylene glycol solutions: reorganization energy and interplay between different solvent modes.

We explore solvent dynamics effects in interfacial bond breaking electron transfer in terms of a multimode approach and make an attempt to interpret challenging recent experimental results (the nonmonotonous behavior of the rate constant of electroreduction of S2O8(2-) from mixed water-EG solutions when increasing the EG fraction; see Zagrebin, P.A. et al. J. Phys. Chem. B 2010, 114, 311). The exact expansion of the solvent correlation function (calculated using experimental dielectric spectra) in a series predicts the splitting of solvent coordinate in three independent modes characterized by different relaxation times. This makes it possible to construct a 5D free-energy surface along three solvent coordinates and one intramolecular degree of freedom describing first electron transfer at the reduction of a peroxodisulphate anion. Classical molecular dynamics simulations were performed to study the solvation of a peroxodisulphate anion (S2O8(2-)) in oxidized and reduced states in pure water and ethylene glycol (EG) as well as mixed H2O-EG solutions. The solvent reorganization energy of the first electron-transfer step at the reduction of S2O8(2-) was calculated for several compositions of the mixed solution. This quantity was found to be significantly asymmetric. (The reorganization energies of reduction and oxidation differ from each other.) The averaged reorganization energy slightly increases with increasing the EG content in solution. This finding clearly indicates that for the reaction under study the static solvent effect no longer competes with solvent dynamics. Brownian dynamics simulations were performed to calculate the electron-transfer rate constants as a function of the solvent composition. The results of the simulations explain the experimental data, at least qualitatively.

Introduction

The electrochemical reduction of a peroxodisulphate anion (S2O82– ) is a well-known reaction, although some details of its mechanism are not properly understood so far. This is a good example of bond-breaking electron transfer; the reaction takes place in the vicinity of activationless discharge and demonstrates a “polarization pit” on the current–voltage curves as well as cation catalysis (see works[1,2] and refs therein). Such qualitatively interesting features of this electrochemical redox process make it attractive to employ modern quantum mechanical theory of charge transfer in condensed media.[3,4] This has been done in works[1,2] with the help of quantum chemical modeling. The reduction of S2O82– on a mercury electrode from solutions with variable viscosity (water–sugar and water–ethylene glycol (EG) mixtures) was investigated experimentally as well.[5,6] The current–solution viscosity dependences built for overvoltage values in the region of the “polarization pit” reveal a challenging nonmonotonous behavior: the current decreases first reaching a plateau. However, it begins to increase from a certain viscosity value.[6] As the reduction of S2O82– is accompanied by the −O–O– bond break, the intramolecular reorganization plays a crucial role. Another important driving force of this electron-transfer process is solvent (outer sphere, or environmental) reorganization. One needs, therefore, at least a 3D reaction free-energy surface (along solvent and intramolecular coordinates) to calculate the activation barrier. The experimental data obtained for water–sugar mixtures[5] were interpreted in ref (7) by considering the interplay between two effects of different nature: the first one originates from a pure solvent dynamics (i.e., slow diffusion along the solvent coordinate) and leads to a decrease in the rate constant when the solvent viscosity increases. In contrast, the second effect is solely of static nature and results from the decrease in solvent reorganization energy with increasing viscosity due to the change of dielectric properties of mixed solutions.

The “static” effect was demonstrated to play a crucial role in describing the experimentally observable effects[7] on the basis of the Sumi–Marcus model,[8] where only one solvent mode is addressed. Therefore, it would be tempting to extend such an approach to elucidate the results obtained in ref (6) for water–EG mixtures as well. It was formerly assumed[6] on the basis of an analysis of the solvent correlation function that the experimental effect might result from the interplay of dynamical properties of different solvent modes. The main aim of this work is to check this hypothesis performing Brownian dynamics simulations.

A preliminary analysis of the Pekar factor (C = 1/ε∞–1/ε0, where ε∞ and ε0 are the optical (fast) and static (slow) medium dielectric constants, respectively) versus EG content in mixed solutions points to a very slight static effect.[6] Because solvent reorganization energy is an important parameter in the above-mentioned calculations, it is therefore worth it to go beyond simplified continuum approaches and to describe a possible “static” solvent effect in water–EG solutions at the reduction of S2O82– at atomistic level. The dependence of solvent reorganization on the composition of various mixed solvents in the regime of homogeneous ET mixed solvents was thoroughly investigated for a number of redox couples in works.[9−13] This quantity was extracted from experimental rate constants as well as estimated using model calculations. The authors put primary attention on the analysis of the Pekar factor as a function of the solution composition (this dependence can be both linear and nonlinear) and concluded that preferential solvation crucially affects the solvent reorganization.

The paper is organized as follows: In Section 2, pertinent model and computational details are reported. The solvation energies of S2O82–, the structure of nearest solvation shell of reactant, and the reorganization energies are discussed in Section 3. The model rate constants calculated for different water–EG mixtures are presented in this Section as well. Some concluding remarks can be found in Section 4.

Computational Details

Solvent Correlation Function and Brownian Dynamics

In our work a correlation function M(τ) = <εA(τ)εA(0)> describing the fluctuating reactant energy level in solution ε(τ) plays an important role. This function can be calculated directly from molecular dynamics simulations. However, following a more general and flexible way we will employ another formalism:[14]where λs is solvent reorganization energy and Q(τ) is a function of the complex dielectric spectrum ε(ω):with ε0 being the static dielectric constant, ε0 = ε(ω → 0). (The frequency-dependent function 1/ε(ω) – 1/ε0 can be treated as a generalized Pekar factor.)

For the S2O82– reduction the solvent reorganization energy in eq 1 refers to the bond-breaking first electron transfer:

This step is rate-controlling; the second step (reduction of the SO4– radical) proceeds significantly faster. It was found in work[6] that the dielectric spectrum for water–EG mixtures reveals a 3-D behavior:where τ are characteristic relaxation times.

For such solvents the correlation function M(τ) can be exactly expanded into the sum:[15,16]where the δ are contributions to the solvent reorganization energy λs from ith solvent mode (∑3δ = 1) and τ* are their relaxation times.

To perform the expansion of M(τ), we used the data on the water −EG dielectric spectra reported in ref (6) and employed an inverse Laplace transform technique using the Mathematica 8 program suite. The results are shown in Figure 1 and reveal a different behavior of the different solvent modes. The lowest relaxation time contributes slightly to the solvent reorganization energy; however, its effect becomes more noticeable at higher EG content. The contribution of the fast mode decreases with the EG content, while the middle mode changes opposite. The most important consequence of expansion 5 is that we can introduce three independent solvent coordinates q⃗ = (q1,q2,q3) to address the solvent contribution to the multidimensional reaction free-energy surface E. In the weak coupling limit (small reactant–electrode orbital overlap) this surface can be recast in the form:andwhere r is the intramolecular coordinate (the O–O bond length in peroxodisulfate); indices i and f refer to the initial and final states, respectively; η is electrode overvoltage. The “solvent” parts of eqs 6 and 7, U(q⃗) and U(q⃗), take the form:and

Figure 1

Three solvent relaxation times τ* (a) and their contributions to the reorganization energy δ. (b) Calculated for EG-water mixtures using experimental dielectric spectra (see eq 5).

The potentials describing the intramolecular reorganization of a peroxodisulphate anion were obtained from DFT calculations in ref (1):where D = 1.86 eV; α = 2.38 Ǻ–1; B = 3.709 eV; and β = 1.95 Ǻ–1.

A theoretical analysis of the outer-sphere electron transfer in solvents with two characteristic relaxation times was performed by the authors in ref (17). They used the Agmond–Hopfield formalism[8] to make numerical estimations of the rate constant. Because we have to deal with a partial differential equation, a similar analysis is hardly possible in our case. (See pertinent discussion in ref (6).) Therefore, we employed Langevin (Brownian) molecular dynamics[18] to calculate the reaction rate (k). A system of four differential equations describes the motion in a 5D reaction space along three solvent coordinates and one intramolecular degree of freedom:where m is the reduced mass of S2O82–; dq1rand, dq2rand, dq3rand are random increments of solvent coordinates; dυrand is a random increment of velocity for the movement along the r coordinate; and γr is the intramolecular friction coefficient.

In eq 11, the first three equations describe the solvent dynamics in an overdamped regime;[4] the solvent coordinates are treated as slow. In opposite, the intramolecular coordinate was found to be fast. The random increments were generated according to the algorithm presented in ref (18). This algorithm presumes constructing Gaussian-like distribution functions for the random increments (see eq 11), which depend on friction coefficients γ. The latter were defined in the following way:[4]

The frequencies ω in eq 12 were estimated from the system of equations:[3]where ωeff = 1013 Hz.

To estimate the intramolecular friction coefficient γr, we took a value of 1 ps for the relaxation time;[7] the corresponding frequency was calculated for the bottom of the Morse potential U*(r). A set of the reactant life times (τlife,) was computed by generating several thousands trajectories. The lifetime for the ith trajectory is defined as time that is required to reach the reaction saddle line starting from the bottom of initial well. The resulting lifetime is their average:where N is the number of trajectories (in our simulations N = 12 000).

Then, the rate constant k can be estimated as 1/<τlife>. The system was integrated by the Verlet method with the help of a Matlab 2009b program with an integration step of 1 ps. The reliability of the computational scheme (eqs 11) was tested first for a simple two-well potential

A crucial point for further estimations of the rate constant is a possible dependence of the solvent reorganization energy on the solution viscosity (EG content), that is, “static” solvent effect (see Introduction). To elucidate this issue, we have performed atomistic molecular dynamics simulations for S2O82–/3– in water–EG mixed solvents. Some details are reported in the next section.

Solvent Reorganization Energy

In general, the solvent reorganization energy should be calculated in the vicinity of the activation barrier, that is, for an elongated O–O bond length. For the sake of simplicity, however, all molecular dynamics simulations were performed for the geometry of the optimized oxidized state, S2O82– (independent of the total charge of reactant, −2 or −3).

The frequently used way to calculate the solvent reorganization energy from molecular dynamics simulations is based on the following equation for the reaction free-energy surface G as function of a collective solvent coordinate ΔE (see, for example, ref (19))where ΔE is the energy of Coulomb interaction of a reactant or product with solvent molecules; f is a distribution function that shows a probability to reach a certain ΔE interval due to statistical fluctuations of solvent molecules; and indexes i and f refer to reactant and product, respectively.

After G(ΔE) and G(ΔE) are built, one can fit both curves (assuming a linear response from the solvent environment) by intersecting parabolas and estimate both the activation barrier and the reorganization energies of reduction (λ⃗s) and oxidation (λ⃖s) as a difference between nonequilibrium and equilibrium solvation energies of reactant and product, respectively. Some results of such calculations obtained for x(EG) = 0.5 are presented as an illustration in the Supporting Information (see Figure 1S, Section A2). This approach is, however, rather computer time demanding, as its accuracy noticeably depends on the ensemble size and simulation time. A special bias sampling technique is employed as well to perform more efficient calculations.[19] In this work we use a more simplified method to estimate λ⃗s and λ⃖s:where Eeq(S2O82–) is the interaction energy of S2O82– with solvent molecules in configurations, which are in equilibrium for S2O82–; Enoneq(S2O82–) is the interaction energy of S2O82– with solvent molecules in configurations that are equilibrium for S2O83–. In turnwhere Eeq(S2O83–) is the interaction energy of S2O83– with solvent molecules in configurations that are in equilibrium for S2O83– and Enoneq(S2O83–) is the interaction energy of S2O83– with solvent molecules in configurations which are in equilibrium for S2O82–.

Of course, the linear response theory (LRT) is sometimes too crude to describe solvent contribution to the activation barrier of electron transfer in a proper way. Deviations from the LRT were thoroughly analyzed in work;[20] a new three-parametric model was developed as well.[20] The behavior of the transfer coefficient in interfacial electron transfer reactions was theoretically investigated in ref (21), assuming a nonlinear solvent response. Redox reactions at the interface of an active protein to water might be another example of systems, where nonlinear effects (nongauss fluctuations of solvent molecules) take place.[22] Nevertheless, the solvent reorganization energy still remains a convenient and commonly adopted language to describe the Franck–Condon barrier in homogeneous electron-transfer reactions in mixed solvents.[9−13] It is important to note that in this work we are interested first of all in addressing the qualitative behavior of the solvent reorganization energy as a function of the EG content in solutions; this justifies to some extent the simplified approach used to estimate the λs values (eqs 16, 17).

Potential Energy Functions and Details of the Simulations

The geometry of peroxodisulfate ion with charge q = −2 was optimized by the Möller–Plesset many-body perturbation theory (MP2) using the standard 6-31++G(d, p) basis set. Atomic partial charges in the S2O82– and S2O83– ions were computed with the restrained electrostatic potential (RESP) model.[23] All interactions are described with Lennard-Jones pair potentials and atomic partial charges:where q and q are the partial charges of atoms i and j; r is the distance between these atoms; σ is the van der Waals diameter, σ = 1/2(σ + σ); and ε is the depth of the potential well, ε = (εε)1/2.

The interaction parameters and partial charges used in this work are listed in Table 1S (see Supporting Information, Section A1). The Amber force fields[24] were adapted to reproduce the solvation free energy of S2O82– by slightly modifying the partial charges of oxygen and sulfur atoms and were employed for S2O82– and S2O83– interacting with water and EG. For water, we used the SPC/E model,[25] which describes fairly well the radial distribution function of the oxygen atoms, other structural parameters, and the water dielectric constant.[26] The SPC/E model is rigid with oxygen–hydrogen distances fixed at 1.0 Å and the valence angle at 109.47°. Partial charges of −0.847q and +0.423q reside on oxygen and hydrogen, respectively (q is the elementary charge). For the EG–EG interactions we used a model developed by Kusalik et al.,[27] where authors used the all-atom AMBER/OPLS force field[28] with MM3[29] torsions, as shown in Table IS in the Supporting Information.

Simulations on peroxodisulfate in liquid water were performed at constant pressure and temperature using the Berendsen algorithm for maintaining constant pressure and temperature in a cubic box (starting dimensions 25 × 25 × 25 Å3) filled with 503 water molecules. The equations of motion were solved using the Verlet algorithm with a 1 fs time step. Long-range Coulomb forces were addressed by the Ewald method. All of the other interactions were calculated within a sphere with a radius of Rcutoff = 9 Å.

Calculations on pure EG were performed in a periodic cubic box (starting size 37.6 × 37.6 × 37.6 Å3) filled with 570 EG molecules. Here we have used the Parrinello–Rahman algorithm with 10.5 Å cutoff for maintaining constant pressure and temperature; the C–H and O–H bonds were constrained using the LINCS algorithm. The value of the static dielectric constant is important to study reorganization energies of pure and binary solutions.[11] The model we employed gives a value of 37 for the dielectric constant of pure EG, which is close to experimental data.[30] The Gibbs energy of solvation of peroxodisulfate ion in EG/water mixtures was calculated by the thermodynamic integration method;[31] pertinent corrections[32−37] were addressed as well. Electrostatic and van der Waals contributions were integrated at 12 discrete steps along 0 ≤ λ ≤ 1. The simulation time for each value of λ was 240 and 40 ps discarded for equilibration. All simulations were performed using the GROMACS 4.5.3 package of computer codes.[38]

Results and Discussion

The solvation free energies ΔGsolv of S2O82– calculated for four values of the EG mole fraction are collected in Table 1.

Table 1

Solvation Free Energy of Peroxodisufate in Pure Water and EG and in Mixed Solutions Calculated from MD Simulationsa

x(EG)	0	0.3	0.5	1
–ΔGsolv, kcal mol–1	200 (176.4)	199	196.5	174.4 (174.2)

Solvation free energies obtained from DFT calculations and the polarized continuum model (PCM)[50] are given in parentheses.

The solvation energy decreases from water to pure EG. Because there are no available experimental data on the solvation energies of a peroxodisulfate anion, additional estimates of ΔGsolv were obtained for pure water and EG from DFT calculations with the polarizable continuum model (PCM). For EG the agreement between this estimate and the one obtained from MD simulations is surprisingly good, while for water the difference between two predictions is larger.

As can be seen from Figure 2, the shape of radial distribution functions calculated for S2O82– and S2O83– in water differ noticeably from each other. The RDF for S2O83– has a shoulder and a main peak, while for S2O83– it reveals two peaks and a small plateau in the vicinity of the first RDF minimum. Similar features can be observed for the RDF with EG (Figure 3). Two smooth maxima on the RDF for the peroxodisulphate anion become sharper for its reduced form. The coordination numbers (ncoord) characterizing the nearest solvation sheath of the anion were found about two times smaller for EG as compared with water. (See Table 2.) These quantities change only slightly going from the oxidized to reduced state. For water, we observe an increase in ncoord by 1; in contrast, ncoord for EG decreases by 1.

Figure 2

Coordination number (dashed line) plotted up to first minimum of the radial distribution function (solid line, left) calculated for S2O82– (a) and S2O83– (b) in water (x(EG) = 0).

Figure 3

Coordination number (dashed line) plotted up to first minimum of the radial distribution function (solid line, left) calculated for S2O82– (a) and S2O83– (b) in ethylene glycol (x(EG) = 1).

Table 2

Coordination Numbers Corresponding to the Nearest Solvation Sheath of S2O82– and S2O83– in Water and EG Calculated from the Molecular Dynamics Simulations

	water		ethylene glycol
	S2O82–	S2O83–	S2O82–	S2O83–
coordination number	25	26	13	12
distance of the first minimum of RDF (nm)	0.577	0.564	0.673	0.638

The solvent reorganization energies are presented in Table 3. The environmental reorganization was found to be strongly asymmetric (the effect of the asymmetric intramolecular reorganization on the activation barrier of interfacial electron transfer reactions was thoroughly investigated in work[51]), that is, λ⃗s and λ⃖s values differ significantly; the asymmetry parameter v = λ⃗s/λ⃖s ranges from 0.2 to 1.2. The average values (λs) were computed using the formula:[39,40]

Table 3

Reorganization Energies (kcal mol–1) of the Reduction (Oxidation) of a Peroxodisulphate Anion in Different Water–EG Mixturesa

x(EG)	λ⃗s	λ⃖s	λs (eq 19)
0	22.4 (15.8)	54.5 (60)	33.2 (27.5)
0.3	18.7 (17.1)	68 (32.3)	32.1 (22.9)
0.5	25.6	54.7	36.3
1.0	38.3	32.2	34.9

Results obtained for T = 350 K are given in parentheses.

It can be seen from Table 3 that the averaged reorganization energy increases slightly and nonmonotonously with increasing EG content. In other words, the activation barrier of electron transfer does not decrease with the growth of the solution viscosity and can be treated nearly constant. This is the most important and qualitatively interesting feature predicted on the basis of atomistic MD simulations. In a previous work[12] a very slight change of the solvent reorganization energy for the thermal electron transfer within a Co(NH3)5(pz)3+/Fe(CN)64– ion pair in water–EG mixtures (x(EG) = 0.09 to 0.58) was assumed using the analysis of the Pekar factor. Note that the reduction of S2O32– at a metal electrode is a heterogeneous electron transfer. To address the effect of a metal surface in terms of a continuum approach, one should correct the solvent reorganization energy by the image term:[41]where x is the distance from the reactant center to the image plane (which is shifted ∼0.1 nm from the position of the first layer of metal atoms toward solution).

If we take a value of 0.4 nm as a reasonable estimate for x, the image term ranges from 9.5 (EG) to 11 kcal mol–1 (water) and reduces the λs values obtained from MD simulations. Note that quantum solvent modes do not contribute to the classical activation barrier and should be cut from the outer-sphere reorganization energy. This effect is normally addressed by a coefficient ξ (<1) reducing λs:[3]where ε(ω) is the complex dielectric spectrum of the medium and the boundary cyclic frequency ω* roughly separates classic and quantum regions of solvent modes.

For pure water, the coefficient ξ was previously estimated to be 0.83. It is difficult to make pertinent estimations for EG and its mixtures with water because the available dielectric spectra poorly describe the high-frequency region.[6] Nevertheless, corrections of λs for the pure water image term (x(EG) = 0, see Table 3) and for quantum modes yield a value of 17.8 kcal mol–1, which is very close to an estimate of the environmental reorganization energy of the reduction of S2O82– made in ref (1) on the basis of a continuum model.

We also tried to investigate the influence of temperature on the reorganization energy for pure water and a certain EG–water mixture. (See Table 3.) The increase in temperature leads to a significant decrease in λs for a pure aqueous solution; in the mixed solution at x(EG) = 0.3 the decrease in the reorganization energy is smaller. This temperature effect cannot be explained in terms of a temperature-dependent Pekar factor.[42−44] According to the results of molecular dynamics simulations reported in works[43,44] the solvent reorganization energy can be recast as a sum of two terms. The first term arises from the orientational structure factors, while the second contribution originates from solvent density fluctuations. The second term was found by the authors[43,44] and is basically responsible for the negative temperature dependence of the solvent reorganization energy.

All Brownian dynamics simulations were performed at room temperature; a value of 17.8 kcal mol–1 was taken for λs, which does nor depend on the EG content (see above). We observed that the reaction system moves slowly in “grooves” on the free-energy surface along the solvent coordinates permanently making fast attempts to overcome the saddle line along the intramolecular degree of freedom. The rate constants versus the EG content calculated at several electrode overvoltages are shown in Figure 4. It can be seen from Figure 4 that starting from η = 2.9 V the k versus x(EG) dependencies increase with increasing EG content. The ascending plots become more pronounced at higher overvoltages. Therefore, the results obtained are qualitatively in good agreement with the experimental data.[6] Another quantity that is useful to describe electrochemical redox systems is transfer coefficient α:

Figure 4

Rate constant of the S2O82– reduction from content in water–EG mixtures versus EG calculated at different electrode overpotentials (2.8 to 3.1 V).

The transfer coefficient α calculated as a function of x(EG) is shown in Figure 5; a linear approximation was used to fit the ln k(η) dependencies. The α values are noticeably <0.5 which is typical for electron-transfer reactions proceeding in the vicinity of activationless discharge, that is, at large electrode overvoltages. The increase in α with the growth of EG content (solvent viscosity) agrees with the experimental findings previously reported[6] and the data obtained for water–glucose (sucrose) solutions.[7] Other examples of a nonmonotonous dependence of transfer coefficient on electrode overpotential calculated on the basis of Sumi–Marcus model as a function of solution viscosity are also presented in work.[45]

Figure 5

Symmetry coefficient describing the reduction of a peroxodisulphate anion at a mercury electrode from mixed water–EG solutions calculated as a function of x(EG).

We also calculated the Kramer’s transmission coefficient (κKr) (one should distinguish between the Kramer’s transmission coefficient and electronic transmission coefficient; the latter is assumed to be nearly 1 in our case), which is a measure of the effect of solvent dynamics, or, in other words, of the deviation from simple transition state theory (TST):where k is the rate constant calculated from the Brownian dynamics simulations and kTST is calculated on the basis of TST.

The κKrvalues were found to be noticeably <1 (see Figure 6), which points to the significant role of solvent dynamics effect.

Figure 6

Kramer’s transmission coefficient (κKr) describing the S2O82– reduction versus the EG content in water–EG mixtures calculated at different electrode overpotentials (η).

The saddle point avoidance (SPA) is one of the most remarkable manifestations of solvent dynamics effects; the larger this effect is, the slower the reaction rate is. Because in our case the intramolecular coordinate is fast, it is reasonable to describe the SPA as a sum of contributions from all solvent modes:where q* are the solvent modes corresponding to the reaction window and qsaddle are the solvent coordinates of the saddle point of 5D free-energy surface.

The dependence of δqsolv on x(EG) calculated at η = 3 V is plotted in Figure 7. As shown in our analysis, the nonmonotonous behavior of the model k versus x(EG) curves can be explained in terms of the SPA; namely, this effect was found to be the largest at the vicinity of the k(x(EG)) minima.

Figure 7

Saddle-point avoidance effect (see eq 23) as a function of the EG content presented at η = 3.1 V.

In former works,[46−48] another version of Brownian dynamics (developed by Kast et al.[49]) was employed to address solvent dynamics effects in molecular dynamics simulations of electron-transfer reactions. In this method, friction effects are addressed by a bath of virtual particles with a Maxwell distribution of velocities that collide with the real particles. These authors[46−48] performed molecular dynamics simulations on 2D[46] and 3D[47,48] reaction free-energy surfaces to obtain the rate constants.

Concluding Remarks

With the help of molecular dynamics simulations we investigated the reduction of peroxodisulfate anion solvated in water–EG mixtures. The outer-sphere reorganization energy of electron transfer at the reduction of S2O82– from mixed solutions was calculated in a simplified way. The averaged reorganization energy increases slightly with increasing EG fraction in the mixed solutions. This means that the static solvent effect does not contribute to the presence of a minimum in the experimental current versus solvent viscosity curves describing the electroreduction of a peroxodisulphate anion in these solvents.[6] This effect must, therefore, result from intrinsic dynamical properties of the key solvent modes that contribute to the observed electron-transfer rate. Our numerical results resting on the Brownian dynamics simulations with a nearly constant λs performed on a 5D reaction free-energy surface confirm this conclusion. The characteristic times of these modes and their contributions to the reaction free-energy surface can be extracted from the analysis of solvent correlation function. Each mode can be attributed to some special motion of EG and water molecules (rotation, the change of torsion angles, the lifetime of intra- or intermolecular H bond, etc). Such a “multi-mode” treatment of solvent dynamics in the Zusman’s sense[4,16,53] seems to be crucial, because our calculations performed in the framework of the “one-mode” Sumi–Marcus model (with a constant λs) fail to describe the experimental effect even in a wide interval of the solution viscosity values (see Figure 2S in Supporting Information, Section A3). One may state that various solvent modes in electron-transfer reactions do not behave chaotically; their interplay resembles somewhat an orchestra where each solvent mode plays the role of one instrument. To judge, however, about the extent of generality of this finding, one needs first of all to investigate systematically solvent dynamics effects in pure “non-Debye” solvents (like glycerol, for example) from both experimental and theoretical viewpoints. The most remarkable manifestation of such an effect might be expected for interfacial electron transfer reactions occurring in the vicinity of activationless discharge because in this case large electrode overpotentials result in a significant decreasing the activation barrier. Thus, for a homogeneous redox process in water–EG mixtures the authors of ref (12) observed only a slight decreasing of the rate constant with the growth of EG content. The significant asymmetry of solvent reorganization is another conspicuous feature that deserves explanation at a microscopic level. A deeper insight into this feature can be gained from the analysis of the power spectra (the Fourier transform of the velocity autocorrelation function) describing the solvation shell of a peroxodisulphate anion in oxidized and reduced states.