Proton-Coupled Electron Transfer: Moving Together and Charging Forward.

Proton-coupled electron transfer (PCET) is ubiquitous throughout chemistry and biology. This Perspective discusses recent advances and current challenges in the field of PCET, with an emphasis on the role of theory and computation. The fundamental theoretical concepts are summarized, and expressions for rate constants and kinetic isotope effects are provided. Computational methods for calculating reduction potentials and pKa's for molecular electrocatalysts, as well as insights into linear correlations and non-innocent ligands, are also described. In addition, computational methods for simulating the nonadiabatic dynamics of photoexcited PCET are discussed. Representative applications to PCET in solution, proteins, electrochemistry, and photoinduced processes are presented, highlighting the interplay between theoretical and experimental studies. The current challenges and suggested future directions are outlined for each type of application, concluding with an overall view to the future.

Proton-Coupled Electron Transfer Is Ubiquitous

Proton-coupled electron transfer (PCET) is broadly defined as any process that involves the transfer of at least one electron and one proton.[1−11] The electron and proton can be transferred between the same sites or between different sites, and they can be transferred in the same direction or in different directions. PCET reactions can be sequential, with a stable intermediate corresponding to electron or proton transfer, or concerted, without such a stable intermediate. The distinction between sequential and concerted PCET reactions is not rigorous because it depends on the definition of a stable intermediate, but nevertheless it serves as a useful guide for discussion. In practice, the identification of a sequential mechanism is clear if the intermediate can be isolated experimentally, but this straightforward distinction is not always possible and may require the specification of a lifetime that depends on the experimental apparatus. A concerted mechanism can be identified if the single electron and single proton transfer reactions lead to intermediates that are known to be much less thermodynamically stable than the product of the concerted mechanism according to the reduction potentials and pKa values.

Traditionally, concerted PCET reactions in which the electron and proton transfer between the same donors and acceptors are denoted hydrogen atom transfer (HAT), whereas concerted PCET reactions in which the electron and proton transfer between different donors and acceptors are denoted electron–proton transfer (EPT). In this context, the donor and acceptor may be defined in terms of a molecular orbital,[12] a chemical bond, or an atom, although this distinction is not rigorous because the quantum mechanical electron and proton are delocalized, and the molecular orbital or chemical bond analysis depends on the level of theory and representation. A more quantitative distinction is provided by the degree of electron–proton nonadiabaticity, as defined below.[13,14] In this Perspective, all of these types of processes fall under the general umbrella of PCET.

PCET occurs in a wide range of biological processes, such as photosynthesis and respiration, as well as chemical processes relevant to solar cells and other energy devices. Many of these processes are quite complex, with multiple electrons and protons transferring in disparate locations in different directions. The investigation of more well-defined, simpler model systems is important for elucidating the fundamental physical principles underlying PCET reactions. Subsequently, these physical principles can be applied to the more complex systems by breaking them down into simpler components that can be understood in terms of the model systems. A combination of experimental and theoretical approaches is required to fully understand both the simple and the more complex systems.

A complete review of the vast field of PCET is not possible within the framework of a Perspective. Instead, this Perspective presents a broad overview of the fundamental concepts of PCET and provides examples that are inherently biased but nevertheless are representative of the various types of PCET. These examples focus on the interplay between theoretical and experimental studies, thereby illustrating the power of such combined efforts. The reader is directed to other sources for more in-depth discussions of the complex issues related to PCET processes and to comprehensive reviews that cover more diverse applications.[1−3,5,6,8−11] The first section herein will cover the fundamental concepts of PCET, and the subsequent sections will present examples of PCET in solution, proteins, electrochemistry, and photoexcited processes. Each of the sections on applications will conclude with a discussion of the current challenges from the theoretical and computational viewpoint. The Perspective will conclude with a brief general outlook of the field and a view to the future.

Fundamental Concepts

Theoretical Framework

PCET reactions involving the transfer of one electron and one proton can be described in terms of the four diabatic electronic states depicted in Figure .[15] In the initial state, the electron and proton are on their donors, and in the final state, the electron and proton are on their acceptors. The other two states correspond to either only the proton or only the electron being transferred. Within this framework, a sequential reaction corresponds to moving along the edges of the rectangle, and a concerted reaction corresponds to moving along the diagonal. The mechanism is determined by the relative energies and couplings among these four diabatic states. When the off-diagonal states are much higher in energy, the reaction will be concerted, which is designated EPT for distinct donors and acceptors.

Figure 1

Four diabatic electronic states used in the PCET theory. The sequential mechanism corresponds to following the edges of the rectangle, either PT followed by ET or ET followed by PT, and the concerted mechanism corresponds to following the diagonal, labeled EPT.

The theory for concerted PCET[2,6] combines concepts from Marcus theory for electron transfer[16] and analogous theories for proton transfer.[17] In the simplest case, the reaction is described in terms of the two diagonal states in Figure , where the electron is localized on the donor for the reactant and on the acceptor for the product. As depicted in Figure , the free energy curves for the reactant and product are approximately parabolic along a collective solvent coordinate, which corresponds to the reorganization of the solvent (or protein) environment associated with the charge transfer reaction. The transferring hydrogen nucleus is represented by a quantum mechanical wavefunction. Thus, these free energy curves correspond to electron–proton vibronic states, rather than electronic states as in Marcus theory. As indicated by the proton potential energy curves and the associated ground proton vibrational states in Figure , the proton donor well is lower in energy for the reactant, with the proton localized on the donor side, and the proton acceptor well is lower in energy for the product, with the proton localized on the acceptor side. The energies of the ground proton vibrational states change along the collective reaction coordinate and are identical at the crossing point. Note that this PCET theory is different from the Bixon-Jortner extension of Marcus theory for electron transfer,[18,19] which includes intramolecular harmonic modes that are not coupled to solvent fluctuations. In this PCET theory, the proton motion is not assumed to be harmonic, is of higher frequency, exhibits a larger change of equilibrium position, and is coupled to solvent fluctuations. The effects of the other intramolecular solute modes have been included in this PCET theory with analogous approaches as those used in electron transfer theory.[20]

Figure 2

Free energy curves for the ground reactant (I) and product (II) diabatic electron–proton vibronic states along the collective solvent coordinate for an EPT reaction. The reactant (blue) and product (red) diabatic states correspond to the electron localized on the donor and acceptor, respectively. The proton potential energy curves along the proton coordinate and the corresponding ground state proton vibrational wavefunctions are depicted for the reactant minimum, the crossing point, and the product minimum of the free energy curves. The energies of these proton vibrational states correspond to the open circles on the free energy curves. Adapted with permission from ref (6). Copyright 2008 American Chemical Society.

Within this framework, the general mechanism for a thermal nonadiabatic PCET reaction is as follows: (1) reorganization of the environment leads to the crossing point; (2) a nonadiabatic transition occurs between the reactant and product degenerate vibronic states, corresponding to the simultaneous tunneling of the electron and proton from their donors to their acceptors; (3) further reorganization of the environment stabilizes the product. The simultaneous tunneling in the second step refers to the characters of the initial and final vibronic states during the nonadiabatic transition, where the initial state corresponds to the electron and proton localized on their donors and the final vibronic state corresponds to the electron and proton localized on their acceptors. In general, the excited proton vibrational states must also be considered. Thus, the system should be described by two sets of stacked parabolas corresponding to the different proton vibrational states for the reactant and product diabatic electronic states. For the excited vibronic states, the second step may not correspond to simultaneous tunneling of the electron and proton from donors to acceptors because the excited proton vibrational states may be delocalized.

Two Types of Nonadiabaticity

The issue of nonadiabaticity in the context of PCET reactions has been discussed extensively elsewhere.[6,21,22] This subsection summarizes the main concepts that are essential for motivating the form of the rate constant expressions given in the next subsection. In the theoretical description of PCET, the system is divided into three subsystems: the solute electrons, the transferring proton(s), and the other nuclei. The electrons and transferring proton are treated quantum mechanically, and the other nuclei are treated classically. The two types of nonadiabaticity are as follows: (1) vibronic nonadiabaticity, which is related to the response of the electron–proton subsystem to motion of the other nuclei, and (2) electron–proton nonadiabaticity, which is related to the response of the electrons to motion of the transferring proton.

The vibronic nonadiabaticity is characterized mainly by the vibronic coupling, which is the Hamiltonian matrix element between the electron–proton vibronic wavefunctions associated with the reactant and product. When this vibronic coupling is much less than the thermal energy, and several other criteria are satisfied,[23] the reaction is vibronically nonadiabatic. In this case, Fermi’s Golden Rule can be used to derive the PCET rate constant, which is proportional to the square of the vibronic coupling for each pair of vibronic states. The form of the vibronic coupling is influenced by the electron–proton nonadiabaticity. When proton transfer is electronically nonadiabatic, the vibronic coupling is the product of the electronic coupling and the overlap integral of the reactant and product proton vibrational wavefunctions. This regime is typically the most relevant to EPT reactions and will be the focus of the next subsection on rate constant expressions.

Electron–proton nonadiabaticity is characterized by the relative time scales of the electrons and the proton or by the nonadiabatic coupling matrix element defined in electronic structure theory.[13,14] A semiclassical formulation[24] can be used to estimate the effective proton tunneling and electronic transition time scales. When the electrons are faster than the proton, they respond instantaneously to the proton motion, and the reaction occurs on the electronic ground state and is electronically adiabatic. In contrast, when the electrons are not able to respond fast enough to the proton motion, the excited electronic states participate in the reaction, which is therefore electronically nonadiabatic. To complement the semiclassical formulation, electronic structure calculations can be used to determine the component of the nonadiabatic coupling between the ground and excited electronic states along the proton transfer coordinate. This nonadiabatic coupling is large when the electronic wavefunction changes significantly and abruptly with respect to the proton motion. Thus, it can be viewed as a measure of the change in electronic charge distribution as the proton transfers. The reaction is electronically nonadiabatic when this nonadiabatic coupling (i.e., the change in electronic charge distribution upon proton transfer) is large. This physical behavior can also be identified by evaluating the change in dipole moment, electrostatic potential, or partial atomic charges along the proton transfer coordinate.

These diagnostics for electron–proton nonadiabaticity have been applied to the phenoxyl-phenol and benzyl-toluene self-exchange reactions[12] and have illustrated that proton transfer is electronically nonadiabatic in the former and electronically adiabatic in the latter reaction.[13,14] Furthermore, the degree of electron–proton nonadiabaticity also enables a quantitative distinction between HAT and EPT reactions: HAT reactions are electronically adiabatic, whereas EPT reactions are electronically nonadiabatic.[13,14] These definitions are consistent with the traditional view that HAT reactions do not involve significant changes in electronic charge distribution (i.e., are associated with the transfer of a neutral hydrogen atom between a single donor and acceptor), in contrast to EPT reactions, which do involve substantial changes in electronic charge distribution because the electron and proton transfer between different donors and acceptors. As discussed above, the donors and acceptors may be defined qualitatively in terms of a molecular orbital, a chemical bond, or an atom, although such definitions are not rigorous. The essential distinction between HAT and EPT is that only EPT is accompanied by significant changes in charge distribution. Because the nonadiabatic coupling along the proton transfer coordinate is a measure of the change in electronic charge distribution as the proton transfers, the electron–proton nonadiabaticity can be used to distinguish between HAT and EPT. Note that multisite concerted PCET reactions, in which the electron and proton clearly transfer between distinct sites and often in different directions, are categorized as EPT rather than HAT.

The difference between the phenoxyl-phenol and benzyl-toluene systems has been interpreted in terms of the lone pair of electrons on the former that is lacking on the latter system.[12] For the transition state geometries analyzed in this manner, the singly occupied molecular orbital is dominated by 2p orbitals perpendicular to the proton donor–acceptor axis for the phenoxyl-phenol system but is dominated by atomic orbitals oriented along the proton donor–acceptor axis for the benzyl-toluene system. These differences in geometry and electronic structure result in the differences in the degree of electron–proton nonadiabaticity. Interestingly, another transition state structure for the phenoxyl-phenol system[25] has a proton transfer interface more similar to that of the benzyl-toluene system and has been characterized as electronically adiabatic, occurring via the HAT mechanism.

Rate Constants and Kinetic Isotope Effects

PCET rate constant expressions have been derived in various well-defined limits. Typically concerted PCET reactions are vibronically nonadiabatic due to the small vibronic coupling. The simplest nonadiabatic PCET rate constant expression is[20]where the summations are over reactant and product vibronic states, Pμ is the Boltzmann probability for the reactant state μ, Vel is the electronic coupling, Sμν is the overlap between the reactant and product proton vibrational wavefunctions for states μ and ν, λ is the reorganization energy, ΔGμν0 is the reaction free energy for states μ and ν, and kB is Boltzmann’s constant. The reorganization energy and reaction free energy are depicted in Figure , and the overlap is depicted at the crossing point. In eq the proton transfer is assumed to be electronically nonadiabatic (i.e., significant electron–proton nonadiabaticity), thereby leading to this specific form of the vibronic coupling. Thus, this expression is valid for EPT but not necessarily for HAT reactions. Alternative expressions have been derived for other regimes.[17,22]

The proton vibrational wavefunction overlap plays an important role in determining the rate constants and kinetic isotope effects (KIEs) of PCET reactions. Moreover, this overlap depends strongly on the proton donor–acceptor distance: the overlap is larger for shorter distances. Thus, the rate constant increases as the proton donor–acceptor distance decreases. The KIE, defined as the ratio of the rate constant for hydrogen transfer to the rate constant for deuterium transfer, is proportional to the square of the ratio of the hydrogen and deuterium overlap integrals for a given pair of vibronic states:where SH and SD are the overlaps of the hydrogen and deuterium wavefunctions, respectively. The overlaps decrease for both hydrogen and deuterium as the proton donor–acceptor distance increases, but the deuterium overlap falls off faster because of its larger mass. As a result, the ratio of the hydrogen to deuterium overlap increases as the proton donor–acceptor distance increases for a given pair of vibronic states. The relation between the KIE and the proton donor–acceptor distance is not always straightforward, however, because of complexities due to different contributions from excited vibronic states.[6,22,26,27] According to this analysis, the rate constant decreases and the KIE often increases as the proton donor–acceptor distance increases. This trend will be illustrated by the PCET reaction between ruthenium polypyridyl complexes discussed in the next section.

Given the importance of the proton vibrational wavefunction overlap and its strong dependence on the proton donor–acceptor distance, the PCET theory has been expanded to include the proton donor–acceptor motion.[28] The proton donor–acceptor mode is characterized by M and Ω, its effective mass and frequency, respectively, where MΩ2 is the force constant. Rate constant expressions have been derived in various well-defined limits. When the energy associated with the proton donor–acceptor motion, ℏΩ, is similar to or lower than the thermal energy kBT, the following rate constant expression is applicable:[28]Here αμν is an attenuation parameter that reflects the exponential decrease of the proton vibrational wavefunction overlap with the proton donor–acceptor distance. Comparing this expression to eq , the only additional term arising from inclusion of the proton donor–acceptor motion in this regime is exp[2kBTαμν2/(MΩ2)].

An alternative form of the rate constant that avoids the approximation of an exponential decrease in the overlap is obtained by a thermal averaging procedure, which assumes that an equilibrium thermal distribution is maintained. In this approach, the rate constant in eq is calculated for different values of the proton donor–acceptor distance R and weighted by the probability P(R) of sampling that value of R, followed by integration over all R. Thus, this rate constant is given byThe rate constants in eqs and 4 have been shown to be mathematically identical in certain well-defined regimes.[29]

When only the ground reactant and product vibronic states are included in the rate constant expression given by eq , the KIE can be approximated as[30]where SH and SD are the overlaps of the hydrogen and deuterium wavefunctions, respectively, at the equilibrium proton donor–acceptor distance, and αH and αD represent the exponential attenuation parameters for hydrogen and deuterium, respectively. Note that the inclusion of the proton donor–acceptor motion leads to the temperature dependence of the KIE, as has been observed experimentally.[31−33] If the proton donor–acceptor distance is assumed to be fixed, as in eq , the KIE is simply the ratio of the squares of the overlaps, as in eq .

PCET in Solution and Proteins

Ruthenium Bipyridyl Complexes in Solution

Numerous PCET reactions in solution have been studied experimentally and theoretically. A particularly illustrative example is PCET between ruthenium polypyridyl complexes, as depicted in Figure for systems labeled CompA and CompB. This reaction is presumed to occur by a concerted EPT mechanism because the single ET and single PT reactions are significantly endoergic, while the EPT mechanism is slightly exoergic.[34] Thus, the concerted mechanism avoids high-energy intermediates. The experimental data indicate that the rate constant is 9.6 times larger for CompB than for CompA, but the KIE is 11.3 for CompB and 16.1 for CompA.[35] Density functional theory (DFT) was used to optimize the geometries of the acceptor complexes for both CompA and CompB.[34] These calculations indicated significantly more steric crowding near the acceptor oxygen atom for CompA than for CompB. In other words, changing the bipyridine ligand to a tripyridine ligand pulled the ligands away from the acceptor oxygen atom in CompB. As a result, the donor complex can get closer to the acceptor complex for CompB, leading to a smaller proton donor–acceptor distance by ∼0.06 Å for CompB than for CompA. According to the analysis in the previous section, a shorter proton donor–acceptor distance tends to lead to a larger rate constant and a smaller KIE, as observed experimentally for these systems.

Figure 3

PCET reaction in Ru bipyridyl complexes, where the electron transfers between the two Ru centers and the proton transfers from the water ligand to the oxygen ligand. Experiments indicate that the rate constant is 9.6 times faster for CompB than for CompA, and the KIE is 16.1 for CompA and 11.3 for CompB. The PCET theory explains these differences in terms of a longer O---O distance for CompA than for CompB. Reproduced with permission from ref (6). Copyright 2008 American Chemical Society.

Soybean Lipoxygenase

PCET also plays a key role in many different types of proteins. A well-studied example is the enzyme soybean lipoxygenase (SLO), which catalyzes the PCET reaction from the linoleic acid substrate to the iron cofactor.[33] Similar to the ruthenium bipyridyl complexes, the single electron and single proton transfer reactions are highly endoergic, whereas the concerted EPT reaction is slightly exoergic, thereby favoring the concerted mechanism.[36] Moreover, analysis of the orbitals obtained from DFT calculations[37] indicates that the electron transfers from the π-backbone of the linoleic acid to the iron, and the proton transfers from C11 of the linoleic acid to the oxygen of the OH ligand, as depicted in Figure . The KIE was measured experimentally to be ∼80 at room temperature and to exhibit relatively weak temperature dependence.[33] Because of this unusually large KIE upon replacement of the transferring hydrogen with deuterium, this system has been the subject of many theoretical studies.[30,33,36−40]

Figure 4

Schematic representation of the net hydrogen atom transfer catalyzed by SLO with the linoleic acid substrate. This process is thought to occur via an EPT mechanism. The red arrow indicates the electron transfer from the π-backbone of the linoleic acid substrate to the iron of the cofactor, and the blue arrow indicates the proton transfer from C11 of the substrate to the iron-bound hydroxide to form water. Reproduced with permission from ref (22). Copyright 2014 American Chemical Society.

An issue of contention among these theoretical studies has been whether this reaction is nonadiabatic.[22] Recently this issue has been studied extensively using the quantitative diagnostics for the two types of nonadiabaticity described above.[41] Constrained DFT was used to calculate the proton potential energy curves associated with the diabatic electronic states, as well as the electronic coupling between these two electronic states, for a model of the SLO active site. The effective electronic transition time was found to be a factor of 85 greater than the effective proton tunneling time, and the dipole moment was shown to change dramatically along the proton transfer coordinate. As discussed in the previous section, these observations are signatures of electron–proton nonadiabaticity. Moreover, the vibronic coupling was found to be significantly less than the thermal energy, signifying vibronic nonadiabaticity. This evidence for these two types of nonadiabaticity, in conjunction with the relatively low frequency of the proton donor–acceptor motion (i.e., the C---O motion), validates the use of the nonadiabatic rate constant expression given by eq .

Several different varieties of the nonadiabatic PCET theory were applied to SLO.[30,33,36,38] In an early study, the protein and solvent environments were described by a dielectric continuum, and two parameters were fit to experimental data.[36] In a later study, the explicit protein and solvent were included, and the input quantities to the rate constant expression were obtained from classical molecular dynamics simulations and DFT calculations on a model system.[30] In both cases, the theoretical calculations reproduced the magnitude and temperature dependence of the KIE. The unusually high KIE was found to arise from the relatively small overlap between the reactant and product proton vibrational wavefunctions and the dominance of the lowest-energy vibronic states. The small overlap is due to the weak C–H---O hydrogen-bonding interaction, leading to a relatively large equilibrium C---O distance. On the basis of eq , the PCET theory predicts that the magnitude of the KIE will increase as the equilibrium C---O distance increases because the ratio of squared overlaps increases. Moreover, the temperature dependence of the KIE is predicted to increase as the frequency of the C---O motion decreases due to the behavior of the temperature-dependent exponential factor.

These predictions were verified by experimental studies in which Ile553, which is ∼15 Å from the iron, was mutated to less bulky residues.[42] The magnitude and temperature dependence of the KIE were found to increase as residue 553 became less bulky. Utilizing the full expression given by eq , the equilibrium proton donor–acceptor distance and associated frequency were fit to the experimental data while all other parameters were held fixed to the values determined for wild-type SLO.[43] According to these calculations, the equilibrium proton donor–acceptor distance increases and the associated frequency decreases as residue 553 becomes less bulky, leading to the increase in the magnitude and temperature dependence of the KIE.

More recently, experimental measurements revealed that the double mutant Leu546Ala/Leu754Ala has an enormous KIE of 500–700 at room temperature.[44] The crystal structure of the mutant showed a predominantly unaltered backbone with a slightly expanded active site cavity. Utilizing eq , this large KIE could be reproduced by increasing the equilibrium proton donor–acceptor distance by 0.1–0.2 Å and retaining a similar frequency for the associated motion as for the wild-type enzyme.[44] These calculations suggest that this enormous KIE is observed because the CH---O interface is constrained to configurations with poor hydrogen vibrational wavefunction overlap. For the Ile553 mutants, the proton donor–acceptor mode frequency decreased as the equilibrium distance increased, thereby allowing the system to effectively sample the shorter distances associated with more moderate KIEs. For the double mutant, however, the frequency did not decrease even though the equilibrium distance increased, and the system was unable to effectively sample the shorter distances, leading to the colossal KIE.

Current Challenges for PCET in Solution and Proteins

The analytical rate constant expressions for PCET provide conceptual understanding and generate predictions that can be tested experimentally. Electronic structure calculations and classical molecular dynamics simulations can be used to obtain the input quantities for these analytical expressions. However, in some cases the level of accuracy provided by classical molecular dynamics simulations is not sufficient for reproducing the subtle changes in bond lengths (∼0.1 Å) and vibrational frequencies (∼50 cm–1) due to chemical modifications or mutations that can lead to substantial changes in the rates, KIEs, and their temperature dependences. The main limitations are the accuracy of the molecular mechanical force fields, even if a portion of the system is treated quantum mechanically, and the extent of conformational sampling. To avoid such limitations in computer simulation methods, the parameters in the analytical rate constant expressions can be fit to the experimental data to provide understanding and enable predictions in trends, as described above. A significant challenge in this field is to develop computational methods that are capable of describing these subtle changes upon relatively minor chemical modifications or distal mutations of large, condensed phase systems, particularly complex biological systems.

Another challenge is the development of dynamical methods that include both electronic and nuclear quantum effects and can directly simulate thermal PCET reactions with the level of required accuracy. Ring polymer molecular dynamics (RPMD) has been applied to PCET in solution[45] and provides useful visualizations of the reactions. These simulations have also reconfirmed information such as the concerted nature of the ruthenium bipyridyl reaction,[35] which was ascertained previously from the thermodynamics of the single ET and PT reactions compared to the EPT reaction.[34] On the other hand, RPMD is also subject to the limitations of accurate potential energy surfaces and adequate conformational sampling, as well as problematic issues related to the inverted Marcus region,[46,47] which is particularly important for PCET reactions. Surface hopping molecular dynamics with the proton treated quantum mechanically using grid-based methods is a viable approach for simulating thermal PCET reactions in the normal and inverted Marcus regimes.[22]

A more general challenge arises in the study of PCET systems that are less well-defined than a system such as lipoxygenase. For example, the donors and acceptors for the proton and electron may be unknown, making it difficult to define the relevant diabatic states given in Figure . Moreover, for systems such as those involved in photosynthesis, many electrons and protons are transferring simultaneously in different directions at various locations, and understanding how all of these reactions are coupled together is challenging. The theoretical framework described above has been extended to processes involving multiple electron and/or proton transfer reactions,[15,48] but the complexity increases, and the analysis becomes significantly more challenging.

Electrochemical PCET

PCET in Molecular Electrocatalysts

PCET also plays a critical role in electrochemical processes, which are relevant to a wide range of energy devices. This section focuses on theoretical studies of molecular electrocatalysts, although analogous methods and principles apply to heterogeneous electrocatalysts. Over the past several years, combined theoretical and experimental studies have guided the design of more active molecular electrocatalysts. The overall objective of these types of studies is to design catalysts with high turnover frequency and low overpotential, preferably composed of environmentally friendly, cost-effective, and earth-abundant materials. The catalytic cycles are comprised of a series of steps that are often proton transfer (PT), electron transfer (ET), or concerted (EPT) reactions. Both the thermodynamics, namely the relative free energies of the intermediates along the cycle, and the kinetics, namely the free energy barriers connecting these intermediates, must be considered in catalyst design. Modifying the catalysts can change the mechanism (i.e., the order of the steps), as well as the thermodynamics and kinetics along the reaction pathway.

A variety of different protocols have been devised to calculate the reduction potentials and pKa’s of molecular electrocatalysts.[49−57] Typically these protocols are based on DFT geometry optimizations combined with a polarizable continuum model to calculate the solvation free energies. The most reliable strategy is to calculate the reduction potentials or pKa’s relative to a related reference reaction for which experimental data are available.[52,53,55,58] This strategy, often discussed in terms of isodesmic reactions, avoids the necessity of calculating the electrode potential, the free energies of the solvated electron and proton, and changes in standard states, thereby eliminating associated systematic errors. Additional errors arising from the DFT functional, basis set, and solvation model also tend to cancel out in this approach. The reduction potentials and pKa’s provide the reaction free energies for the associated ET and PT reactions, respectively, and therefore enable the generation of the free energy pathway for any proposed mechanism. Standard transition state calculations provide the free energy barriers for the PT steps, and various approaches have been developed for calculating the inner-sphere (solute) and outer-sphere (solvent) reorganization energies for the ET and EPT steps.[59,60]

Many theoretical studies of molecular electrocatalysts have focused on sequential mechanisms, where ET and PT occur in a series of separate steps. However, efforts to reduce the overpotential requirement have also focused on designing catalysts that favor the concerted mechanism, which avoids high-energy intermediates and therefore tends to be associated with a lower overpotential. The concerted PCET (EPT) theory described above has been extended to electrochemical PCET.[6,61] In this case, the electron transfers between the electrode and a molecule or hydrogen-bonded complex in solution, and the proton is assumed to transfer within this complex. The electrochemical rate constants are of a similar form as the homogeneous rate constants given above, except they require integration over the electronic energy levels of the electrode with appropriate weighting by the Fermi distribution and density of states. This theory for electrochemical EPT has been applied to various molecular electrocatalysts.[62] Note that this theory differs from that developed by Savéant and co-workers[63,64] in terms of the prefactor and the reorganization energy in the derived rate constant expressions, as well as the treatment of the proton donor–acceptor motion. A more detailed comparison of these two theoretical treatments is provided in ref (61).

Cobaloximes and Cobalt Dithiolenes

Cobaloximes, depicted in Figure , have been shown to produce H2 from protic solutions at modest overpotentials.[65−67] Various reaction pathways by which H2 is evolved monometallically or bimetallically from a Co(III)H or a Co(II)H intermediate have been proposed.[53,67−69] The reduction potentials and pKa’s for each step of the proposed mechanisms have been calculated with a well-defined protocol.[53] These calculations allowed the construction of the free energy diagrams for proposed mechanisms and the determination of the thermodynamically favored pathways.

Figure 5

Top: Structure of a cobaloxime, where L is a solvent molecule and R is a substituent that was varied to investigate linear correlations. Bottom: Calculated reduction potentials and pKa’s as functions of the Hammett constants for a series of substituents: R = −NH2, −OH, −OCH3, −CH3, −C6H5, −H, −Cl, −CF3, −CN. Reproduced in part with permission from ref (70). Copyright 2011 American Chemical Society.

The impact of altering the substituents on the cobaloximes was investigated by performing these calculations for a series of substituents characterized by the Hammett constant, which reflects the electron-donating or electron-withdrawing character.[70] All of the reduction potentials and pKa’s were found to be linearly correlated with the Hammett constant, as depicted in Figure . Such linear correlations are valuable because all of these quantities can be determined for any substituent if the Hammett constant is known, and if the Hammett constant is not known, then only one of these quantities must be calculated or measured to obtain the rest of them. Knowledge of these quantities enables the construction of the free energy diagram associated with any proposed mechanism.

Figure indicates that the linear plots corresponding to the reduction potentials for the Co(II/I) and Co(III/II)H couples intersect at the Hammett constant associated with the methyl substituent, suggesting that the cyclic voltammogram (CV) peaks associated with these two couples would overlap. This observation led to the reassignment of the CV peak at ca. −1.0 V vs SCE in acetonitrile to the Co(II/I)H couple[53,68] rather than the Co(III/II)H couple[66] because the calculations suggested that the Co(III/II)H peak was obscured by the Co(II/I) peak at −0.55 V vs SCE in acetonitrile. Moreover, Figure leads to the prediction that these overlapping peaks will separate as the Hammett constant becomes more positive.[70] This example illustrates that DFT can assist in the assignment of CV peaks. Moreover, these assignments have mechanistic implications in terms of whether the monometallic or bimetallic pathway is thermodynamically favored.[57,69] In some cases, the formation of nanoparticles on the electrodes may complicate the mechanistic interpretation.[71]

In addition to the substituents, the ligands of the cobaloximes can be modified. In particular, the BF2 bridge has been replaced by a H bridge.[72,73] In these complexes, the H bridge can be protonated, which reduces the required overpotential. Although the BF2 bridge is not as easily protonated, the BF2 bridge is more electron-withdrawing than the H bridge, also reducing the required overpotential. Thus, minimizing the required overpotential entails a balance of these two effects, and the outcome depends on the identity of the metal center. According to DFT calculations,[74] a single H bridge is favorable for cobalt and nickel centers, whereas two BF2 bridges are more favorable for an iron center. Ligand protonation could also occur at the nitrogen and, in some cases, could lead to decomposition of the catalyst.

Ligand protonation has also been shown to play an important role in cobalt dithiolenes, as depicted in Figure . Experimental studies showed that the electrocatalytic overpotentials do not behave as expected in terms of the electron-withdrawing character of a series of ligands.[75] DFT calculations illustrated that one or two sulfur atoms can become protonated, as shown in Figure , thereby explaining the anomalous trend in the catalytic potentials.[76] The mechanism for H2 evolution in these catalysts has been proposed to involve proton transfer from the sulfur atom to the cobalt center, forming an active cobalt-hydride species that could produce hydrogen with an acid or with another protonated sulfur ligand.[76] This example further highlights the importance of ligand protonation in electrocatalysis.

Figure 6

Structure of two different cobalt dithiolene catalysts with one or two of the sulfur atoms protonated. The color scheme is as follows: magenta, cobalt; yellow, sulfur; cyan, carbon; blue, nitrogen; white, hydrogen. In the mechanism for hydrogen evolution, the proton may transfer from the sulfur to the cobalt, possibly after a thermodynamically accessible isomerization to orient the proton toward the cobalt. Reproduced in part with permission from ref (76). Copyright 2012 American Chemical Society.

Hydrogenase Models

A variety of biomimetic catalysts modeled after hydrogenase enzymes, which catalyze both H2 oxidation and production, have been studied with theoretical methods.[52,55,62,77−81] Often a key step in the catalytic cycle is protonation to form a metal-hydride species. A combined experimental and theoretical study of the Ni–Fe molecular electrocatalyst depicted in Figure illustrated the significance of isomerization at the metal centers in these types of catalysts.[81] Initially the Ni(I)Fe(I) complex[82] and the oxidized[83] and protonated[84] species were characterized with X-ray crystallography, as well as infrared spectroscopy to determine the CO vibrational mode frequencies. These experiments indicated that the geometry at the Ni site is tetrahedral for the Ni(I)Fe(I) species but square planar for the protonated Ni(II)HFe(II) and oxidized Ni(II)Fe(I) species. Similar results were obtained for complexes with Pd substituted for Ni; however, the infrared spectroscopy on complexes with Pt substituted for Ni indicated that the geometry at the Pt site of the neutral species is square planar with properties of Pt(II)Fe(0).

Figure 7

Structure of the Ni–Fe catalyst and illustration of the isomerization process at the Ni center that is proposed to occur prior to protonation on the basis of theoretical calculations and experimental measurements. The thermodynamically accessible isomerization from the tetrahedral to square planar geometry at the Ni site is accompanied by electron transfer from Ni to Fe, thereby enhancing the basicity of the Fe site by ∼108 and facilitating the protonation step.

DFT calculations provided further insights into these systems.[81] For the neutral Pt catalyst, DFT calculations confirmed that the isomer with a square planar geometry at the Pt site is ∼30 kcal/mol lower than the isomer with a tetrahedral geometry at the Pt site. Interestingly, DFT calculations on the Ni and Pd catalysts revealed a previously undetected isomer with a square planar geometry at the Ni or Pd site. In these systems, the minimum associated with the square planar isomer is nearly isoergic with the minimum associated with the tetrahedral geometry at the metal center, and the free energy barrier for isomerization from the tetrahedral to the square planar geometry was ∼7 kcal/mol for the Ni center and ∼3 kcal/mol for the Pd center. Supporting these calculations, 31P NMR experiments implicated transient or intermediate square planar species interconverting the 31P sites for the Ni and Pd catalysts with barriers that are consistent with the DFT results. Thus, both the DFT calculations and the 31P NMR experiments indicate that the previously unobserved square planar Ni(II)Fe(0) isomer is thermodynamically accessible.

The isomerization from tetrahedral to square planar geometry at the Ni site is associated with electron transfer from the Ni to the Fe center, thereby weakening the metal–metal bond and enhancing the basicity of the Fe metal site. Quantitatively, DFT calculations illustrated that the square planar species is vastly more basic than the tetrahedral species by a factor of ∼108 (i.e., the pKa is ∼8 units higher for the square planar than for the tetrahedral isomer). The results of this combined experimental and theoretical study represented a paradigm shift in the mechanistic interpretation of these catalysts.[81] As depicted in Figure , this work suggested that protonation occurs via a previously undetected square planar Ni(II)Fe(0) isomer with enhanced basicity to facilitate protonation at the Fe site. The catalytic importance of two-electron mixed-valence species that may be slightly higher in free energy but nevertheless are the active species is relevant to other molecular electrocatalysts as well.

The Ni(P2N2)2 catalysts depicted in Figure exhibit significantly higher turnover frequencies for hydrogen evolution due to the presence of the pendant amine, which serves as a proton relay.[85−87] Similar to the cobaloximes, the reduction potentials and pKa’s of the possible ET and PT steps have been calculated and used to generate free energy pathways, Pourbaix diagrams, and in some cases the complete thermodynamic cycle.[52,55,79] The effects of modifying the substituents on the nitrogen and phosphorus have also been investigated, and linear correlations between various properties have been identified.[80]

Figure 8

Structure of the Ni(P2N2)2 catalyst. The substituents on the nitrogen and phosphorus groups are not shown but were varied to investigate linear correlations. The PCET step shown here involves proton transfer (PT) from the nitrogen of the pendant amine to the Ni and electron transfer (ET) from the electrode to the molecule.

Furthermore, to reduce the overpotential requirement, efforts were aimed at designing catalysts that favor the concerted PCET mechanism (i.e., the EPT mechanism).[62] These efforts focused on the steps in the catalytic cycle involving PT between Ni and N and ET between the complex and the electrode, as depicted in Figure . The expression given in eq was used to calculated the EPT rate constant for the Ni(P2N2)2 catalyst with methyl substituents. In this system, the proton transfers between the Ni and the N of the pendant amine. The equilibrium Ni---N distance was calculated to be 3.25 Å, which is unfavorably long for a PT reaction. However, the catalyst undergoes thermal fluctuations that decrease this distance to facilitate PT. The calculations predicted that the EPT rate constant will increase as the equilibrium Ni---N distance decreases and as the amine ligand becomes more flexible to facilitate contraction of this distance with a lower energy penalty.

On the basis of this prediction, several Ni catalysts with more flexible pendant amines were examined.[88] However, often the catalysts with more flexible amines (i.e., with a lower frequency associated with the Ni---N motion) exhibited larger equilibrium Ni---N distances, thereby counteracting the advantage, although other factors such as inner-sphere reorganization energy could favor the flexible amine ligands. Thus, the EPT mechanism requires a balance between a well-positioned pendant amine with a short equilibrium Ni---N distance and a flexible pendant amine that facilitates a further decrease in this distance. The overall guiding design principle generated by these calculations was that the pendant amines should be flexible enough to allow motion toward the Ni center but still reasonably well-positioned.[88]

Current Challenges for Electrochemical PCET

A major challenge for calculations of electrochemical PCET is the description of the explicit ions, solvent molecules, and catalysts at the electrode surface. Most current methods for calculating reduction potentials are based on a dielectric continuum description of the solvent, thereby neglecting the effects of explicit solvent molecules interacting with the molecule, although a small number of explicit solvent molecules can be included in the electronic structure calculations. Furthermore, most of these dielectric continuum methods for calculating reduction potentials neglect the effects of the electrode and the electrolyte ions. Recently developed methods for calculating the electrochemical solvent reorganization energy[60] include the effects of the electrode but neglect the effects of ions at the interface and in the bulk solvent. The calculation of accurate electronic couplings between a molecular and an electrode, particularly in the presence of solvent and electrolyte ions, is also still an unsolved problem. Thus, the development of computational methods that provide an accurate description of the ions, solvent molecules, and catalysts at the interface in the context of calculating reduction potentials, reorganization energies, and couplings is a key challenge in the field.

Another significant challenge is the investigation of heterogeneous catalysis, in which the surface plays a direct role in a chemical step of the catalytic cycle. For example, surface atoms could accept and/or donate protons in a catalytic step. The description of such heterogeneous processes requires an accurate, atomic-level, quantum mechanical description of the electrode surface as well as the molecular catalyst. Moreover, the challenges associated with describing the solvent molecules and ions at the interface are also significant in heterogeneous catalysis. The optimal methods may combine a quantum mechanical treatment of the atoms directly participating in chemistry, an explicit molecular mechanical treatment of nearby molecules, a continuum treatment of the bulk solvent and/or electrode, and analytical expressions to provide understanding and predictive power.

Photoinduced PCET

In systems such as photosynthetic reaction centers and solar cells, PCET is induced by light. Understanding the nonequilibrium dynamics of a condensed phase PCET system following photoexcitation is another important direction of research. Photoinduced PCET is different from the more common excited state proton transfer, in which photoexcitation alters the electronic charge distribution predominantly in the proton transfer interface region to induce proton transfer. In contrast, photoinduced PCET involves a veritable electron transfer reaction as well as a proton transfer reaction upon photoexcitation. Moreover, the photoinduced PCET discussed in this section is also distinct from excited state PCET, in which photoexcitation induces charge separation to prepare the system for a subsequent PCET process.[89−93] The analytical rate constant expressions discussed above are not applicable to nonequilibrium photoinduced PCET because such expressions are based on the assumption that the system is initially at equilibrium. Typically photoexcitation induces an instantaneous change in the electronic charge distribution of the solute; consequently, the system is no longer at equilibrium. Simulation of the real-time nonequilibrium dynamics of the solute and the solvent following photoexcitation requires the use of nonadiabatic molecular dynamics methods to allow relaxation of the system from the excited electronic state down to the ground electronic state.

The surface hopping molecular dynamics with quantum transitions method was developed by Tully for electronic surfaces[94] and subsequently extended for proton vibrational surfaces.[95] More recently it was extended for electron–proton vibronic surfaces in the context of photoinduced PCET reactions.[96,97] In this approach, the classical nuclei move on a single surface except for instantaneous transitions incorporated according to Tully’s fewest switches algorithm.[94] For an ensemble of trajectories, the fraction of trajectories on each surface at each time is approximately equivalent to the quantum probability obtained by integrating the time-dependent Schrödinger equation. This approach has been used to study photoinduced PCET in model systems that have exhibited concerted, sequential, and complex branching pathways.[96,97]

Recently, this nonadiabatic dynamics method was used to simulate photoinduced PCET in a hydrogen-bonded phenol-amine complex that had been studied experimentally.[98,99] The experiments implicated two different mechanisms in this system, as depicted in Figure . The first mechanism was sequential, namely ET followed by PT, where photoexcitation to the intramolecular charge transfer (ICT) state was followed by PT on a lower excited state. The second mechanism was concerted, where photoexcitation directly to the EPT state was followed by relaxation within this state. For the EPT mechanism, Raman experiments[98] indicated that the electronic charge distribution at the hydrogen-bonding interface shifted upon photoexcitation: the electronic density at the O–H bond shifted to the N–H bond, although the hydrogen nucleus did not move on this time scale, thereby leading to an elongated N–H bond that subsequently relaxed.

Figure 9

Top left: Experimentally studied hydrogen-bonded complex composed of p-nitrophenylphenol and tert-butylamine. For the nonadiabatic dynamics simulations, tert-butylamine was replaced with ammonia. Bottom left: Schematic depiction of EPT (S1) and ICT (S2) excited state potential energy surfaces as functions of the H transfer coordinate and a collective reaction coordinate. The ground state (S0) is not shown. Right: Representative MDQT trajectory initiated on the ICT (S2) state. The decay from the S2 to the S1 state at ∼164 fs is followed by proton transfer from O to N on the S1 state at ∼600 fs. Upon decay to the S0 state at ∼745 fs, the proton transfers back to O. The O–H and N–H distances are calculated for the proton that is involved in the hydrogen bond. Reproduced in part with permission from ref (99). Copyright 2015 American Chemical Society.

In the initial simulations of this photoinduced PCET process, the transferring hydrogen nucleus was treated classically, and surface hopping molecular dynamics trajectories were propagated on the S0, S1, and S2 electronic state surfaces.[99] The potential energy surfaces were obtained with a quantum mechanical/molecular mechanical (QM/MM) method, where the solute was treated quantum mechanically and the surrounding solvent was treated with a molecular mechanical force field. The solute electronic states were generated on-the-fly with a semiempirical implementation of the floating occupation molecular orbital complete active space configuration interaction (FOMO-CASCI) method.[100,101] This multiconfigurational method inherently includes non-dynamical electron correlation but also includes dynamical electron correlation by fitting some of the semiempirical parameters to data from complete active space second-order perturbation theory (CASPT2) calculations.

Prior to the nonadiabatic dynamics simulations, the electronic states were characterized for this system using CASPT2 in the gas phase and QM/MM FOMO-CASCI molecular dynamics in solution.[99] In the gas phase, the minimum energy structure of the phenol-amine complex for all three electronic states corresponds to the proton bonded to the oxygen. In solution, however, the free energy profiles along the proton transfer coordinate exhibit different behavior. For the S0 and S2 states, the most thermodynamically stable configuration still corresponds to the proton bonded to the oxygen, with a substantial free energy barrier for PT to the nitrogen. For the S1 state, however, the most thermodynamically stable configuration corresponds to the proton bonded to the nitrogen, and the free energy barrier for PT from the oxygen to the nitrogen is only ∼4 kcal/mol, which is easily surmountable through zero point energy effects. Thus, these calculations are consistent with the experimental interpretation of S1 as an EPT state and S2 as an ICT state.

In the nonadiabatic dynamics simulations, ∼230 trajectories were initiated on the S1 state and on the S2 state. The system decayed from the S2 to the S1 state in ∼100 fs, and the decay time scale from the S1 to the S0 state was ∼0.9 ps. These time scales are in qualitative agreement with the experimental transient absorption data.[98] In addition, 54% of the trajectories exhibited PT on the S1 state. Figure depicts this PT for a representative trajectory following photoexcitation to the S2 state. After initial fast decay from S2 to S1, the system exhibited PT from the oxygen to the nitrogen on the S1 state, followed by PT from the nitrogen back to the oxygen upon decay from S1 to S0. Thus, the simulations provided atomic-level evidence of PT on the EPT state, which was suggested by the experiments but could not be detected experimentally. Recent analysis of these simulations highlighted the significant role of solvent dynamics in this photoinduced PCET process. Current simulations are treating the transferring hydrogen nucleus quantum mechanically to incorporate vibrational relaxation effects. Overall, these types of simulations enable the investigation of nonequilibrium solute, solvent, charge transfer, and vibrational relaxation dynamics.

Current Challenges for Photoinduced PCET

A critical challenge for simulating photoinduced PCET is the efficient generation of accurate excited state potential energy surfaces. Time-dependent DFT (TDDFT) is reasonably efficient but can be problematic for describing charge transfer states[102] and conical intersections,[103] although active research is aimed at addressing these issues.[104−106] On the other hand, ab initio multiconfigurational methods that include dynamical correlation are not computationally practical for on-the-fly dynamics. The semiempirical implementation of FOMO-CASCI[99−101] includes both non-dynamical and dynamical correlation to some extent but requires a fitting procedure that is cumbersome and could lead to unreliable potential energy surfaces in regions that are not part of the fitting procedure. Thus, each of these methods has limitations that need to be considered seriously. Moreover, QM/MM methods are necessary to include the effects of the solvent and/or protein.

Another significant challenge is the incorporation of nuclear quantum effects in a computationally efficient manner. Path integral methods such as RPMD have been applied to thermal PCET[46] but are problematic for photoinduced PCET because of difficulties in describing nonadiabatic or nonequilibrium dynamics on excited electronic states. The grid-based surface hopping methods that have been applied to photoinduced PCET[96,97,99] have been shown to be computationally tractable for the quantum mechanical treatment of a single proton but are not easily extended to a three-dimensional quantum treatment of many protons, which may be required for more complex systems. In addition, surface hopping is not a rigorous method for nonadiabatic dynamics, although many of the key issues, such as decoherence, are not expected to be important for photoinduced processes that decay to the ground state relatively quickly. Moreover, the surface hopping algorithm requires the propagation of a large number of trajectories to ensure convergence, thereby restricting the level of theory used to generate the potential energy surfaces and incorporate nuclear quantum effects. Thus, further developments in all of these areas will be essential for future progress.

General Outlook

The examples described above illustrate that the current theories and computational methods are able to provide useful mechanistic insights and predictions that have been experimentally validated. All of these examples have emphasized the necessity of combining experimental and theoretical approaches to fully understand a given system or process. The current challenges for each type of PCET process have also been discussed in each section. While the existing PCET theories provide a conceptual framework for understanding these processes, further advances in the field will require the development of new computational methods or innovative combinations of existing methods. The collaboration between experimentalists and theoreticians will be critical for unraveling the mysteries of more complex PCET processes.