Empirical Valence Bond Research Articles

A classical reactive potential for molecular clusters of sulphuric acid and water

ABSTRACTWe present a two-state empirical valence bond (EVB) potential describing interactions between sulphuric acid and water molecules and designed to model proton transfer between them within a classical dynamical framework. The potential has been developed in order to study the properties of molecular clusters of these species, which are thought to be relevant to atmospheric aerosol nucleation. The particle swarm optimisation method has been used to fit the parameters of the EVB model to density functional theory (DFT) calculations. Features of the parametrised model and DFT data are compared and found to be in satisfactory agreement. In particular, it is found that a single sulphuric acid molecule will donate a proton when clustered with four water molecules at 300 K and that this threshold is temperature dependent.

Molecular Dynamics Study of Proton Transport in Modeled Water Cluster Structure of Polymer Electrolyte Membrane

Molecular dynamics simulations have been performed to study the effects of water cluster structure on proton transport (PT) properties. The systems were built to closely approximate the proposed hydrophilic cluster structures (e.g. the cylinder model and the lamellar model) that are the most typical morphological models of Nafion. The PT properties were estimated in terms of diffusion coefficient and proton pathway using a modified two-state empirical valence bond (EVB) model. The diffusion coefficients in each dimension were calculated and correlated with the cluster size as well as the type of cluster models. It was found that proton diffusion is strongly affected by the cluster model and its size, which can be associated with the proton pathway in each cluster model. The influence of the sulfonate groups on the PT properties were also investigated. Our simulation results provide insight into quantitative information about PT properties in atomic level.

Effects of Water Structure on Proton Transport in Nafion Thin Films with Molecular Dynamics Simulations

Proton transport properties in Nafion thin films have been investigated using molecular dynamics simulations to characterize the nanoscopic flow phenomena in the electrode of polymer electrolyte fuel cell. Protons transfer in the ionomer by a combination of the Vehicle mechanism and the Grotthuss mechanism. In this study, the Grotthuss mechanism is described using the anharmonic two-state empirical valence bond method to evaluate quantitative proton transport within the framework of classical molecular dynamics simulations. The electrode consists of platinum-supported carbon and ionomers, and it is well known that the wettability of supported carbon depends on materials and operating environment. Our results suggested that the wettability affects morphology of the ionomer and proton transport.

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.

Empirical Valence Bond Theory Studies of the CH4 + Cl → CH3 + HCl Reaction

We report a theoretical investigation of the CH4 + Cl hydrogen abstraction reaction in the framework of empirical valence bond (EVB) theory. The main purpose of this study is to benchmark the EVB method against previous experimental and theoretical work. Analytical potential energy surfaces for the reaction have been developed on which quasi-classical trajectory calculations were carried out. The surfaces agree well with ab initio calculations at stationary points along the reaction path and dynamically relevant regions outside the reaction path. The analysis of dynamical data obtained using the EVB method, such as vibrational, rotational, and angular distribution functions, shows that this method compares well to both experimental measurements and higher-level theoretical calculations, with the additional benefit of low computational cost.

Expanding the Catalytic Triad in Epoxide Hydrolases and Related Enzymes.

Potato epoxide hydrolase 1 exhibits rich enantio- and regioselectivity in the hydrolysis of a broad range of substrates. The enzyme can be engineered to increase the yield of optically pure products as a result of changes in both enantio- and regioselectivity. It is thus highly attractive in biocatalysis, particularly for the generation of enantiopure fine chemicals and pharmaceuticals. The present work aims to establish the principles underlying the activity and selectivity of the enzyme through a combined computational, structural, and kinetic study using the substrate trans-stilbene oxide as a model system. Extensive empirical valence bond simulations have been performed on the wild-type enzyme together with several experimentally characterized mutants. We are able to computationally reproduce the differences between the activities of different stereoisomers of the substrate and the effects of mutations of several active-site residues. In addition, our results indicate the involvement of a previously neglected residue, H104, which is electrostatically linked to the general base H300. We find that this residue, which is highly conserved in epoxide hydrolases and related hydrolytic enzymes, needs to be in its protonated form in order to provide charge balance in an otherwise negatively charged active site. Our data show that unless the active-site charge balance is correctly treated in simulations, it is not possible to generate a physically meaningful model for the enzyme that can accurately reproduce activity and selectivity trends. We also expand our understanding of other catalytic residues, demonstrating in particular the role of a noncanonical residue, E35, as a “backup base” in the absence of H300. Our results provide a detailed view of the main factors driving catalysis and regioselectivity in this enzyme and identify targets for subsequent enzyme design efforts.

Computational Study of Competition between Direct Abstraction and Addition-Elimination in the Reaction of Cl Atoms with Propene.

Quasi-classical trajectory calculations on a newly constructed and full-dimensionality potential energy surface (PES) examine the dynamics of the reaction of Cl atoms with propene. The PES is an empirical valence bond (EVB) fit to high-level ab initio energies and incorporates deep potential energy wells for the 1-chloropropyl and 2-chloropropyl radicals, a direct H atom abstraction route to HCl + allyl radical (CH2CHCH2(•)) products (Δ(r)H(298K)(⊖) = −63.1 kJ mol(-1)), and a pathway connecting these regions. In total, 94 000 successful reactive trajectories were used to compute distributions of angular scattering and HCl vibrational and rotational level populations. These measures of the reaction dynamics agree satisfactorily with available experimental data. The dominant reaction pathway is direct abstraction of a hydrogen atom from the methyl group of propene occurring in under 500 fs. Less than 10% of trajectories follow an addition–elimination route via the two isomeric chloropropyl radicals. Large amplitude motions of the Cl about the propene molecular framework couple the addition intermediates to the direct abstraction pathway. The EVB method provides a good description of the complicated PES for the Cl + propene reaction despite fitting to a limited number of ab initio points, with the further advantage that dynamics specific to certain mechanisms can be studied in isolation by switching off coupling terms in the EVB matrix connecting different regions of the PES.

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.

Molecular Dynamics Study of Proton Transport in Modeled Water Cluster Structure of Polymer Electrolyte Membrane

The polymer electrolyte membrane (PEM) is one of the principal components for polymer electrolyte fuel cells. In the nanoscopic structure of the membrane, ion clusters are formed by water molecules gathered vicinity of sulfonate groups which are hydrophilic parts of Nafion membrane, and proton transport (PT) are largely attributed to the nanoscopic structure of polymer membranes and water aggregations. Thus, it is critical to understand the PT mechanisms through polymer membrane and clarify an important link between the membrane nanostructure and the PT properties. In the bulk aqueous solution, it is observed that proton mobility is roughly 5 times higher than other cations such as K+ which have an ionic radius similar to the hydronium ions. This property of PT is attributed to a combination of proton hopping, known as the Grotthuss mechanism, and vehicular mechanism where protons diffuse in solution as a hydrated form (e.g. hydronium ions). Many research efforts have been focused on understanding the PT mechanism using theoretical basis for the proton hopping models such as the empirical valence bond (EVB) approaches which allow PT phenomena to be simulated within a molecular dynamics (MD) framework. In PEM, it is considered that proton hopping transport still takes place although the water dynamic property in the bulk aqueous solution is inhibited by a lack of bulk-like water structure due to the phase separation in PEM. Since the phase separation is considered one of the primary factors affecting its performance, the morphological properties of PEM have been studied experimentally and been characterized by proposing the cluster models (e.g., the cylinder model, the lamellar model, and the cluster-rod model) which reasonably fit the experimental scattering spectra. However, a detailed relationship between the morphological features and PT is still under debate. In this study, atomistic MD simulations have been performed to study the effects of water cluster structure on PT properties. The systems were built to closely approximate the proposed hydrophilic cluster structures (the cylinder model and the lamellar model) that are the most major morphological models of Nafion. The PT properties were estimated in terms of diffusion coefficient and proton pathway using a modified EVB model developed based on the previous study of the two-state EVB model by Walbran et al. The diffusion coefficients in each dimension were calculated and correlated with the cluster size as well as the type of cluster models. It was found that proton diffusion is strongly affected by the cluster model and its size, which can be associated with the proton pathway in each cluster model. The influence of the sulfonate groups on the PT properties were also investigated. Our simulation results provide insight into quantitative information about PT properties in atomic level.

Effects of Water Structure on Proton Transport in Nafion Thin Films with Molecular Dynamics Simulations

Polymer electrolyte fuel cell (PEFC) is expected to play important roles for next generation power systems in terms of high efficiency and low environmental burden and capability for emergency use, but there are still many hurdles to be largely commercialized. One of them results from problems in catalyst layers (CL). What required for CL is, for example, reducing the amount of Pt particles with keeping high current density. To accomplish this, knowledge about optimum structures of CL based on material transport properties is needed because there is a strong correlation between current density and transport properties. Although there are many reports about material transport properties in the polymer electrolyte membrane (PEM), it is predicted that material transport properties in CL show some unique characteristics because of effects of interfaces in CL. It means that the knowledge obtained from PEM should not be applied to CL. In addition, the researches about CL have been tried by means of experiments and simulations, but especially the transport properties have not been clarified yet. There are nanoscopic structures in Nafion thin films and it is hard to capture the transport phenomena in such a small area by experiments and macroscopic simulations. In this study, we analyze proton transport properties in the cathode CL using molecular dynamics simulations, focusing on the dependences of the diffusion coefficient on water content. We consider the proton transport mechanism as a combination of the Vehicular mechanism and the Grotthuss mechanism. The Grotthuss mechanism is represented using empirical valence bond methods. It has been shown by previous researches in PEM that the Grotthuss mechanism is a dominant factor for proton diffusion in relatively high water contents. Thus, it is essential to consider the Grotthuss mechanism for analysis of proton transport in CL as well. The system consists of graphite surface and Nafion thin films. Pt particles are neglected because of a restriction of calculation cost. In PEM analysis, the diffusion coefficient of protons increases with increasing water contents. First, we report whether the same trend as observed in PEM can be seen in CL and the correlation between diffusibility and structures. Next, unique characteristics of proton transport properties in CL compared to the PEM analysis are reported. They are evaluated in terms of mass distribution of the thin films and solvent molecules. From this analysis, it is found that the diffusibility is not much affected by surface because of a lack of solvent molecules on surface, which is attributed to the hydrophobic properties of graphite.

A modified two-state empirical valence bond model for proton transport in aqueous solutions.

A detailed analysis of the proton solvation structure and transport properties in aqueous solutions is performed using classical molecular dynamics simulations. A refined two-state empirical valence bond (aTS-EVB) method, which is based on the EVB model of Walbran and Kornyshev and the anharmonic water force field, is developed in order to describe efficiently excess proton transport via the Grotthuss mechanism. The new aTS-EVB model clearly satisfies the requirement for simpler and faster calculation, because of the simplicity of the two-state EVB algorithm, while providing a better description of diffusive dynamics of the excess proton and water in comparison with the previous two-state EVB models, which significantly improves agreement with the available experimental data. The results of activation energies for the excess proton and water calculated between 300 and 340 K (the temperature range used in this study) are also found to be in good agreement with the corresponding experimental data.

Chemical reaction mechanisms in solution from brute force computational Arrhenius plots.

Decomposition of activation free energies of chemical reactions, into enthalpic and entropic components, can provide invaluable signatures of mechanistic pathways both in solution and in enzymes. Owing to the large number of degrees of freedom involved in such condensed-phase reactions, the extensive configurational sampling needed for reliable entropy estimates is still beyond the scope of quantum chemical calculations. Here we show, for the hydrolytic deamination of cytidine and dihydrocytidine in water, how direct computer simulations of the temperature dependence of free energy profiles can be used to extract very accurate thermodynamic activation parameters. The simulations are based on empirical valence bond models, and we demonstrate that the energetics obtained is insensitive to whether these are calibrated by quantum mechanical calculations or experimental data. The thermodynamic activation parameters are in remarkable agreement with experiment results and allow discrimination among alternative mechanisms, as well as rationalization of their different activation enthalpies and entropies.

Qgui: A high-throughput interface for automated setup and analysis of free energy calculations and empirical valence bond simulations in biological systems

Structural information and activity data has increased rapidly for many protein targets during the last decades. In this paper, we present a high-throughput interface (Qgui) for automated free energy and empirical valence bond (EVB) calculations that use molecular dynamics (MD) simulations for conformational sampling. Applications to ligand binding using both the linear interaction energy (LIE) method and the free energy perturbation (FEP) technique are given using the estrogen receptor (ERα) as a model system. Examples of free energy profiles obtained using the EVB method for the rate-limiting step of the enzymatic reaction catalyzed by trypsin are also shown. In addition, we present calculation of high-precision Arrhenius plots to obtain the thermodynamic activation enthalpy and entropy with Qgui from running a large number of EVB simulations.

Molecular Dynamics Simulation of Proton Transport in Polymer Electrolyte Membrane.

We have performed a detailed analysis of proton solvation and transport properties in hydrated Nafion using molecular dynamics simulation. The revised empirical valence bond (EVB) method was developed in order to treat the excess proton transport through the Grotthuss mechanism. The new EVB model predicts a significantly enhanced transport in comparison with previous hopping models as well as the classical hydronium diffusion, which largely improves the agreement with the available experimental data. Our results suggest that a proton hopping mechanism has a small effect on the proton dissociation from the first solvation shell of sulfonate groups, namely that protons are not enhanced to separate from the sulfonate groups by the hopping mechanisms. From diffusion comparison between the Grotthuss and vehicular mechanism, the Grotthuss mechanism dominates the proton diffusion at the studied hydration levels including a hydration level of 3. It was also found that the vehicular mechanism dominates the electroosmotic transport of water molecules at the studied hydration levels.

On the challenge of exploring the evolutionary trajectory from phosphotriesterase to arylesterase using computer simulations.

The ability to design effective enzymes presents a fundamental challenge in biotechnology and also in biochemistry. Unfortunately, most of the progress on this field has been accomplished by bringing the reactants to a reasonable orientation relative to each other, rather than by rational optimization of the polar preorganization of the environment, which is the most important catalytic factor. True computer based enzyme design would require the ability to evaluate the catalytic power of designed active sites. This work considers the evolution from a phosphotriesterase (with the paraoxon substrate) to arylesterase (with the 2-naphthylhexanoate (2NH) substrate) catalysis. Both the original and the evolved enzymes involve two zinc ions and their ligands, making it hard to obtain a reliable quantum mechanical description and then to obtain an effective free energy sampling. Furthermore, the options for the reaction path are quite complicated. To progress in this direction we started with DFT calculations of the energetics of different mechanistic options of cluster models and then used the results to calibrate empirical valence bond (EVB) models and to generate properly sampled free energy surfaces for different mechanisms in the enzyme. Interestingly, it is found that the catalytic effect depends on the Zn-Zn distance making the mechanistic analysis somewhat complicated. Comparing the activation barriers of paraoxon and the 2NH ester at the beginning and end of the evolutionary path reproduced the observed evolutionary trend. However, although our findings provide an advance in exploring the nature of promiscuous enzymes, they also indicate that modeling the reaction mechanism in the case of enzymes with a binuclear zinc center is far from trivial and presents a challenge for computer-aided enzyme design.

Catalytic Stimulation by Restrained Active-Site Floppiness—The Case of High Density Lipoprotein-Bound Serum Paraoxonase-1

Despite the abundance of membrane-associated enzymes, the mechanism by which membrane binding stabilizes these enzymes and stimulates their catalysis remains largely unknown. Serum paraoxonase-1 (PON1) is a lipophilic lactonase whose stability and enzymatic activity are dramatically stimulated when associated with high-density lipoprotein (HDL) particles. Our mutational and structural analyses, combined with empirical valence bond simulations, reveal a network of hydrogen bonds that connect HDL binding residues with Asn168—a key catalytic residue residing >15Å from the HDL contacting interface. This network ensures precise alignment of N168, which, in turn, ligates PON1's catalytic calcium and aligns the lactone substrate for catalysis. HDL binding restrains the overall motion of the active site and particularly of N168, thus reducing the catalytic activation energy barrier. We demonstrate herein that disturbance of this network, even at its most far-reaching periphery, undermines PON1's activity. Membrane binding thus immobilizes long-range interactions via second- and third-shell residues that reduce the active site's floppiness and pre-organize the catalytic residues. Although this network is critical for efficient catalysis, as demonstrated here, unraveling these long-rage interaction networks is challenging, let alone their implementation in artificial enzyme design.

Methyltransferases do not work by compression, cratic, or desolvation effects, but by electrostatic preorganization.

The enzyme catechol O-methyltransferase (COMT) catalyzes the transfer of a methyl group from S-adenosylmethionine to dopamine and related catechols. The search for the origin of COMT catalysis has led to different proposals and hypothesis, including the entropic, the NAC, and the compression proposals as well as the more reasonable electrostatic idea. Thus, it is important to understand the catalytic power of this enzyme and to examine the validity of different proposals and in particular the repeated recent implication of the compression idea. The corresponding analysis should be done by well-defined physically-based considerations that involve computations rather than circular interpretations of experimental results. Thus, we explore here the origin of the catalytic efficiency of COMT by using the empirical valence bond and the linear response approximation approaches. The results demonstrate that the catalytic effect of COMT is mainly due to electrostatic preorganization effects. It is also shown that the compression, NAC and entropic proposals do not account for the catalytic effect.

Prediction of enhanced solvent-induced enantioselectivity for a ring opening with a bifurcating reaction path.

Classical molecular dynamics simulations are reported for the deazetisation and ring opening of meso-2,3-difluoro-2,3-dimethyldiazocyclopropane in three solvents: CHCl3, CHFClBr and CH3CH(OH)CF3 (TFIPA). The achiral reactant leads to enantiomeric allene products, and the question addressed in the study is whether either of the chiral, enantiomerically pure solvents can induce significant enantiomeric excess in the products. The direct dynamics calculations use an empirical valence bond potential for the solute, with empirical parameters optimised against M06-2X/cc-pVTZ density functional results. The results reveal that the exothermic N2 loss and ring opening promote transient strong solvent-solute interactions within the first ∼100 fs of the reaction. Because of the bifurcating reaction path, these interactions occur at time when the "decision" about which enantiomer of the product to form has yet to be made (at least for many of the trajectories). Hence, it is possible in principle that the solvent could exert a larger-than-normal influence on the course of the reaction. In fact, the results reveal no such effect for CHFClBr but do predict that TFIPA should induce 15.2 ± 2.1% enantiomeric excess. This is roughly an order of magnitude larger than solvent-induced enantiomeric excesses found experimentally in reactions where the conversion of reactant(s) to enantiomeric products occur over separate transition states.

J0560104 分子動力学法を用いた高分子電解質膜内におけるプロトン輸送現象の解析

The effects of water cluster structures on the PT properties have been investigated using reactive molecular dynamics (MD) simulations. The two-state empirical valence bond (aTS-EVB) model was used to allow proton transport (PT) phenomena including the Grotthuss mechanism to be simulated within the simplicity of the theoretical framework in MD simulations. The systems were constructed to reproduce the experimentally proposed hydrophilic cluster structures (the cylinder model, the lamellar model, and the sphere-rod model) in the Nafion systems. The diffusion coefficients in each dimension were estimated to characterize PT properties. It was found that proton diffusion strongly correlates with the cluster size as well as the type of cluster models. The proton diffusion in the lamellar model is higher than that in the other two types of models when the water volumes are the same. In addition, the rod radius in the sphere-rod model, i.e., the bottleneck of the cluster domain, strongly affects the proton diffusion, suggesting that the rod radius, regardless of the sphere radius, is a dominant factor for the proton diffusion in the sphere-rod model.