A Hamiltonian electrostatic coupling scheme for hybrid Car–Parrinello molecular dynamics simulations

Exploring reaction pathways for O-GlcNAc transferase catalysis. A string method study

The inverting O-GlcNAc glycosyltransferase (OGT) is an important post-translation enzyme, which catalyzes the transfer of N-acetylglucosamine from UDP-N-acetylglucosamine (UDP-GlcNAc) to the hydroxyl group of the Ser/Thr of cytoplasmic, nuclear, and mitochondrial proteins. In the past, three different catalytic bases were proposed for the reaction: His498, α-phosphate, and Asp554. In this study, we used hybrid quantum mechanics/molecular mechanics (QM/MM) Car-Parrinello molecular dynamics to investigate reaction paths using α-phosphate and Asp554 as the catalytic bases. The string method was used to calculate the free-energy reaction profiles of the tested mechanisms. During the investigations, an additional mechanism was observed. In this mechanism, a proton is transferred to α-phosphate via a water molecule. Our calculations show that the mechanism with α-phosphate acting as the base is favorable. This reaction has a rate-limiting free-energy barrier of 23.5 kcal/mol, whereas reactions utilizing Asp554 and water-assisted α-phosphate have barriers of 41.7 and 40.9 kcal/mol, respectively. Our simulations provide a new insight into the catalysis of OGT and may thus guide rational drug design of transition-state analogue inhibitors with potential therapeutic use.

Nonadiabatic molecular dynamics simulations: synergies between theory and experiments

Recent developments in nonadiabatic dynamics enabled ab inito simulations of complex ultrafast processes in the condensed phase. These advances have opened new avenues in the study of many photophysical and photochemical reactions triggered by the absorption of electromagnetic radiation. In particular, theoretical investigations can be combined with the most sophisticated femtosecond experimental techniques to guide the interpretation of measured time-resolved observables. At the same time, the availability of experimental data at high (spatial and time) resolution offers a unique opportunity for the benchmarking and the improvement of those theoretical models used to describe complex molecular systems in their natural environment. The established synergy between theory and experiments can produce a better understanding of new ultrafast physical and chemical processes at atomistic scale resolution. Furthermore, reliable ab inito molecular dynamics simulations can already be successfully employed as predictive tools to guide new experiments as well as the design of novel and better performing materials. In this paper, I will give a concise account on the state of the art of molecular dynamics simulations of complex molecular systems in their excited states. The principal aim of this approach is the description of a given system of interest under the most realistic ambient conditions including all environmental effects that influence experiments, for instance, the interaction with the solvent and with external time-dependent electric fields, temperature, and pressure. To this end, time-dependent density functional theory (TDDFT) is among the most efficient and accurate methods for the representation of the electronic dynamics, while trajectory surface hopping gives a valuable representation of the nuclear quantum dynamics in the excited states (including nonadiabatic effects). Concerning the environment and its effects on the dynamics, the quantum mechanics/molecular mechanics (QM/MM) approach has the advantage of providing an atomistic (even though approximated) description of the solvent molecules, which is crucial for the characterization of all ultrafast relaxation phenomena that depend on the geometrical arrangement at the interface between a molecule and the solvent, for example, the hydrogen bond network. After a short description of the TDDFT-based implementation of Ehrenfest and trajectory surface hopping dynamics, I will present applications in different domains of molecular chemistry and physics: the analysis and the understanding of (time-resolved) X-ray absorption spectra, the interpretation of the ultrafast relaxation dynamics of photoexcited dyes in solution, and the design of specific laser pulses (capable of inducing desired chemical reactions) using local control theory.

Utilizing fast multipole expansions for efficient and accurate quantum-classical molecular dynamics simulations

Recently, a novel approach to hybrid quantum mechanics/molecular mechanics (QM/MM) molecular dynamics (MD) simulations has been suggested [Schwörer et al., J. Chem. Phys. 138, 244103 (2013)]. Here, the forces acting on the atoms are calculated by grid-based density functional theory (DFT) for a solute molecule and by a polarizable molecular mechanics (PMM) force field for a large solvent environment composed of several 10(3)-10(5) molecules as negative gradients of a DFT/PMM hybrid Hamiltonian. The electrostatic interactions are efficiently described by a hierarchical fast multipole method (FMM). Adopting recent progress of this FMM technique [Lorenzen et al., J. Chem. Theory Comput. 10, 3244 (2014)], which particularly entails a strictly linear scaling of the computational effort with the system size, and adapting this revised FMM approach to the computation of the interactions between the DFT and PMM fragments of a simulation system, here, we show how one can further enhance the efficiency and accuracy of such DFT/PMM-MD simulations. The resulting gain of total performance, as measured for alanine dipeptide (DFT) embedded in water (PMM) by the product of the gains in efficiency and accuracy, amounts to about one order of magnitude. We also demonstrate that the jointly parallelized implementation of the DFT and PMM-MD parts of the computation enables the efficient use of high-performance computing systems. The associated software is available online.

The complete conformational free energy landscape of β-xylose reveals a two-fold catalytic itinerary for β-xylanases

Unraveling the conformational catalytic itinerary of glycoside hydrolases (GHs) is a growing topic of interest in glycobiology, with major impact in the design of GH inhibitors. β-xylanases are responsible for the hydrolysis of glycosidic bonds in β-xylans, a group of hemicelluloses of high biotechnological interest that are found in plant cell walls. The precise conformations followed by the substrate during catalysis in β-xylanases have not been unambiguously resolved, with three different pathways being proposed from structural analyses. In this work, we compute the conformational free energy landscape (FEL) of β-xylose to predict the most likely catalytic itineraries followed by β-xylanases. The calculations are performed by means of ab initio metadynamics, using the Cremer-Pople puckering coordinates as collective variables. The computed FEL supports only two of the previously proposed itineraries, 2SO → [2,5B]ǂ → 5S1 and 1S3 → [4H3]ǂ → 4C1, which clearly appear in low energy regions of the FEL. Consistently, 2SO and 1S3 are conformations preactivated for catalysis in terms of free energy/anomeric charge and bond distances. The results however exclude the OE → [OS2]ǂ → B2,5 itinerary that has been recently proposed for a family 11 xylanase. Classical and ab initio QM/MM molecular dynamics simulations reveal that, in this case, the observed OE conformation has been enforced by enzyme mutation. These results add a word of caution on using modified enzymes to inform on catalytic conformational itineraries of glycoside hydrolases.

Unraveling the degradation of artificial amide bonds in nylon oligomer hydrolase: from induced-fit to acylation processes

To unravel the factor that provides the ability to degrade non-biological amide bond with nylon oligomer hydrolase, we investigated the process from induced-fit to acylation by a combination of different theoretical methods.

Human infrared vision is triggered by two-photon chromophore isomerization

Vision relies on photoactivation of visual pigments in rod and cone photoreceptor cells of the retina. The human eye structure and the absorption spectra of pigments limit our visual perception of light. Our visual perception is most responsive to stimulating light in the 400- to 720-nm (visible) range. First, we demonstrate by psychophysical experiments that humans can perceive infrared laser emission as visible light. Moreover, we show that mammalian photoreceptors can be directly activated by near infrared light with a sensitivity that paradoxically increases at wavelengths above 900 nm, and display quadratic dependence on laser power, indicating a nonlinear optical process. Biochemical experiments with rhodopsin, cone visual pigments, and a chromophore model compound 11-cis-retinyl-propylamine Schiff base demonstrate the direct isomerization of visual chromophore by a two-photon chromophore isomerization. Indeed, quantum mechanics modeling indicates the feasibility of this mechanism. Together, these findings clearly show that human visual perception of near infrared light occurs by two-photon isomerization of visual pigments.

Asymmetric environmental effects on the structure and vibrations of cis-[Pt(NH3)2Cl2] in condensed phases

We report the structural and vibrational properties of anticancer drug cisplatin (cis-[Pt(NH3)2Cl2]) in gas phase, in solid phase, and in aqueous solution using density functional theory (DFT) calculations, quantum mechanical/molecular mechanical (QM/MM) molecular dynamics, and effective normal modes analysis. In contrast with the gas-phase case, asymmetric hydrogen bonding environments are found in both solid phase and aqueous solution. It is shown that the discrepancy of the molecular geometry between previous gas phase calculations and the X-ray crystal structure can be resolved by considering intermolecular hydrogen bonds in the calculations of solid phase. In addition, our simulations in solid phase and aqueous solution reveal that asymmetric environmental effects lead to several spectral features observed in experiments, such as the blue-shift in the N-H stretching region and the frequency splitting of NH3 symmetric deformation modes. Furthermore, a similar decoupling and localization of several vibrational modes of cisplatin is found in solid phase and aqueous solution, in comparison to those of O-H stretching modes of water molecules in liquid water [ J. Phys. Chem. Lett. 2013 , 4 ( 19 ), 3245 - 3250 ].

Molecular simulation of water and hydration effects in different environments: challenges and developments for DFTB based models

We discuss the description of water and hydration effects that employs an approximate density functional theory, DFTB3, in either a full QM or QM/MM framework. The goal is to explore, with the current formulation of DFTB3, the performance of this method for treating water in different chemical environments, the magnitude and nature of changes required to improve its performance, and factors that dictate its applicability to reactions in the condensed phase in a QM/MM framework. A relatively minor change (on the scale of kBT) in the O-H repulsive potential is observed to substantially improve the structural properties of bulk water under ambient conditions; modest improvements are also seen in dynamic properties of bulk water. This simple change also improves the description of protonated water clusters, a solvated proton, and to a more limited degree, a solvated hydroxide. By comparing results from DFTB3 models that differ in the description of water, we confirm that proton transfer energetics are adequately described by the standard DFTB3/3OB model for meaningful mechanistic analyses. For QM/MM applications, a robust parametrization of QM-MM interactions requires an explicit consideration of condensed phase properties, for which an efficient sampling technique was developed recently and is reviewed here. The discussions help make clear the value and limitations of DFTB3 based simulations, as well as the developments needed to further improve the accuracy and transferability of the methodology.

Computational insights into function and inhibition of fatty acid amide hydrolase

The Fatty Acid Amide Hydrolase (FAAH) enzyme is a membrane-bound serine hydrolase responsible for the deactivating hydrolysis of a family of naturally occurring fatty acid amides. FAAH is a critical enzyme of the endocannabinoid system, being mainly responsible for regulating the level of its main cannabinoid substrate anandamide. For this reason, pharmacological inhibition of FAAH, which increases the level of endogenous anandamide, is a promising strategy to cure a variety of diseases including pain, inflammation, and cancer. Much structural, mutagenesis, and kinetic data on FAAH has been generated over the last couple of decades. This has prompted several informative computational investigations to elucidate, at the atomic-level, mechanistic details on catalysis and inhibition of this pharmaceutically relevant enzyme. Here, we review how these computational studies - based on classical molecular dynamics, full quantum mechanics, and hybrid QM/MM methods - have clarified the binding and reactivity of some relevant substrates and inhibitors of FAAH. We also discuss the experimental implications of these computational insights, which have provided a thoughtful elucidation of the complex physical and chemical steps of the enzymatic mechanism of FAAH. Finally, we discuss how computations have been helpful for building structure-activity relationships of potent FAAH inhibitors.

The CHARMM-TURBOMOLE interface for efficient and accurate QM/MM molecular dynamics, free energies, and excited state properties

The quantum mechanical (QM)/molecular mechanical (MM) interface between Chemistry at HARvard Molecular Mechanics (CHARMM) and TURBOMOLE is described. CHARMM provides an extensive set of simulation algorithms, like molecular dynamics (MD) and free energy perturbation, and support for mature nonpolarizable and Drude polarizable force fields. TURBOMOLE provides fast QM calculations using density functional theory or wave function methods and excited state properties. CHARMM-TURBOMOLE is well-suited for extended QM/MM MD simulations using first principles methods with large (triple-ζ) basis sets. We demonstrate these capabilities with a QM/MM simulation of Mg(2+) (aq), where the MM outer sphere water molecules are represented using the SWM4-NDP Drude polarizable force field and the ion and inner coordination sphere are represented using QM PBE, PBE0, and MP2 methods. The relative solvation free energies of Mg(2+) and Zn(2+) were calculated using thermodynamic integration. We also demonstrate the features for excited state properties. We calculate the time-averaged solution absorption spectrum of indole, the emission spectrum of the indole 1La excited state, and the electronic circular dichroism spectrum of an oxacepham.

Molecular Recognition of Platinated DNA from Chromosomal HMGB1

Cisplatin cures testicular and ovarian cancers with unprecedented potency. It induces its beneficial activity by covalently binding to DNA. Repair enzymes, which remove the platinated lesions from DNA, cause drug resistance. Chromosomal High Mobility Group Box proteins (HMGB) may interfere with this process by binding to platinated DNA. Using 8 μs multiple-walker well-tempered metadynamics simulations, here, we investigated the structural and the energetic determinants of one of the HMGB proteins (HMGB1A) in complex with the platinated oligonucleotide [Pt(NH3)2](2+)-d(CCUCTCTG*G*ACCTTCC)-d(GGAGAGACCTGGAAGG) (*G are platinated guanines), for which experimental structural information is available. The calculated affinity is in good agreement with experiment. The process is predicted to be enthalpy-driven, as found for other protein/DNA complexes. The Lys7 residue, whose side-chain was not resolved in the X-ray structure, is found to interact with the C4 5'-phosphate and this interaction emerges as a key facet for the molecular recognition process. In addition, our calculations provide a molecular basis for the experimentally measured decreased affinity of HMGB1A for platinated DNA, as a consequence of Cys22-Cys44 S-S bridge formation (such an oxidation cannot take place in some members of this protein family present in the testis, where the drug is particularly effective). This decrease is likely to be caused by a small yet significant rearrangement of helices H1 and H2 with consequent alteration of the Phe37 juxtaposition.

QM/MM molecular dynamics studies of metal binding proteins

Mixed quantum-classical (quantum mechanical/molecular mechanical (QM/MM)) simulations have strongly contributed to providing insights into the understanding of several structural and mechanistic aspects of biological molecules. They played a particularly important role in metal binding proteins, where the electronic effects of transition metals have to be explicitly taken into account for the correct representation of the underlying biochemical process. In this review, after a brief description of the basic concepts of the QM/MM method, we provide an overview of its capabilities using selected examples taken from our work. Specifically, we will focus on heme peroxidases, metallo-β-lactamases, α-synuclein and ligase ribozymes to show how this approach is capable of describing the catalytic and/or structural role played by transition (Fe, Zn or Cu) and main group (Mg) metals. Applications will reveal how metal ions influence the formation and reduction of high redox intermediates in catalytic cycles and enhance drug metabolism, amyloidogenic aggregate formation and nucleic acid synthesis. In turn, it will become manifest that the protein frame directs and modulates the properties and reactivity of the metal ions.

Are dinucleoside monophosphates relevant models for the study of DNA intrastrand cross-link lesions? The example of g[8-5m]T

Oxidatively generated tandem lesions such as G[8-5m]T pose a potent threat to genome integrity. Direct experimental studies of the kinetics and thermodynamics of a specific lesion within DNA are very challenging, mostly due to the variety of products that can be formed in oxidative conditions. Dinucleoside monophosphates (DM) involving only the reactive nucleobases in water represent appealing alternative models on which most physical chemistry and structural techniques can be applied. However, it is not yet clear how relevant these models are. Here, we present QM/MM MD simulations of the cyclization step involved in the formation of G[8-5m]T from the guanine-thymine (GpT) DM in water, with the aim of comparing our results to our previous investigation of the same reaction in DNA ( Garrec , J. , Patel , C. , Rothlisberger , U. , and Dumont , E. ( 2012 ) J. Am. Chem. Soc. 134 , 2111 - 2119 ). We show that, despite the different levels of preorganization of the two systems, the corresponding reactions share many energetic and structural characteristics. The main difference lies in the angle between the G and T bases, which is slightly higher in the transition state (TS) and product of the reaction in water than in the reaction in DNA. This effect is due to the Watson-Crick H-bonds, which are absent in the {GpT+water} system and restrain the relative positioning of the reactive nucleobases in DNA. However, since the lesion is accommodated easily in the DNA macromolecule, the induced energetic penalty is relatively small. The high similarity between the two reactions strongly supports the use of GpT in water as a model of the corresponding reaction in DNA.

A mechanochemical switch to control radical intermediates

B₁₂-dependent enzymes employ radical species with exceptional prowess to catalyze some of the most chemically challenging, thermodynamically unfavorable reactions. However, dealing with highly reactive intermediates is an extremely demanding task, requiring sophisticated control strategies to prevent unwanted side reactions. Using hybrid quantum mechanical/molecular mechanical simulations, we follow the full catalytic cycle of an AdoB₁₂-dependent enzyme and present the details of a mechanism that utilizes a highly effective mechanochemical switch. When the switch is "off", the 5'-deoxyadenosyl radical moiety is stabilized by releasing the internal strain of an enzyme-imposed conformation. Turning the switch "on," the enzyme environment becomes the driving force to impose a distinct conformation of the 5'-deoxyadenosyl radical to avoid deleterious radical transfer. This mechanochemical switch illustrates the elaborate way in which enzymes attain selectivity of extremely chemically challenging reactions.

Reaction Path Force Matching: A New Strategy of Fitting Specific Reaction Parameters for Semiempirical Methods in Combined QM/MM Simulations

We present a general strategy of reparametrizing semiempirical (SE) methods against ab initio (AI) methods for combined quantum mechanical and molecular mechanical (QM/MM) simulations of specific chemical reactions in condensed phases. The resulting approach, designated Reaction Path Force Matching (RP-FM), features cycles of sampling configurations along a reaction path on an efficient SE/MM potential energy surface (PES) and adjusting specific reaction parameters (SRPs) in the SE method such that the atomic forces computed at the target AI/MM level are reproduced. Iterative applications of the RP-FM cycle make possible achieving the accuracy of AI/MM simulations without explicitly sampling the computationally expensive AI/MM PES. The bypassed sampling, nevertheless, is implicitly accomplished through the aid of the efficient SE-SRP/MM PES, on which the target-level reaction path is expected to be obtained upon convergence. We demonstrate the effectiveness of the RP-FM procedure for a symmetric proton transfer reaction in the gas phase and in solution. The remarkable agreements between the RP-FM optimized SE-SRP methods and the target AI method on various properties, including energy profiles, potential of mean force free energy profiles, atomic forces, charge populations, and solvation effects, suggest that RP-FM can be used as an efficient and reliable strategy for simulating condensed-phase chemical reactions.

Molecular mechanism of a hotdog-fold acyl-CoA thioesterase

Thioesterases are enzymes that hydrolyze thioester bonds between a carbonyl group and a sulfur atom. They catalyze key steps in fatty acid biosynthesis and metabolism, as well as polyketide biosynthesis. The reaction molecular mechanism of most hotdog-fold acyl-CoA thioesterases remains unknown, but several hypotheses have been put forward in structural and biochemical investigations. The reaction of a human thioesterase (hTHEM2), representing a thioesterase family with a hotdog fold where a coenzyme A moiety is cleaved, was simulated by quantum mechanics/molecular mechanics metadynamics techniques to elucidate atomic and electronic details of its mechanism, its transition-state conformation, and the free energy landscape of the process. A single-displacement acid-base-like mechanism, in which a nucleophilic water molecule is activated by an aspartate residue acting as a base, was found, confirming previous experimental proposals. The results provide unambiguous evidence of the formation of a tetrahedral-like transition state. They also explain the roles of other conserved active-site residues during the reaction, especially that of a nearby histidine/serine pair that protonates the thioester sulfur atom, the participation of which could not be elucidated from mutation analyses alone.

Study of the redox properties of singlet and triplet Tris(2,2'-bipyridine)ruthenium(II) ([Ru(bpy)3]2+) in aqueous solution by full quantum and mixed quantum/classical molecular dynamics simulations

The oxidation of ground-state (singlet) and triplet [Ru(bpy)3](2+) were studied by full quantum-mechanical (QM) and mixed quantum/classical (QM/MM) molecular dynamics simulations. Both approaches provide reliable results for the redox potentials of the two spin states. The two redox reactions closely obey Marcus theory for electron transfer. The free energy difference between the two [Ru(bpy)3](2+) states amounts to 1.78 eV from both QM and QM/MM simulations. The two methods also provide similar results for the reorganization free energy associated with the transition from singlet to triplet [Ru(bpy)3](2+) (0.06 eV for QM and 0.07 eV for QM/MM). On the basis of single-point calculations, we estimate the entropic contribution to the free energy difference between singlet and triplet [Ru(bpy)3](2+) to be 0.27 eV, which is significantly greater than previously assumed (0.03 eV) and in contradiction with the assumption that the transition between these two states can be accurately described using purely energetic considerations. Employing a thermodynamic cycle involving singlet [Ru(bpy)3](2+), triplet [Ru(bpy)3](2+), and [Ru(bpy)3](3+), we calculated the triplet oxidation potential to be -0.62 V vs the standard hydrogen electrode, which is significantly different from a previous experimental estimate based on energetic considerations only (-0.86 V).

Integrated Hamiltonian sampling: a simple and versatile method for free energy simulations and conformational sampling

Motivated by specific applications and the recent work of Gao and co-workers on integrated tempering sampling (ITS), we have developed a novel sampling approach referred to as integrated Hamiltonian sampling (IHS). IHS is straightforward to implement and complementary to existing methods for free energy simulation and enhanced configurational sampling. The method carries out sampling using an effective Hamiltonian constructed by integrating the Boltzmann distributions of a series of Hamiltonians. By judiciously selecting the weights of the different Hamiltonians, one achieves rapid transitions among the energy landscapes that underlie different Hamiltonians and therefore an efficient sampling of important regions of the conformational space. Along this line, IHS shares similar motivations as the enveloping distribution sampling (EDS) approach of van Gunsteren and co-workers, although the ways that distributions of different Hamiltonians are integrated are rather different in IHS and EDS. Specifically, we report efficient ways for determining the weights using a combination of histogram flattening and weighted histogram analysis approaches, which make it straightforward to include many end-state and intermediate Hamiltonians in IHS so as to enhance its flexibility. Using several relatively simple condensed phase examples, we illustrate the implementation and application of IHS as well as potential developments for the near future. The relation of IHS to several related sampling methods such as Hamiltonian replica exchange molecular dynamics and λ-dynamics is also briefly discussed.

Ligand substitutions between ruthenium-cymene compounds can control protein versus DNA targeting and anticancer activity

Ruthenium compounds have become promising alternatives to platinum drugs by displaying specific activities against different cancers and favourable toxicity and clearance properties. Nonetheless, their molecular targeting and mechanism of action are poorly understood. Here we study two prototypical ruthenium-arene agents-the cytotoxic antiprimary tumour compound [(η(6)-p-cymene)Ru(ethylene-diamine)Cl]PF6 and the relatively non-cytotoxic antimetastasis compound [(η(6)-p-cymene)Ru(1,3,5-triaza-7-phosphaadamantane)Cl2]-and discover that the former targets the DNA of chromatin, while the latter preferentially forms adducts on the histone proteins. Using a novel 'atom-to-cell' approach, we establish the basis for the surprisingly site-selective adduct formation behaviour and distinct cellular impact of these two chemically similar anticancer agents, which suggests that the cytotoxic effects arise largely from DNA lesions, whereas the protein adducts may be linked to the other therapeutic activities. Our study shows promise for developing new ruthenium drugs, via ligand-based modulation of DNA versus protein binding and thus cytotoxic potential, to target distinguishing epigenetic features of cancer cells.

Thiol redox biochemistry: insights from computer simulations

Thiol redox chemical reactions play a key role in a variety of physiological processes, mainly due to the presence of low-molecular-weight thiols and cysteine residues in proteins involved in catalysis and regulation. Specifically, the subtle sensitivity of thiol reactivity to the environment makes the use of simulation techniques extremely valuable for obtaining microscopic insights. In this work we review the application of classical and quantum-mechanical atomistic simulation tools to the investigation of selected relevant issues in thiol redox biochemistry, such as investigations on (1) the protonation state of cysteine in protein, (2) two-electron oxidation of thiols by hydroperoxides, chloramines, and hypochlorous acid, (3) mechanistic and kinetics aspects of the de novo formation of disulfide bonds and thiol-disulfide exchange, (4) formation of sulfenamides, (5) formation of nitrosothiols and transnitrosation reactions, and (6) one-electron oxidation pathways.

Coupling density functional theory to polarizable force fields for efficient and accurate Hamiltonian molecular dynamics simulations

Hybrid molecular dynamics (MD) simulations, in which the forces acting on the atoms are calculated by grid-based density functional theory (DFT) for a solute molecule and by a polarizable molecular mechanics (PMM) force field for a large solvent environment composed of several 10(3)-10(5) molecules, pose a challenge. A corresponding computational approach should guarantee energy conservation, exclude artificial distortions of the electron density at the interface between the DFT and PMM fragments, and should treat the long-range electrostatic interactions within the hybrid simulation system in a linearly scaling fashion. Here we describe a corresponding Hamiltonian DFT/(P)MM implementation, which accounts for inducible atomic dipoles of a PMM environment in a joint DFT/PMM self-consistency iteration. The long-range parts of the electrostatics are treated by hierarchically nested fast multipole expansions up to a maximum distance dictated by the minimum image convention of toroidal boundary conditions and, beyond that distance, by a reaction field approach such that the computation scales linearly with the number of PMM atoms. Short-range over-polarization artifacts are excluded by using Gaussian inducible dipoles throughout the system and Gaussian partial charges in the PMM region close to the DFT fragment. The Hamiltonian character, the stability, and efficiency of the implementation are investigated by hybrid DFT/PMM-MD simulations treating one molecule of the water dimer and of bulk water by DFT and the respective remainder by PMM.

Accurate basis set truncation for wavefunction embedding

Density functional theory (DFT) provides a formally exact framework for performing embedded subsystem electronic structure calculations, including DFT-in-DFT and wavefunction theory-in-DFT descriptions. In the interest of efficiency, it is desirable to truncate the atomic orbital basis set in which the subsystem calculation is performed, thus avoiding high-order scaling with respect to the size of the MO virtual space. In this study, we extend a recently introduced projection-based embedding method [F. R. Manby, M. Stella, J. D. Goodpaster, and T. F. Miller III, J. Chem. Theory Comput. 8, 2564 (2012)] to allow for the systematic and accurate truncation of the embedded subsystem basis set. The approach is applied to both covalently and non-covalently bound test cases, including water clusters and polypeptide chains, and it is demonstrated that errors associated with basis set truncation are controllable to well within chemical accuracy. Furthermore, we show that this approach allows for switching between accurate projection-based embedding and DFT embedding with approximate kinetic energy (KE) functionals; in this sense, the approach provides a means of systematically improving upon the use of approximate KE functionals in DFT embedding.

Perspectives for hybrid ab initio/molecular mechanical simulations of solutions: from complex chemistry to proton-transfer reactions and interfaces

As a consequence of the ongoing development of enhanced computational resources, theoretical chemistry has become an increasingly valuable field for the investigation of a variety of chemical systems. Simulations employing a hybrid quantum mechanical/molecular mechanical (QM/MM) molecular dynamics (MD) technique have been shown to be a particularly promising approach, whenever ultrafast (i.e., picosecond) dynamical properties are to be studied, which are in many cases difficult to access via experimental techniques. Details of the quantum mechanical charge field (QMCF) ansatz, an advanced QM/MM protocol, are discussed and simulation results for various systems ranging from simple ionic hydrates to solvated organic molecules and coordination complexes in solution are presented. A particularly challenging application is the description of proton-transfer reactions in chemical simulations, which is a prerequisite to study acidified and basic systems. The methodical requirements for a combination of the QMCF methodology with a dissociative potential model for the description of the solvent are discussed. Furthermore, the possible extension of QM/MM approaches to solid/liquid interfaces is outlined.

Assigning the EPR fine structure parameters of the Mn(II) centers in Bacillus subtilis oxalate decarboxylase by site-directed mutagenesis and DFT/MM calculations

Oxalate decarboxylase (OxDC) catalyzes the Mn-dependent conversion of the oxalate monoanion into CO2 and formate. EPR-based strategies for investigating the catalytic mechanism of decarboxylation are complicated by the difficulty of assigning the signals associated with the two Mn(II) centers located in the N- and C-terminal cupin domains of the enzyme. We now report a mutational strategy that has established the assignment of EPR fine structure parameters to each of these Mn(II) centers at pH 8.5. These experimental findings are also used to assess the performance of a multistep strategy for calculating the zero-field splitting parameters of protein-bound Mn(II) ions. Despite the known sensitivity of calculated D and E values to the computational approach, we demonstrate that good estimates of these parameters can be obtained using cluster models taken from carefully optimized DFT/MM structures. Overall, our results provide new insights into the strengths and limitations of theoretical methods for understanding electronic properties of protein-bound Mn(II) ions, thereby setting the stage for future EPR studies on the electronic properties of the Mn(II) centers in OxDC and site-specific variants.

An extensible interface for QM/MM molecular dynamics simulations with AMBER

We present an extensible interface between the AMBER molecular dynamics (MD) software package and electronic structure software packages for quantum mechanical (QM) and mixed QM and classical molecular mechanical (MM) MD simulations within both mechanical and electronic embedding schemes. With this interface, ab initio wave function theory and density functional theory methods, as available in the supported electronic structure software packages, become available for QM/MM MD simulations with AMBER. The interface has been written in a modular fashion that allows straight forward extensions to support additional QM software packages and can easily be ported to other MD software. Data exchange between the MD and QM software is implemented by means of files and system calls or the message passing interface standard. Based on extensive tests, default settings for the supported QM packages are provided such that energy is conserved for typical QM/MM MD simulations in the microcanonical ensemble. Results for the free energy of binding of calcium ions to aspartate in aqueous solution comparing semiempirical and density functional Hamiltonians are shown to demonstrate features of this interface.

Electron flow in multiheme bacterial cytochromes is a balancing act between heme electronic interaction and redox potentials

The naturally widespread process of electron transfer from metal reducing bacteria to extracellular solid metal oxides entails unique biomolecular machinery optimized for long-range electron transport. To perform this function efficiently, microorganisms have adapted multiheme c-type cytochromes to arrange heme cofactors into wires that cooperatively span the cellular envelope, transmitting electrons along distances greater than 100 Å. Implications and opportunities for bionanotechnological device design are self-evident. However, at the molecular level, how these proteins shuttle electrons along their heme wires, navigating intraprotein intersections and interprotein interfaces efficiently, remains a mystery thus far inaccessible to experiment. To shed light on this critical topic, we carried out extensive quantum mechanics/molecular mechanics simulations to calculate stepwise heme-to-heme electron transfer rates in the recently crystallized outer membrane deca-heme cytochrome MtrF. By solving a master equation for electron hopping, we estimate an intrinsic, maximum possible electron flux through solvated MtrF of 10(4)-10(5) s(-1), consistent with recently measured rates for the related multiheme protein complex MtrCAB. Intriguingly, our calculations show that the rapid electron transport through MtrF is the result of a clear correlation between heme redox potential and the strength of electronic coupling along the wire: thermodynamically uphill steps occur only between electronically well-connected stacked heme pairs. This observation suggests that the protein evolved to harbor low-potential hemes without slowing down electron flow. These findings are particularly profound in light of the apparently well-conserved staggered cross-heme wire structural motif in functionally related outer membrane proteins.