Dependence of Vibronic Coupling on Molecular Geometry and Environment: Bridging Hydrogen Atom Transfer and Electron-Proton Transfer.

The rate constants for typical concerted proton-coupled electron transfer (PCET) reactions depend on the vibronic coupling between the diabatic reactant and product states. The form of the vibronic coupling is different for electronically adiabatic and nonadiabatic reactions, which are associated with hydrogen atom transfer (HAT) and electron-proton transfer (EPT) mechanisms, respectively. Most PCET rate constant expressions rely on the Condon approximation, which assumes that the vibronic coupling is independent of the nuclear coordinates of the solute and the solvent or protein. Herein we test the Condon approximation for PCET vibronic couplings. The dependence of the vibronic coupling on molecular geometry is investigated for an open and a stacked transition state geometry of the phenoxyl-phenol self-exchange reaction. The calculations indicate that the open geometry is electronically nonadiabatic, corresponding to an EPT mechanism that involves significant electronic charge redistribution, while the stacked geometry is predominantly electronically adiabatic, corresponding primarily to an HAT mechanism. Consequently, a single molecular system can exhibit both HAT and EPT character. The dependence of the vibronic coupling on the solvent or protein configuration is examined for the soybean lipoxygenase enzyme. The calculations indicate that this PCET reaction is electronically nonadiabatic with a vibronic coupling that does not depend significantly on the protein environment. Thus, the Condon approximation is shown to be valid for the solvent and protein nuclear coordinates but invalid for the solute nuclear coordinates in certain PCET systems. These results have significant implications for the calculation of rate constants, as well as mechanistic interpretations, of PCET reactions.

Introduction

Proton-coupled electron transfer (PCET) reactions encompass a broad spectrum of mechanisms.[1−6] For example, such reactions may be sequential or concerted, depending on the existence of a stable intermediate. If a stable intermediate is observed, then the reaction is clearly sequential, but detection of the intermediate may depend on the experimental apparatus. Thus, the absence of an observable intermediate does not definitively imply a concerted reaction. In some cases, the reaction can be determined to be concerted if the products of single electron transfer (ET) and single proton transfer (PT) are much less thermodynamically favorable than the product of the concerted mechanism. This information can often be obtained from the pKa’s and reduction potentials. This paper will focus on PCET reactions that have been determined to be concerted through such an analysis.

Concerted PCET reactions may be further broken down into hydrogen atom transfer (HAT) and electron–proton transfer (EPT).[7] Traditionally, HAT reactions are characterized by the electron and proton transferring between the same donor and acceptor and hence do not involve a significant change in the electronic charge distribution. In contrast, EPT reactions are characterized by the electron and proton transferring between different donors and acceptors and thus result in a significant change in the electronic charge distribution. According to these traditional definitions, the donor and acceptor could be an atom, a molecular orbital (MO), or a chemical bond, although such definitions are not rigorous because the quantum mechanical electron and proton tend to be delocalized, and the MO or chemical bond analysis depends on the representation and level of theory. To provide a more quantitative and well-defined distinction, previously HAT and EPT reactions were shown to be associated with electronically adiabatic and nonadiabatic proton transfer, respectively.[8,9] The distinction between HAT and EPT by the degree of electron–proton nonadiabaticity is consistent with the traditional characterizations mentioned above because the nonadiabatic coupling along the proton transfer coordinate reflects the change in electronic charge distribution as the proton transfers. Thus, a significant change in charge distribution is associated with the electronically nonadiabatic EPT reaction but not with the electronically adiabatic HAT reaction.

Several diagnostics have been devised for distinguishing between HAT and EPT reactions in terms of electron–proton nonadiabaticity.[8−10] A semiclassical formalism developed by Georgievskii and Stuchebrukhov[11] can be used to calculate an effective proton tunneling time τp and an electronic transition time τe as well as an adiabaticity parameter that is defined to be the ratio of these two quantities, p = τp/τe. The reaction is electronically adiabatic if p ≫ 1 because the electrons respond instantaneously to the proton motion, and the system remains on the electronic ground state. The reaction is electronically nonadiabatic if p ≪ 1 because the electrons are unable to respond instantaneously to the proton motion, and excited electronic states are involved. Another diagnostic of electron–proton nonadiabaticity is the magnitude of the first-derivative nonadiabatic coupling between the ground and first excited electronic states along the proton coordinate. A related diagnostic is the magnitude of the change in the electronic charge distribution along the proton coordinate, as reflected by the dipole moment or partial atomic charges.

These diagnostics have been applied to several molecular systems and, more recently, to an enzymatic system. The prototypical examples are the benzyl-toluene and phenoxyl-phenol self-exchange reactions. The former has been shown to be electronically adiabatic (HAT), while the latter has been shown to be electronically nonadiabatic (EPT) according to the diagnostics for electron–proton nonadiabaticity.[8,9] These systems, as well as related systems, have also been studied with other theoretical methods.[12−15] More recently, the PCET reaction catalyzed by the soybean lipoxygenase (SLO) enzyme was shown to be electronically nonadiabatic (EPT) by applying these diagnostics to a gas-phase model system.[16] All of these systems were shown to be vibronically nonadiabatic in that the overall vibronic coupling is small compared to the thermal energy kBT, thereby validating the use of a golden rule rate constant expression. This vibronic nonadiabaticity is related to the response of the solvent or protein environment to the electron–proton subsystem and is determined by different criteria that have been discussed elsewhere.[3,17] Theoretical calculations based on the nonadiabatic treatment of the SLO enzyme have reproduced the experimentally observed hydrogen/deuterium kinetic isotope effect of ∼80 in the wild-type enzyme and up to ∼500 in mutant enzymes.[18−23]

The nonadiabatic PCET rate constant expressions rely on the Condon approximation for the vibronic coupling.[3,24,25] In nonadiabatic electron transfer theory, the Condon approximation is based on the assumption that the electronic coupling is independent of the nuclear configuration.[26−29] For vibronically nonadiabatic PCET reactions, the Condon approximation is based on the assumption that the vibronic coupling is independent of the nuclear configuration, including both the molecular geometry of the PCET solute complex and the solvent or protein environment. An exception is that the dependence of the vibronic coupling on the proton donor–acceptor distance is included explicitly in the nonadiabatic PCET rate constant expressions.[25] For situations in which the Condon approximation breaks down, a given PCET system could span the electronically adiabatic and nonadiabatic regimes. For both electron transfer and PCET systems, these two regimes may be spanned as the donor–acceptor distances are varied, but the dependence of the coupling and the degree of nonadiabaticity on other geometrical coordinates, as well as the environmental configuration, is less obvious.

In this paper, we examine the dependence of the magnitude of the vibronic coupling and the degree of electron–proton nonadiabaticity on the molecular geometry and on the solvent or protein configuration. The dependence on molecular geometry is investigated for the phenoxyl-phenol self-exchange system because two transition states corresponding to either an open or a stacked geometry have been identified for this type of system.[13] In particular, for the related benzyl-toluene system, transition states have been optimized with either an open geometry[12] or a stacked geometry,[13] and the stacked geometry was found to be lower in energy by 4.0 kcal/mol at the level of theory used in ref (13). The present study focuses on the phenoxyl-phenol system, which also exhibits both types of transition state geometries. The open geometry of the phenoxyl-phenol system has already been shown to be electronically nonadiabatic,[8,9] but herein we also study the stacked geometry. The dependence on the solvent and protein environment is investigated for the PCET reaction catalyzed by the SLO enzyme. The previous studies focused on a fully quantum mechanical gas-phase model of the substrate–cofactor complex.[16] Herein we use mixed quantum mechanical/molecular mechanical (QM/MM) methods to include the protein and solvent environment. These calculations test the validity of the Condon approximation for PCET vibronic couplings in solution and proteins. The results have significant implications for applications to PCET in chemical and biological processes.

Theory and Computational Methods

General Theory

The semiclassical formalism for calculating the effective proton tunneling time and electronic transition time is described in detail elsewhere.[8,11] In this section, we only provide the expressions that are used to calculate the quantities necessary to determine the adiabaticity parameter and the semiclassical vibronic coupling. The effective proton tunneling time isand the effective electronic transition time isHere Vel is the electronic coupling between the two diabatic electronic states, |ΔF| is the difference between the slopes of the diabatic proton potential energy curves at the crossing point, and v = (2(Vc – E)/mp)1/2, where Vc is the energy at which the potential energy curves cross, mp is the proton mass, and E is the tunneling energy, which is defined as the energy of the degenerate proton vibrational levels in the reactant and product potential wells. As mentioned above, the adiabaticity parameter is defined as p = τp/τe.

The vibronic coupling can be calculated in several different ways, including a full basis set Hamiltonian matrix diagonalization or a semiclassical approach, which have been shown to provide numerically equivalent results for these types of systems.[9] For systems that are known to be in the electronically adiabatic or nonadiabatic regime, expressions derived in these limits can be utilized. Specifically, the vibronic coupling in the electronically adiabatic regime, denoted V(ad), is half the tunneling splitting associated with the ground electronic state. The vibronic coupling in the electronically nonadiabatic regime, denoted V(nad), is given by the following expression:where S is the overlap integral between the proton vibrational wave functions calculated for the reactant and product diabatic potentials. In principle, this overlap can be calculated for any pair of reactant and product proton vibrational wave functions, but in this paper we calculate it for only the ground proton vibrational states.

The semiclassical coupling spans the electronically adiabatic and nonadiabatic regimes and is expressed aswhereHere Γ(x) is the gamma function, and p is the adiabaticity parameter defined above. The derivation of these equations is based on the general semiclassical tunneling model and is given in ref (11). In this paper, we calculated the vibronic coupling with all of these methods for the phenoxyl-phenol system to enable a comparison. For the SLO system, however, the vibronic coupling was calculated only with the electronically nonadiabatic expression given in eq .

The diabatic proton potential energy curves can be calculated in two different ways. For relatively small molecular systems, the adiabatic proton potential energy curves associated with the ground and first excited adiabatic electronic states can be calculated along the proton coordinate with the complete active space self-consistent field (CASSCF) method.[30,31] Subsequently, a diabatization procedure that is exact for two states along a single coordinate (i.e., the nonadiabatic coupling between the two states is identically zero along this coordinate)[9,32] can be used to generate the diabatic proton potential energy curves from these adiabatic electronic states as well as the nonadiabatic coupling between them.[9] Alternatively, constrained density functional theory-configuration interaction (CDFT-CI)[33−35] can be used to generate the diabatic proton potential energy curves and the corresponding electronic couplings. In CDFT-CI, the coupling between the two constrained states is approximated as the off-diagonal Hamiltonian matrix element between the two Slater determinants comprised of the Kohn–Sham orbitals for the constrained states.[33−35] We applied both the CASSCF and CDFT-CI methods to the phenoxyl-phenol molecular system and showed that these two methods lead to qualitatively similar results. Due to computational limitations, we applied only the CDFT-CI method to the SLO enzymatic system and included the protein and solvent environment using a QM/MM approach.

The other two diagnostics for electron–proton nonadiabaticity are the nonadiabatic coupling and the change in electronic charge distribution along the proton coordinate. We calculated the nonadiabatic coupling between the lowest two adiabatic electronic states along the one-dimensional proton coordinate with the CASSCF method for the phenoxyl-phenol system. This scalar coupling is defined aswhere Ψlel(re;rp) and Ψ2el(re;rp) are the ground and first excited adiabatic electronic states, respectively, along the proton coordinate rp. In addition, the dipole moment of the ground electronic state as the proton moves along the proton donor–acceptor axis was calculated using CASSCF for the phenoxyl-phenol system and using ground state DFT for the SLO system.

Computational Details for Phenoxyl-Phenol System

For the phenoxyl-phenol calculations, two different transition state structures, denoted the “open” and “stacked” geometries, were optimized with density functional theory (DFT) at the M06-2X/6-311+G** level of theory[36−38] in the gas phase using Gaussian09.[39] These structures are depicted in Figure . At this level of theory, the stacked transition state structure is 4.1 kcal/mol lower in energy than the open transition state structure. For each structure, the adiabatic proton potential energy curves associated with the ground and first excited adiabatic electronic states were obtained by calculating the state-averaged CASSCF energies in the gas phase as the transferring proton was moved along a grid spanning the proton donor–acceptor axis with all other atoms fixed. On the basis of careful analysis of the active spaces over the range of proton coordinates, CAS(3,6) calculations state-averaged over two states were performed for the open structure, and CAS(7,8) calculations state-averaged over four states were performed for the stacked structure to ensure that the active space was conserved along the proton transfer coordinate. Note that these structures are not transition states at the CASSCF level but are useful symmetric structures for the analysis described below.

Figure 1

Open (a) and stacked (b) transition state structures for the self-exchange phenoxyl-phenol reaction calculated at the DFT/M06-2X/6-311+G** level of theory. The proton is transferring between the two red oxygen atoms.

Additional types of calculations were performed for comparison to these CASSCF results. The effects of dynamic correlation were investigated by performing second-order perturbation theory CASSCF (CASPT2) calculations and comparing the results with those obtained from the CASSCF calculations. The CASSCF energies and first-derivative nonadiabatic couplings reported in the main paper were performed with Molpro,[40] but the comparison of CASSCF and CASPT2 results provided in the Supporting Information was performed with Molcas.[41−43] In addition, CDFT-CI calculations with the ωB97X functional[44] were performed using Q-Chem[45] to obtain the diabatic proton potential energy curves for comparison to the CASSCF results. For the CDFT-CI calculations, the spin density was constrained to be zero on the phenol (left side) and unity on the phenoxyl (right side) fragment for the reactant diabatic state and the reverse for the product diabatic state. The 6-31G** basis set was used for the CASSCF and CDFT-CI calculations. Note that the M06-2X functional was used for the geometry optimizations because it includes dispersion effects, whereas the ωB97X functional was used for calculating the diabatic states and couplings with CDFT-CI because it includes long-range corrections, which are important for describing charge transfer states. A previous benchmarking study[16] illustrated that CDFT-CI calculations with the ωB97X functional resulted in similar diabatic states and couplings as those obtained with the CASSCF method for the phenoxyl-phenol system.

The adiabatic electronic states obtained from the CASSCF calculations were diabatized using the method described in previous work.[9] The proton vibrational wave functions for each diabatic proton potential energy curve were calculated using the Fourier Grid Hamiltonian method.[46] The double adiabatic states defined as products of diabatic electronic states and associated proton vibrational states were used as a basis set to construct a Hamiltonian matrix. Diagonalization of this full basis set Hamiltonian matrix provides the vibronic eigenfunctions and eigenvalues. The vibronic coupling calculated from this full basis set Hamiltonian diagonalization, denoted V(full), is half the energy difference between the two lowest-energy vibronic states. This vibronic coupling was compared to the vibronic coupling calculated with the semiclassical approach, as given in eq , and to the vibronic couplings calculated in the adiabatic and nonadiabatic limits.

Computational Details for Soybean Lipoxygenase System

The SLO calculations included the protein environment based on snapshots from an equilibrated molecular dynamics (MD) trajectory of the protein solvated with explicit TIP3P[47] water molecules. The crystal structure of SLO with PDB code 3PZW(48) was used as the initial protein structure, and the linoleic acid (LA) substrate was docked to the active site of this enzyme. A description of the initial preparation of the system and the equilibration procedure, as well as other computational details related to the classical MD simulations, is provided in the Supporting Information. Three snapshots were obtained from the 20 ns production MD trajectory for the subsequent QM/MM calculations. Because of the large size of this system, the CASSCF calculations performed on the phenoxyl-phenol system were not computationally tractable, and standard QM/MM methods that are typically applied to protein systems[49,50] were used instead.

For each of the three snapshots obtained from the MD simulation, a QM/MM geometry optimization was performed for a system comprised of the SLO–LA complex and all solvent molecules and ions within 5 Å of at least one atom in the SLO–LA complex. The QM region contained 87 atoms, including Fe–OH, partial side chains of residues 499, 504, 690, 694, and 839, and a portion of the linoleic acid. The QM atoms are depicted in the ball and stick representation in Figure , where the transferring hydrogen atom is highlighted in yellow. The boundary between the QM and MM regions was treated with the hydrogen capping method.[51] In the geometry optimization, only the atoms within 20 Å of Fe were allowed to move. The QM region was treated with DFT using with the B3LYP functional,[52,53] in conjunction with the 6-31G* basis set for all nonmetal atoms and the LANL2DZ[54,55] basis set for the Fe atom. The MM region was described by the OPLS2005 force field.[47,49,56] To obtain a reasonable QM/MM configuration corresponding to the crossing point of the diabatic free energy curves in PCET theory,[3,24] constrained QM/MM geometry optimizations were performed. In one set of optimizations, the donor-proton and acceptor-proton distances were constrained to be 1.32 and 1.38 Å, respectively. In another set of optimizations, both of these distances were constrained to be 1.35 Å. All QM/MM geometry optimizations were performed using the QSite module in the Schrödinger package.[49,50,57]

Figure 2

(a) Depiction of the SLO-LA complex used for the QM/MM calculations, where the QM region is indicated by the ball-and-stick representation, the transferring hydrogen is highlighted in yellow, and the MM region is indicated by the ribbon representation. The proton (yellow) is transferring from the carbon (cyan) to the oxygen (red), and the electron is effectively transferring from the π backbone of the LA substrate (cyan) to the iron (mauve). (b) Depiction of the system used for the QM/MM CDFT-CI calculations, where the QM region is indicated by the ball-and-stick representation, and the surrounding MM charges associated with the protein and solvent environment are indicated by purple spheres.

The diabatic proton potential energy curves were obtained by performing QM/MM CDFT-CI calculations for each of the six optimized geometries.[33−35] In these calculations, the proton was moved along a grid spanning the linear axis connecting the donor carbon and the acceptor oxygen atoms with all other nuclei fixed. The MM atoms within 10 Å of at least one of the QM atoms were included in the MM region for the QM/MM CDFT-CI calculations. These MM atoms were treated as external charges with magnitudes defined by the OPLS2005 force field. Each MM charge was represented as a Gaussian blurred charge with a width of 3 Å to describe their electrostatic interactions with the QM electrons and nuclei for the CDFT-CI calculations. The ωB97X functional with the 6-31G** basis set was used for the QM region. The diabatic states were obtained by constraining the spin densities on the LA and Fe-cofactor fragments to be 0 and 5, respectively, for the reactant state and 1 and 4, respectively, for the product state.[16] The QM/MM CDFT-CI calculations were performed with Q-Chem.[45] The ground proton vibrational wave functions were calculated for each diabatic state using the FGH method, and the vibronic coupling was calculated using the electronically nonadiabatic expression given in eq .

Results

Phenoxyl-Phenol System

The open and stacked transition state geometries of the phenoxyl-phenol system are depicted in Figure . The proton donor–acceptor O–O distance is 2.4 Å in both optimized structures. The O–H–O angle is 180° for the open and 166° for the stacked geometry. These structural properties indicate the presence of a reasonably strong hydrogen bond in both geometries. The vibrational mode associated with the transition state imaginary frequency corresponds to proton transfer between the two oxygen atoms for both geometries. Previous studies indicated that the self-exchange reaction for the open geometry of the phenoxyl-phenol system corresponds to the EPT mechanism and is electronically nonadiabatic, while the self-exchange reaction for the open geometry of the benzyl-toluene system corresponds to the HAT mechanism and is electronically adiabatic.[8,12] For the benzyl-toluene system, a stacked transition state geometry has been found to be lower in energy than the open geometry of this system.[13]

Figure depicts the two highest-energy occupied MOs for both the open and stacked geometries of the phenoxyl-phenol system, where the lower MO is doubly occupied and the higher MO is singly occupied. For the open geometry, the PT interface region of these MOs is dominated by 2p orbitals perpendicular to the proton donor–acceptor axis with a π-bonding interaction in the doubly occupied MO. In contrast, for the stacked geometry, the PT interface region of these MOs is dominated by atomic orbitals oriented along the proton donor–acceptor axis with a σ-bonding interaction in the doubly occupied MO. The character of the MOs in the PT interface region for the stacked geometry is similar to that observed for both the open[12] and stacked[13] geometries of the benzyl-toluene system, which was determined to be electronically adiabatic.[8,9] Another significant difference between the MOs for the open and stacked geometries of the phenoxyl-phenol system is that the stacked geometry exhibits a π-stacking interaction between the two aromatic rings. As depicted in Figure , the doubly occupied and singly occupied MOs exhibit bonding and antibonding interactions, respectively, between the aromatic ring orbitals, resulting in a net bonding interaction between the ring moieties for the stacked geometry. This π-stacking bonding interaction increases the electronic coupling, thereby decreasing the electronic transition time relative to the effective proton tunneling time. This analysis suggests that the stacked geometry of the phenoxyl-phenol system may be associated with electronically adiabatic HAT, in contrast to the previously studied open geometry of this system, which was determined to be associated with electronically nonadiabatic EPT.

Figure 3

Highest occupied MOs for the open (left) and stacked (right) geometries of the phenoxyl-phenol system for the dominant configuration obtained from the CASSCF/6-31G** calculations. The lower MO is doubly occupied, and the higher MO is singly occupied.

The CASSCF and CDFT-CI proton potential energy curves for both the open and stacked geometries of the phenoxyl-phenol system are depicted Figure . For each geometry, the CASSCF and CDFT-CI proton potential energy curves are qualitatively similar, thereby supporting the use of the CDFT-CI method with the ωB97X functional for other PCET systems, including the SLO system discussed below. Additional CASSCF calculations that included four adiabatic electronic states were also performed, but the second and third excited states were found to be much higher in energy than the first excited state (Figure S1), providing validation for the use of a two-state model. Moreover, the adiabatic proton potential energy curves were also calculated with CASPT2 to examine the effects of dynamical correlation, and the CASSCF and CASPT2 curves were found to be similar (Figure S2).

Figure 4

Adiabatic (black dashed lines) and diabatic (blue and red solid lines) proton potential energy curves for the open (left panels) and stacked (right panels) geometries of the phenoxyl-phenol system obtained using the CASSCF/6-31G** (upper panels) and CDFT-CI/ωB97X/6-31G** (lower panels) methods.

Figure depicts both the adiabatic (black dashed lines) and diabatic (blue and red solid lines) proton potential energy curves. These curves are qualitatively different for the open and stacked geometries. In particular, the splitting between the ground and first excited adiabatic states is much greater for the stacked geometry. Moreover, the adiabatic and diabatic proton potential energy curves are virtually indistinguishable except in the crossing region for the open geometry but differ significantly for the entire range of proton coordinates for the stacked geometry. These differences are consistent with electronically nonadiabatic self-exchange for the open geometry but electronically adiabatic self-exchange for the stacked geometry.

The degree of electron–proton nonadiabaticity was quantified within the semiclassical formalism by calculating the effective proton tunneling time τp and the electronic transition time τe, as well as the adiabaticity parameter p, which is the ratio of these two quantities. The values of these parameters are given in Table . For the open structure, the effective proton tunneling time is much smaller than the electronic transition time, with p ≪ 1, indicating that the reaction is electronically nonadiabatic. For the stacked structure, the effective proton tunneling time is larger than the electronic transition time, with p > 1, indicating that the reaction is predominantly electronically adiabatic. As discussed below, however, the proton tunneling time and electronic transition time are similar for the stacked structure, with a ratio of p = 1.4, so the self-exchange reaction for this geometry can be viewed as being in the intermediate regime between electronically adiabatic and nonadiabatic.

Table 1

Electronic Couplings, Semiclassical Parameters, and Vibronic Couplings Calculated with Various Methods for Open and Stacked Geometries of Phenoxyl-Phenol System

geometry	Vel (cm–1)	τp (fs)	τe (fs)	p = τp/τe
open	376	0.076	14.12	0.00535
stacked	5735	1.272	0.926	1.374

Vibronic couplings given in cm−1.

The first value for V(na) is the matrix element of the product of Vel(rp) and the ground reactant and product proton vibrational wave functions, and the value in parentheses is obtained from eq with Vel calculated at rp = 0 (i.e., the product of Vel and the overlap integral between the ground reactant and product proton vibrational wave functions). The similarity between these two values indicates that Vel does not depend strongly on rp.

Table also provides the vibronic couplings calculated with the full basis set diagonalization method, the semiclassical approach, and the methods that are valid in the adiabatic and nonadiabatic regimes. For both geometries, the full basis set diagonalization and semiclassical couplings are similar to each other because both of these approaches are valid in the adiabatic and nonadiabatic limits as well as the intermediate regime. Thus, the coupling calculated with these two approaches will be denoted the “general” coupling. For the open geometry, the nonadiabatic coupling agrees well with the general coupling, whereas for the stacked geometry, the adiabatic coupling agrees better with the general coupling. These calculations provide further evidence that the open and stacked geometries correspond to electronically nonadiabatic and predominantly electronically adiabatic reactions, respectively.

Figure depicts the first-derivative nonadiabatic coupling vector and the dipole moment vector projected along the proton donor–acceptor axis as the proton moves from the donor to the acceptor. The open geometry exhibits a substantial peak in the nonadiabatic coupling and a drastic change in the dipole moment as the proton moves across the midpoint of the proton donor–acceptor axis, whereas the stacked geometry does not exhibit any significant nonadiabatic coupling and only relatively minor and more gradual changes in dipole moment as the proton moves along this axis. Furthermore, the electrostatic potential maps shown in Figure exhibit significant electronic charge transfer between the two aromatic rings as the proton transfers for the open geometry and a much smaller degree of electronic charge transfer between the two rings as the proton transfers for the stacked geometry. In addition, the spin densities depicted in Figure illustrate that the unpaired spin density shifts from one ring to the other as the proton transfers in the open geometry but remains delocalized over both rings during proton transfer for the stacked geometry. These observations are consistent with electronically nonadiabatic behavior for the open geometry, corresponding to an EPT mechanism, and predominantly electronically adiabatic behavior for the stacked geometry, corresponding more closely to an HAT mechanism.

Figure 5

(a) Component of the first-order nonadiabatic coupling vector, as defined in eq , between the CASSCF/6-31G** ground and first excited adiabatic electronic states as the proton moves along the proton donor–acceptor axis for the open (solid) and stacked (dashed) geometries of the phenoxyl-phenol system. (b) Component of the dipole moment vector as the proton moves along the proton donor–acceptor axis for the CASSCF/6-31G** ground adiabatic electronic state for the open (solid) and stacked (dashed) geometries of the phenoxyl-phenol system.

Figure 6

Electrostatic potential maps for the ground adiabatic electronic states generated with DFT/ωB97X/6-31G** for the reactant (top), transition state (middle), and product (bottom) positions of the transferring hydrogen for the open (left) and stacked (right) geometries of the phenoxyl-phenol system. The density isosurface value is 0.005, and negatively and positively charged regions are indicated by red and blue coloring, respectively.

Figure 7

Spin densities for the open (top) and stacked (bottom) geometries of the phenoxyl-phenol system obtained from CASSCF/6-31G** ground state calculations for the reactant (left) and product (right) positions of the transferring hydrogen.

We emphasize that the stacked geometry does exhibit a small amount of electronic charge redistribution between the two aromatic rings during proton transfer, as indicated by the changes in dipole moment and electrostatic potential, and therefore is not a pure HAT reaction. In other words, the self-exchange reaction for the stacked geometry is not fully electronically adiabatic, as also indicated by the adiabaticity parameter, which is greater than unity but not as large as was observed for the open geometry of the benzyl-toluene system, which is considered to be a pure HAT reaction. The adiabaticity parameter is 1.4 for the stacked geometry of the phenoxyl-phenol system and 3.5 for the open geometry of the benzyl-toluene system.[8] (For further comparison, calculations on the stacked geometry of the benzyl-toluene system[13] are provided in Figure S3, indicating an adiabaticity parameter of 10.0.) On the basis of this analysis, we classify the stacked geometry of the phenoxyl-phenol system as a predominantly electronically adiabatic reaction that can be described as an HAT mechanism with a small amount of EPT character. Consequently, the stacked geometry represents an example of a system that is in the intermediate regime between electronically adiabatic and nonadiabatic, or between HAT and EPT, although it is closer toward the electronically adiabatic HAT limit.

Thus, all of these analyses indicate that the open and stacked geometries of the phenoxyl-phenol system are in different regimes. Specifically, the open geometry is electronically nonadiabatic, corresponding to an EPT reaction, while the stacked geometry is in the intermediate regime but predominantly electronically adiabatic, corresponding to an HAT reaction. As given in Table , the electronic coupling is significantly greater for the stacked geometry than for the open geometry because of the π-stacking interaction between the rings, as indicated by the MOs in Figure . This stacking interaction decreases the electronic transition time to the extent that the electrons are able to respond instantaneously to the proton motion, thereby leading to an electronically adiabatic proton transfer that remains on the electronic ground state. This reaction involves only a small amount of electronic charge redistribution between the rings, supporting the designation of a primarily HAT mechanism for this geometry. In contrast, the open geometry involves a significant shift in electronic charge distribution from one ring to the other, and these rings are further apart with weaker interactions, supporting the designation of EPT for this geometry. These calculations illustrate that a single molecular system can span the electronically adiabatic and nonadiabatic limits as it explores configurational space via thermal fluctuations. Moreover, at certain geometries the reaction may lie in the intermediate regime between the electronically adiabatic and nonadiabatic limits and therefore can no longer be designated as either HAT or EPT.

Soybean Lipoxygenase System

To investigate the impact of the protein and solvent environment on the vibronic coupling, we performed QM/MM calculations on the SLO-LA system depicted in Figure . The classical MD simulations and QM/MM geometry optimizations were performed for the solvated enzyme system depicted in Figure a, where the QM region is shown in the ball-and-stick representation. The QM/MM CDFT-CI calculations were performed for the somewhat truncated system depicted in Figure b. The atomic charges in the MM region, depicted as purple spheres in Figure b, were included in the QM/MM CDFT-CI calculations using an electrostatic embedding method. Two types of constrained QM/MM geometry optimizations were performed for each of three different configurations along the classical MD trajectory, leading to a total of six SLO configurations. For each configuration, the reactant and product diabatic proton potential energy curves, as well as the electronic coupling between these two states, were calculated using the QM/MM CDFT-CI approach.

The reactant and product diabatic proton potential energy curves obtained for one of these configurations are depicted in Figure . These curves were shifted to ensure that the ground proton vibrational energy levels are degenerate. The analogous curves for the other five configurations are provided in the Supporting Information. These proton potential energy curves are similar to each other and to those obtained previously for a gas-phase model of the SLO-LA system.[16] These results demonstrate that the electrostatic effects from the protein environment in this system do not significantly perturb the shape of the diabatic proton potential energy curves.

Figure 8

Diabatic proton potential energy curves for the SLO system calculated with QM/MM CDFT-CI/ωB97X/6-31G** for a configuration obtained by QM/MM geometry optimization with the donor-proton and acceptor-proton distances constrained to be 1.35 Å. The diabatic states were shifted to ensure that the ground proton vibrational energy levels are degenerate. The analogous diabatic proton potential energy curves for five other configurations are provided in Figure S4.

Table provides the electronic couplings, as well as the effective proton tunneling time τp, the electronic transition time τe, and the adiabaticity parameter p, for three different SLO configurations. For this system, the effective proton tunneling time is much smaller than the electronic transition time, with p ≪ 1, demonstrating that this reaction is electronically nonadiabatic. Further evidence of nonadiabaticity is provided by Figure , which depicts the dipole moment of the QM region in the field of MM point charges as the proton moves along the donor–acceptor axis for the ground electronic state obtained with QM/MM DFT/ωB97X/6-31G**. The drastic change in dipole moment as the proton passes through the middle of the proton donor–acceptor axis indicates a substantial amount of electronic charge transfer from the LA substrate to the Fe-cofactor as the proton transfers. In addition, the spin densities depicted in Figure illustrate that the electron effectively transfers from the π backbone of the LA substrate to the iron of the cofactor, resulting in unpaired spin density along the π backbone of the LA substrate after the proton has transferred.

Table 2

Electronic Couplings, Semiclassical Parameters, and Nonadiabatic Vibronic Couplings for the SLO System

geometry	Vel (cm–1)	τp (fs)	τe (fs)	p = τp/τe	V(na) (cm–1)
QM/MM Ia	1500	0.23	3.53	0.07	2.6
QM/MM II	1798	0.27	2.95	0.09	3.6
QM/MM III	2273	0.31	2.34	0.13	3.9
gas-phase modelb	1607	0.22	3.3	0.07	1.8

Geometries I, II, and III were obtained from QM/MM geometry optimizations of three different snapshots from a classical MD trajectory. The QM/MM optimizations were conducted with the donor-proton and acceptor-proton distances constrained to be 1.35 Å. Results for three other geometries are provided in Table S2.

This result was obtained by previous studies on a gas-phase SLO-LA model system.[16]

Figure 9

Dipole moment as the proton transfers from the donor to the acceptor for an SLO configuration obtained by QM/MM geometry optimization with the donor-proton and acceptor-proton distances constrained to be 1.35 Å. The magnitude of the total dipole moment vector (solid blue line) and the dipole moment vector projected onto the axis connecting the donor carbon and the Fe (dashed red line) were calculated for the QM region in the field of MM point charges for the QM/MM DFT/ωB97X/6-31G** ground adiabatic electronic state as the proton moves along the proton donor–acceptor axis. Analogous figures for five other SLO configurations are provided in Figure S5.

Figure 10

Spin densities for the reactant (left) and product (right) positions of the transferring hydrogen corresponding to the minima of the reactant and product diabatic curves, respectively, in Figure for the SLO-LA system. The spin densities were calculated for the QM region in the field of MM point charges for the QM/MM DFT/ωB97X/6-31G** ground adiabatic electronic state.

Because this reaction was determined to be electronically nonadiabatic, the vibronic coupling was calculated using eq , which is valid in the electronically nonadiabatic regime. As shown in Table , the semiclassical parameters and vibronic couplings are similar for the three different protein configurations as well as for the gas-phase model studied previously. This agreement for the four different environments illustrates that the protein environment does not significantly impact the vibronic coupling, thereby providing validation for the Condon approximation.[26−29] As further validation of these findings, calculations using nonadiabatic rate constant expressions relying on the Condon approximation have reproduced the experimentally observed hydrogen/deuterium kinetic isotope effects and their temperature dependencies.[18−23]

Comparison of Phenoxyl-Phenol and SLO-LA Systems

A comparison of the analyses of electron–proton nonadiabaticity for the phenoxyl-phenol and SLO-LA systems indicates that the SLO-LA system is more similar to the open structure than to the stacked structure of the phenoxyl-phenol system. Both the open phenoxyl-phenol and the SLO-LA systems exhibit a substantial change in the electronic charge distribution during proton transfer, as illustrated in Figure b and Figure for the phenoxyl-phenol system and in Figure for the SLO-LA system. The spin densities for these two systems also indicate a significant shift of the unpaired spin density during proton transfer, as depicted by Figure for the phenoxyl-phenol system and Figure for the SLO-LA system. Moreover, the adiabaticity parameter is less than unity for both systems, namely 0.005 for the open phenoxyl-phenol system and ∼0.1 for the SLO-LA system, compared to the value of 1.4 for the stacked phenoxyl-phenol system. On the basis of these analyses, both the open structure of the phenoxyl-phenol system and the SLO-LA system are in the electronically nonadiabatic regime and therefore represent EPT, although the SLO-LA system could be viewed as less nonadiabatic in terms of the adiabaticity parameter.

Conclusions

In this paper, we tested the Condon approximation for PCET vibronic couplings, which strongly impact the rate constants. Calculations of the vibronic coupling for the phenoxyl-phenol self-exchange reaction illustrate that the open geometry is electronically nonadiabatic, while the stacked geometry is in the intermediate regime but is predominantly electronically adiabatic. The electronic coupling is significantly greater for the stacked geometry than for the open geometry because of the π-stacking bonding interaction between the rings in the stacked geometry. Moreover, the self-exchange reaction involves substantially more electronic charge redistribution in the open geometry than in the stacked geometry. On the basis of this analysis, the reaction is identified as EPT in the open geometry but primarily HAT, with a small amount of EPT character, in the stacked geometry. These calculations demonstrate a breakdown in the Condon approximation in that the vibronic coupling depends strongly on the geometry of the PCET complex. This geometric dependence is not simply a dependence on the donor–acceptor distance, which is a well-known phenomenon, but rather is a more interesting dependence on the intramolecular angle between the planes of two aromatic rings. Calculations of the vibronic coupling for the SLO-LA enzyme indicate that this reaction is electronically nonadiabatic, corresponding to EPT, similar to the open geometry of the phenoxyl-phenol system. Furthermore, these calculations illustrate that the vibronic coupling for the SLO-LA enzyme does not depend significantly on the solvent or protein environment.

Thus, these calculations suggest that the Condon approximation is valid for the solvent and protein nuclear coordinates but could be invalid for the solute nuclear coordinates, particularly intramolecular coordinates that influence the π-stacking interactions for systems with aromatic rings. Moreover, a single molecular system can span the electronically adiabatic and nonadiabatic limits through thermal fluctuations that lead to conformational changes. The form of the vibronic coupling is different in these two regimes, as well as in the intermediate regime. The mechanistic interpretation is also different in these two regimes, resulting in the possibility that a single system can exhibit both HAT and EPT character. Thus, simulations of PCET reactions, calculations of PCET rate constants, and mechanistic interpretations should account for the possibility of spanning both regimes, as well as the potential breakdown of the Condon approximation, for certain types of systems.