Multistate Empirical Valence Bond Research Articles

Proton transfer in liquid water confined inside graphene slabs.

The microscopic structure and dynamics of an excess proton in water constrained in narrow graphene slabs between 0.7 and 3.1 nm wide has been studied by means of a series of molecular dynamics simulations. Interaction of water and carbon with the proton species was modeled using a multistate empirical valence bond Hamiltonian model. The analysis of the effects of confinement on proton solvation structure and on its dynamical properties has been considered for varying densities. The system is organized in one interfacial and a bulk-like region, both of variable size. In the widest interplate separations, the lone proton shows a marked tendency to place itself in the bulk phase of the system, due to the repulsive interaction with the carbon atoms. However, as the system is compressed and the proton is forced to move to the vicinity of graphene walls it moves closer to the interface, producing a neat enhancement of the local structure. We found a marked slowdown of proton transfer when the separation of the two graphene plates is reduced. In the case of lowest distances between graphene plates (0.7 and 0.9 nm), only one or two water layers persist and the two-dimensional character of water structure becomes evident. By means of spectroscopical analysis, we observed the persistence of Zundel and Eigen structures in all cases, although at low interplate separations a signature frequency band around 2500 cm^{-1} suffers a blue shift and moves to characteristic values of asymmetric hydronium ion vibrations, indicating some unstability of the typical Zundel-Eigen moieties and their eventual conversion to a single hydronium species solvated by water.

Non-equilibrium reaction and relaxation dynamics in a strongly interacting explicit solvent: F + CD3CN treated with a parallel multi-state EVB model.

We describe a parallelized linear-scaling computational framework developed to implement arbitrarily large multi-state empirical valence bond (MS-EVB) calculations within CHARMM and TINKER. Forces are obtained using the Hellmann-Feynman relationship, giving continuous gradients, and good energy conservation. Utilizing multi-dimensional Gaussian coupling elements fit to explicitly correlated coupled cluster theory, we built a 64-state MS-EVB model designed to study the F + CD3CN → DF + CD2CN reaction in CD3CN solvent (recently reported in Dunning et al. [Science 347(6221), 530 (2015)]). This approach allows us to build a reactive potential energy surface whose balanced accuracy and efficiency considerably surpass what we could achieve otherwise. We ran molecular dynamics simulations to examine a range of observables which follow in the wake of the reactive event: energy deposition in the nascent reaction products, vibrational relaxation rates of excited DF in CD3CN solvent, equilibrium power spectra of DF in CD3CN, and time dependent spectral shifts associated with relaxation of the nascent DF. Many of our results are in good agreement with time-resolved experimental observations, providing evidence for the accuracy of our MS-EVB framework in treating both the solute and solute/solvent interactions. The simulations provide additional insight into the dynamics at sub-picosecond time scales that are difficult to resolve experimentally. In particular, the simulations show that (immediately following deuterium abstraction) the nascent DF finds itself in a non-equilibrium regime in two different respects: (1) it is highly vibrationally excited, with ∼23 kcal mol(-1) localized in the stretch and (2) its post-reaction solvation environment, in which it is not yet hydrogen-bonded to CD3CN solvent molecules, is intermediate between the non-interacting gas-phase limit and the solution-phase equilibrium limit. Vibrational relaxation of the nascent DF results in a spectral blue shift, while relaxation of the post-reaction solvation environment results in a red shift. These two competing effects mean that the post-reaction relaxation profile is distinct from what is observed when Franck-Condon vibrational excitation of DF occurs within a microsolvation environment initially at equilibrium. Our conclusions, along with the theoretical and parallel software framework presented in this paper, should be more broadly applicable to a range of complex reactive systems.

Hydrated proton and hydroxide charge transfer at the liquid/vapor interface of water.

The role of the solvated excess proton and hydroxide ions in interfacial properties is an interesting scientific question with applications in a variety of aqueous behaviors. The role that charge transfer (CT) plays in interfacial behavior is also an unsettled question. Quantum calculations are carried out on clusters of water with an excess proton or a missing proton (hydroxide) to determine their CT. The quantum results are applied to analysis of multi-state empirical valence bond trajectories. The polyatomic nature of the solvated excess proton and hydroxide ion results in directionally dependent CT, depending on whether a water molecule is a hydrogen bond donor or acceptor in relation to the ion. With polyatomic molecules, CT also depends on the intramolecular bond distances in addition to intermolecular distances. The hydrated proton and hydroxide affect water's liquid/vapor interface in a manner similar to monatomic ions, in that they induce a hydrogen-bonding imbalance at the surface, which results in charged surface waters. This hydrogen bond imbalance, and thus the charged waters at the surface, persists until the ion is at least 10 Å away from the interface.

Electron transfer activation of a second water channel for proton transport in [FeFe]-hydrogenase.

Hydrogenase enzymes are important because they can reversibly catalyze the production of molecular hydrogen. Proton transport mechanisms have been previously studied in residue pathways that lead to the active site of the enzyme via residues Cys299 and Ser319. The importance of this pathway and these residues has been previously exhibited through site-specific mutations, which were shown to interrupt the enzyme activity. It has been shown recently that a separate water channel (WC2) is coupled with electron transport to the active site of the [FeFe]-hydrogenase. The water-mediated proton transport mechanisms of the enzyme in different electronic states have been studied using the multistate empirical valence bond reactive molecular dynamics method, in order to understand any role WC2 may have in facilitating the residue pathway in bringing an additional proton to the enzyme active site. In a single electronic state A(2-), a water wire was formed through which protons can be transported with a low free energy barrier. The remaining electronic states were shown, however, to be highly unfavorable to proton transport in WC2. A double amino acid substitution is predicted to obstruct proton transport in electronic state A(2-) by closing a cavity that could otherwise fill with water near the proximal Fe of the active site.

Bringing Reactivity to the Aggregation-Volume-Bias Monte Carlo Based Simulation Framework: Water Nucleation Induced by a Reactive Proton.

The development of the aggregation-volume-bias Monte Carlo based simulation technique has led to recent success in studying rare nucleation events, but thus far, this simulation method has been limited to nonreactive systems. This work presents the first application of this technique to study a reactive system of relevance to atmospheric chemistry, i.e., formation of water droplets in the presence of a reactive proton, by combining this approach with a multistate empirical valence bond (MSEVB) description of the excess proton (or the hydronium). It was shown that the ability for the hydronium to share its charge with adjacent water molecules changes dramatically with the cluster size, especially when clusters are small and the distribution of the charge is affected by the presence of an interface, emphasizing the need to use this more sophisticated MSEVB model for such a reactive system. In addition, the simulation results obtained from this system are compared to those with nonreactive hard-sphere ions of different sizes. Overall, the presence of a hydronium or ions appeared to dramatically change the free energy landscape of nucleation compared to the pure water system, leading to the formation of a stable precritical cluster. Although the free energy change due to the addition of the first few water molecules was shown to be very sensitive to the ionic details, the later portion of the free energy profile was found to be nearly independent of the nature of the ion.

Network analysis of proton transfer in liquid water.

Proton transfer in macromolecular systems is a fascinating yet elusive process. In the last ten years, molecular simulations have shown to be a useful tool to unveil the atomistic mechanism. Notwithstanding, the large number of degrees of freedom involved make the accurate description of the process very hard even for the case of proton diffusion in bulk water. Here, multi-state empirical valence bond molecular dynamics simulations in conjunction with complex network analysis are applied to study proton transfer in liquid water. Making use of a transition network formalism, this approach takes into account the time evolution of several coordinates simultaneously. Our results provide evidence for a strong dependence of proton transfer on the length of the hydrogen bond solvating the Zundel complex, with proton transfer enhancement as shorter bonds are formed at the acceptor site. We identify six major states (nodes) on the network leading from the "special pair" to a more symmetric Zundel complex required for transferring the proton. Moreover, the second solvation shell specifically rearranges to promote the transfer, reiterating the idea that solvation beyond the first shell of the Zundel complex plays a crucial role in the process.

Dynamical aspects of intermolecular proton transfer in liquid water and low-density amorphous ices.

The microscopic dynamics of an excess proton in water and in low-density amorphous ices has been studied by means of a series of molecular dynamics simulations. Interaction of water with the proton species was modelled using a multistate empirical valence bond Hamiltonian model. The analysis of the effects of low temperatures on proton diffusion and transfer rates has been considered for a temperature range between 100 and 298 K at the constant density of 1 g cm(-3). We observed a marked slowdown of proton transfer rates at low temperatures, but some episodes are still seen at 100 K. In a similar fashion, mobility of the lone proton gets significantly reduced when temperature decreases below 273 K. The proton transfer in low-density amorphous ice is an activated process with energy barriers between 1-10 kJ/mol depending of the temperature range considered and eventually showing Arrhenius-like behavior. Spectroscopic data indicated the survival of both Zundel and Eigen structures along the whole temperature range, revealed by significant spectral frequency shifts.

Local and Collective Reaction Coordinates in the Transport of the Aqueous Hydroxide Ion

We investigate local and collective reaction coordinates for the structural diffusion of the hydroxide ion in dilute aqueous NaOH solution using a multistate empirical valence bond (MS-EVB) simulation. We characterize a 15 fs time scale associated with shifting of the equally shared proton within a Zundel-like H3O2(-) ion to form a water molecule, a 550 fs relaxation from this transition state largely guided by electrostatic fluctuations of the surrounding environment, and a 9.6 ps time scale that corresponds to the solvation of the water molecule formed by the proton transfer event. When individual proton transfer events are examined, we are unable to identify a unique transition state solely on the basis of a decrease in the hydroxide ion's coordination number. Instead, we find that the collective electric field along the proton transfer direction is better suited to describe the creation and relaxation of Zundel-like transition states that allow structural diffusion of the hydroxide ion.

On the Origin of Proton Mobility Suppression in Aqueous Solutions of Amphiphiles

Recent experiments reported that proton mobility in tetramethylurea (TMU) solutions is much slower than in urea solutions of the same molarity, and this (as well as the significantly retarded water reorientation) was ascribed to hydrophopic effects. In order to further explore the mechanism of proton transport in these solutions, reactive molecular dynamics simulations using a multistate empirical valence bond model were conducted. The simulations showed that the hydrophobic effect of the TMU methyl groups is weaker than believed. Rather, water concentration is the dominant factor determining proton diffusion. This contrasts with water reorientation and self-diffusion in these samples, which are mutually correlated and depend on the identity of the solute. Interestingly, we find that the mean squared displacements (MSDs) of both water and proton grow nonlinearly in time up to at least 1 ns ("transient subdiffusion"). Subdiffusion is more pronounced for the proton, with an exponent as low as 0.85 that depends, again, on water concentration. Hence, the experimentally relevant long-time diffusivity is observably smaller than what is usually deduced from short simulation runs. It exhibits, for both water and proton, a universal dependence on the power-law exponent.

Molecular Dynamics Simulations of Proton Transport in 3M and Nafion Perfluorosulfonic Acid Membranes

Proton transfer and local structures in 3M (EW 825) and Nafion (EW 890) membranes are investigated in this study by both standard nonreactive molecular dynamics and the self-consistent iterative multistate empirical valence bond method, which is capable of simulating multiple reactive protons and accounting for the Grotthuss mechanism of proton transport. The Nafion and 3M systems have the same backbone, so we can isolate and compare the effect of the different side chains by calculating the radial distribution functions (RDFs), self-diffusion constants, and other properties for three hydration levels at 5, 9, and 14 at 300 and 353 K. The conformations of the 3M and Nafion side chains are also compared. We found that even though many results are similar for both F3C and SPC/Fw water models, certain trends such as the sulfonate clustering can depend on the water model selected. The relationship between the different RDFs for the sulfonate, water, and hydronium is discussed. The self-diffusion constants of water for both membranes are found to be close with respect to each water model selected, even though the experimental values for 3M at 300 K are higher. The calculated self-diffusion constants of the excess protons are found to be higher for 3M than Nafion for hydration levels 9 and 14 at 300 K but statistically the same at 353 K.

Structure and Dynamics of Proton Transfer in Liquid Imidazole. A Molecular Dynamics Simulation

Proton transfer (PT) via the Grotthuss mechanism in liquid imidazole (im) at 393 K is studied with molecular dynamics simulation using a reactive multistate empirical valence bond (MS-EVB) model. It is found that the proton is tightly binded to an imidazole to form an imidazolium (imH(+)), which is solvated in a distorted Eigen-like complex (im-imH(+)-im), whereas the Zundel-like complex (im-H(+)-im) is rare. PT occurs via an Eigen-Zundel-Eigen scenario for switching the identity of imH(+) from an Eigen-like complex to another, intermediated by a Zundel-like complex. Structural and dynamical analyses demonstrate that PT in imidazole can be considered as a local event with very short spatial/temporal correlation, characterized by a few "rattling" or recurrent PT events. At long time scale, the trend of the PT correlation function may be recast with the diffusion model of reversible geminate recombination toward the power-law decay. The formation of the hydrogen bonds (HBs) for the imidazole molecules between the first and second solvation shell of imH(+) is crucial to pave the PT pathway. The above features may be understood by the flexibility of the HBs in liquid imidazole, as a stable HB network is essential for the Grotthuss mechanism.

Proton Conduction in Exchange Membranes across Multiple Length Scales

Concerns over global climate change associated with fossil-fuel consumption continue to drive the development of electrochemical alternatives for energy technology. Proton exchange fuel cells are a particularly promising technology for stationary power generation, mobile electronics, and hybrid engines in automobiles. For these devices to work efficiently, direct electrical contacts between the anode and cathode must be avoided; hence, the separator material must be electronically insulating but highly proton conductive. As a result, researchers have examined a variety of polymer electrolyte materials for use as membranes in these systems. In the optimization of the membrane, researchers are seeking high proton conductivity, low electronic conduction, and mechanical stability with the inclusion of water in the polymer matrix. A considerable number of potential polymer backbone and side chain combinations have been synthesized to meet these requirements, and computational studies can assist in the challenge of designing the next generation of technologically relevant membranes. Such studies can also be integrated in a feedback loop with experiment to improve fuel cell performance. However, to accurately simulate the currently favored class of membranes, perfluorosulfonic acid containing moieties, several difficulties must be addressed including a proper treatment of the proton-hopping mechanism through the membrane and the formation of nanophase-separated water networks. We discuss our recent efforts to address these difficulties using methods that push the limits of computer simulation and expand on previous theoretical developments. We describe recent advances in the multistate empirical valence bond (MS-EVB) method that can probe proton diffusion at the nanometer-length scale and accurately model the so-called Grotthuss shuttling mechanism for proton diffusion in water. Using both classical molecular dynamics and coarse-grained descriptions that replace atomistic representations with collective coordinates, we investigated the proton conductivity of polymer membrane structure as a function of hydration level. Nanometer-sized water channels form torturous pathways that are traversed by the charges during fuel cell operation. Using a combination of coarse-grained membrane structure and novel multiscale methods, we demonstrate emerging approaches to treat proton motion at the mesoscale in these complex materials.

Proton Transport Pathways in [NiFe]-Hydrogenase

Hydrogenases reversibly catalyze the production of molecular hydrogen. Current interest in these enzymes is focused on understanding the catalysis, since this may prove useful for hydrogen-based fuel cell and photosynthetic hydrogen production cell technologies. A key step in the hydrogenase catalytic cycle and the focus of this work is proton transport (PT) to and from the active site. The PT mechanism of the enzyme is studied using reactive molecular dynamics simulations of the full protein and the excess proton transfers via the multistate empirical valence bond (MS-EVB) method. Pathways connecting the bulk and the active site are located that suggest possible participation by several protonatable residues. PT free energy surfaces are calculated to differentiate the pathways.

Hydronium Behavior at the Air–Water Interface with a Polarizable Multistate Empirical Valence Bond Model

Molecular dynamics simulations were carried out to understand the propensity of the hydronium ion for the air–water interface with a polarizable multistate empirical valence bond (MS-EVB) model. Reasonable agreement with experiment for radial distribution functions and very good agreement for hydronium diffusion were found for the model. The polarizable MS-EVB model had no free energy minimum at the air–water interface. However, when polarizability on the hydronium ion alone was removed, a free energy of around −1.5 kcal/mol was calculated at the air–water interface. This discrepancy was found to be due to the behavior of water molecules in the first solvation shell of a hydronium ion. These water molecules contained a moderate amount of hydronium character, resulting in the delocalization of the hydronium ion. For the system with polarizable hydronium ions, this delocalization was the same at the interface as in the bulk, but for the system without polarizable hydronium ions, the delocalization increased...

A Refined MS-EVB Model for Proton Transport in Aqueous Environments

In order to improve the description of proton mobility in aqueous environments, a revised multistate empirical valence bond model (aMS-EVB3) is developed. The new aMS-EVB3 model is built upon an anharmonic water force field (aSPC/Fw) in which the OH bond potential is described through a quartic approximation to a Morse potential. First, it is shown that the aSPC/Fw anharmonic water model provides an accurate description of water at ambient conditions and reproduces the available experimental data for several structural, thermodynamic, and dynamical properties. Second, it is shown that, when applied to the study of proton solvation and transport in bulk water, the new aMS-EVB3 model accurately describes the solvation structure around the excess proton. Importantly, the new aMS-EVB3 model predicts a significantly larger proton diffusion coefficient than previous models, which largely improves the agreement with the available experimental data.

Expanding the view of proton pumping in cytochrome c oxidase through computer simulation

In cytochrome c oxidase (CcO), a redox-driven proton pump, protons are transported by the Grotthuss shuttling via hydrogen-bonded water molecules and protonatable residues. Proton transport through the D-pathway is a complicated process that is highly sensitive to alterations in the amino acids or the solvation structure in the channel, both of which can inhibit proton pumping and enzymatic activity. Simulations of proton transport in the hydrophobic cavity showed a clear redox state dependence. To study the mechanism of proton pumping in CcO, multi-state empirical valence bond (MS-EVB) simulations have been conducted, focusing on the proton transport through the D-pathway and the hydrophobic cavity next to the binuclear center. The hydration structures, transport pathways, effects of residues, and free energy surfaces of proton transport were revealed in these MS-EVB simulations. The mechanistic insight gained from them is herein reviewed and placed in context for future studies. This article is part of a Special Issue entitled: Respiratory Oxidases.

The Curious Case of the Hydrated Proton

Understanding the hydrated proton is a critically important problem that continues to engage the research efforts of chemists, physicists, and biologists because of its involvement in a wide array of phenomena. Only recently have several unique properties of the hydrated proton been unraveled through computer simulations. One such process is the detailed molecular mechanism by which protons hop between neighboring water molecules, thus giving rise to the anomalously high diffusion of protons relative to other simple cations. Termed Grotthuss shuttling, this process occurs over multiple time and length scales, presenting unique challenges for computer modeling and simulation. Because the hydrated proton is in reality a dynamical electronic charge defect that spans multiple water molecules, the simulation methodology must be able to dynamically readjust the chemical bonding topology. This reactive nature of the chemical process is automatically captured with ab initio molecular dynamics (AIMD) simulation methods, where the electronic degrees of freedom are treated explicitly. Unfortunately, these calculations can be prohibitively expensive for more complex proton solvation and transport phenomena in the condensed phase. These AIMD simulations remain extremely valuable, however, in validating empirical models, verifying results, and providing insight into molecular mechanisms. In this Account, we discuss recent progress in understanding the solvation and transport properties of the hydrated excess proton. The advances are based on results obtained from reactive molecular dynamics simulations using the multistate empirical valence bond (MS-EVB) methodology. This approach relies on a dynamic linear combination of chemical bond topologies to model charge delocalization and dynamic bonding environments. When parametrized via a variational force-matching algorithm from AIMD trajectories, the MS-EVB method can be viewed as a multiscale bridging of ab initio simulation results to a simpler and more efficient representation. The process allows sampling of longer time and length scales, which would normally be too computationally expensive with AIMD alone. With the MS-EVB methodology, the statistically important components of the excess proton solvation and hopping mechanisms in liquid water have been identified. The most likely solvation structure for the hydrated proton is a distorted Eigen-type complex (H(9)O(4)(+)). In this state, the excess proton charge defect rapidly resonates between three possible distorted Eigen-type structures until a successful proton hop occurs. This process, termed the "special-pair dance", serves as a kind of preparatory phase for the proton hopping while the neighboring water hydrogen-bonding network fluctuates and ultimately rearranges to facilitate a proton hop. The modifications of the solvation structure and transport properties of the excess proton in concentrated acid solutions were further investigated. The Eigen-type solvation structure also possesses both "hydrophilic" and "hydrophobic" sides, which accounts for the affinity of the hydrated proton for the air-water interface. This unusual "amphiphilic" character of the hydrated proton further leads to the metastable formation of contact ion pairs between two hydrated protons. It also engenders a surprisingly constant degree of solubility of hydrophobic species as a function of acid concentration, which contrasts with a markedly variable solubility as a function of salt (such as NaCl or KCl) concentration.

Proton Solvation and Transport in Hydrated Nafion

Proton solvation properties and transport mechanisms have been studied in hydrated Nafion using the self-consistent multistate empirical valence bond (SCI-MS-EVB) method that includes the effects excess proton charge defect delocalization and Grotthuss proton hopping. It was found that sulfonate groups influence excess proton solvation, as well as the proton hydration structure, by stabilizing a more Zundel-like (H(5)O(2)(+)) structure in their first solvation shells. Hydrate proton-related hydrogen bond networks were observed to be more stable than networks with water alone. Diffusion rates, Arrhenius activation energies, and transport pathways were calculated and analyzed to characterize the nature of the proton transport. Diffusion rate analysis suggests that a proton-hopping mechanism dominates the proton transport for the studied water loading levels and that there is a clear degree of anticorrelation with the vehicular transport. The activation energy drops quickly with an increasing water content when the water loading level is smaller than ∼10 H(2)O/SO(3)(-), which is consistent with experimental observations. The sulfonate groups were also found to affect the proton hopping directions. The temperature and water content effects on the proton transport pathways were also investigated.

Chemical rescue of enzymes: proton transfer in mutants of human carbonic anhydrase II.

In human carbonic anhydrase II (HCA II), the mutation of position 64 from histidine to alanine (H64A) disrupts the rate limiting proton transfer (PT) event, resulting in a reduction of the catalytic activity of the enzyme as compared to the wild-type. Potential of mean force (PMF) calculations utilizing the multistate empirical valence bond (MS-EVB) methodology for H64A HCA II yields a PT free energy barrier significantly higher than that found in the wild-type enzyme. This high barrier, determined in the absence of exogenous buffer and assuming no additional ionizable residues in the PT pathway, indicates the likelihood of alternate enzyme pathways that utilize either ionizable enzyme residues (self-rescue) and/or exogenous buffers (chemical rescue). It has been shown experimentally that the catalytic activity of H64A HCA II can be chemically rescued to near wild-type levels by the addition of the exogenous buffer 4-methylimidazole (4MI). Crystallographic studies have identified two 4MI binding sites, yet site-specific mutations intended to disrupt 4MI binding have demonstrated these sites to be nonproductive. In the present work, MS-EVB simulations show that binding of 4MI near Thr199 in the H64A HCA II mutant, a binding site determined by NMR spectroscopy, results in a viable chemical rescue pathway. Additional viable rescue pathways are also identified where 4MI acts as a proton transport intermediary from the active site to ionizable residues on the rim of the active site, revealing a probable mode of action for the chemical rescue pathway.

Infrared spectroscopy of small protonated water clusters at room temperature: An effective modes analysis

We perform infrared vibrational analysis on small protonated water clusters H(+)(H(2)O)(n), with n = 2, [ellipsis (horizontal)], 6, at room temperature. The absorption spectra are calculated based on classical trajectories obtained by the multistate empirical valence bond method. The analysis is carried out based on the effective modes analysis, which has been recently developed [Martinez et al., J. Chem. Phys. 125, 144106 (2006)] as generalization of the normal modes analysis. This technique enables us to decompose the full spectrum in maximally localized bands which are obtained by accounting for temperature and anharmonic effects. These effects are especially considered in the determination of the modes coupling. The spectra of the small clusters are interpreted by identifying the behavior of the excess charge, by understanding the role of hydrogen bonds, and by considering the effect of (micro-)solvation. Our results are presented by showing comparisons with other numerical methods and experimental measurements which are available in the literature.