An interoperable implementation of collective-variable based enhanced sampling methods in extended phase space within the OpenMM package.

Collective variable (CV)-based enhanced sampling techniques are widely used today for accelerating barrier-crossing events in molecular simulations. A class of these methods, which includes temperature accelerated molecular dynamics (TAMD)/driven-adiabatic free energy dynamics (d-AFED), unified free energy dynamics (UFED), and temperature accelerated sliced sampling (TASS), uses an extended variable formalism to achieve quick exploration of conformational space. These techniques are powerful, as they enhance the sampling of a large number of CVs simultaneously compared to other techniques. Extended variables are kept at a much higher temperature than the physical temperature by ensuring adiabatic separation between the extended and physical subsystems and employing rigorous thermostatting. In this work, we present a computational platform to perform extended phase space enhanced sampling simulations using the open-source molecular dynamics engine OpenMM. The implementation allows users to have interoperability of sampling techniques, as well as employ state-of-the-art thermostats and multiple time-stepping. This work also presents protocols for determining the critical parameters and procedures for reconstructing high-dimensional free energy surfaces. As a demonstration, we present simulation results on the high dimensional conformational landscapes of the alanine tripeptide in vacuo, tetra-N-methylglycine (tetra-sarcosine) peptoid in implicit solvent, and the Trp-cage mini protein in explicit water.

Hybrid Functional and Plane Waves based Ab Initio Molecular Dynamics Study of the Aqueous Fe2+ /Fe3+ Redox Reaction.

Kohn-Sham density functional theory and plane wave basis set based ab initio molecular dynamics (AIMD) simulation is a powerful tool for studying complex reactions in solutions, such as electron transfer (ET) reactions involving Fe2+ /Fe3+ ions in water. In most cases, such simulations are performed using density functionals at the level of Generalized Gradient Approximation (GGA). The challenge in modelling ET reactions is the poor quality of GGA functionals in predicting properties of such open-shell systems due to the inevitable self-interaction error (SIE). While hybrid functionals can minimize SIE, standard plane-wave based AIMD at that level of theory is typically 150 times slower than GGA for systems containing ∼100 atoms. Among several approaches reported to speed-up AIMD simulations with hybrid functionals, the noise-stabilized MD (NSMD) procedure, together with the use of localized orbitals to compute the required exchange integrals, is an attractive option. In this work, we demonstrate the application of the NSMD approach for studying the Fe2+ /Fe3+ redox reaction in water. It is shown here that long AIMD trajectories at the level of hybrid density functionals can be obtained using this approach. Redox properties of the aqueous Fe2+ /Fe3+ system computed from these simulations are compared with the available experimental data for validation.

Temperature Accelerated Sliced Sampling to Probe Ligand Dissociation from Protein.

Modeling ligand unbinding in proteins to estimate the free energy of binding and probing the mechanism presents several challenges. They primarily pertain to the entropic bottlenecks resulting from protein and solvent conformations. While exploring the unbinding processes using enhanced sampling techniques, very long simulations are required to sample all of the conformational states as the system gets trapped in local free energy minima along transverse coordinates. Here, we demonstrate that temperature accelerated sliced sampling (TASS) is an ideal approach to overcome some of the difficulties faced by conventional sampling methods in studying ligand unbinding. Using TASS, we study the unbinding of avibactam inhibitor molecules from the Class C ß-lactamase (CBL) active site. Extracting CBL-avibactam unbinding free energetics, unbinding pathways, and identifying critical interactions from the TASS simulations are demonstrated.

Probing the general base for DNA polymerization in telomerase: a molecular dynamics investigation.

Telomerase is an RNA-dependent DNA polymerase that plays a role in the maintenance of the 3' end of the eukaryotic chromosome, known as a telomere, by catalyzing the DNA polymerization reaction in cancer and embryonic stem cells. The detailed molecular details of the DNA polymerization by telomerase, especially the general base for deprotonating the terminal 3'-hydroxyl, which triggers the chemical reaction, remain elusive. We conducted a computational investigation using hybrid quantum mechanical/molecular mechanical (QM/MM) molecular dynamics (MD) simulations to probe the detailed mechanism of the reaction. Our simulations started with the telomerase:RNA:DNA:dNTP ternary complex, and by using enhanced sampling QM/MM MD simulations, we probed the general base involved directly in the polymerization. We report the participation of an aspartate (Asp344) coordinated to Mg and an active site water molecule, jointly acting as a base during nucleic acid addition. The Asp344 residue remains transiently protonated during the course of the reaction, and later it deprotonates by transferring its proton to the water at the end of the reaction.

Exploration of high dimensional free energy landscapes by a combination of temperature-accelerated sliced sampling and parallel biasing.

Temperature-accelerated sliced sampling (TASS) is an enhanced sampling method for achieving accelerated and controlled exploration of high-dimensional free energy landscapes in molecular dynamics simulations. With the aid of umbrella bias potentials, the TASS method realizes a controlled exploration and divide-and-conquer strategy for computing high-dimensional free energy surfaces. In TASS, diffusion of the system in the collective variable (CV) space is enhanced with the help of metadynamics bias and elevated-temperature of the auxiliary degrees of freedom (DOF) that are coupled to the CVs. Usually, a low-dimensional metadynamics bias is applied in TASS. In order to further improve the performance of TASS, we propose here to use a highdimensional metadynamics bias, in the same form as in a parallel bias metadynamics scheme. Here, a modified reweighting scheme, in combination with artificial neural network is used for computing unbiased probability distribution of CVs and projections of high-dimensional free energy surfaces. We first validate the accuracy and efficiency of our method in computing the four-dimensional free energy landscape for alanine tripeptide in vacuo. Subsequently, we employ the approach to calculate the eight-dimensional free energy landscape of alanine pentapeptide in vacuo. Finally, the method is applied to a more realistic problem wherein we compute the broad four-dimensional free energy surface corresponding to the deacylation of a drug molecule which is covalently complexed with a ß-lactamase enzyme. We demonstrate that using parallel bias in TASS improves the efficiency of exploration of high-dimensional free energy landscapes.

Improving the scaling and performance of multiple time stepping-based molecular dynamics with hybrid density functionals.

Density functionals at the level of the generalized gradient approximation (GGA) and a plane-wave basis set are widely used today to perform ab initio molecular dynamics (AIMD) simulations. Going up in the ladder of accuracy of density functionals from GGA (second rung) to hybrid density functionals (fourth rung) is much desired pertaining to the accuracy of the latter in describing structure, dynamics, and energetics of molecular and condensed matter systems. On the other hand, hybrid density functional based AIMD simulations are about two orders of magnitude slower than GGA based AIMD for systems containing ~100 atoms using ~100 compute cores. Two methods, namely MTACE and s-MTACE, based on a multiple time step integrator and adaptively compressed exchange operator formalism are able to provide a speed-up of about 7-9 in performing hybrid density functional based AIMD. In this work, we report an implementation of these methods using a task-group based parallelization within the CPMD program package, with the intention to take advantage of the large number of compute cores available on modern high-performance computing platforms. We present here the boost in performance achieved through this algorithm. This work also identifies the computational bottleneck in the s-MTACE method and proposes a way to overcome it.

Boosting the conformational sampling by combining replica exchange with solute tempering and well-sliced metadynamics.

Methods that combine collective variable (CV) based enhanced sampling and global tempering approaches are used in speeding-up the conformational sampling and free energy calculation of large and soft systems with a plethora of energy minima. In this paper, a new method of this kind is proposed in which the well-sliced metadynamics approach (WSMTD) is united with replica exchange with solute tempering (REST2) method. WSMTD employs a divide-and-conquer strategy wherein high-dimensional slices of a free energy surface are independently sampled and combined. The method enables one to accomplish a controlled exploration of the CV-space with a restraining bias as in umbrella sampling, and enhance-sampling of one or more orthogonal CVs using a metadynamics like bias. The new hybrid method proposed here enables boosting the sampling of more slow degrees of freedom in WSMTD simulations, without the need to specify associated CVs, through a replica exchange scheme within the framework of REST2. The high-dimensional slices of the probability distributions of CVs computed from the united WSMTD and REST2 simulations are subsequently combined using the weighted histogram analysis method to obtain the free energy surface. We show that the new method proposed here is accurate, improves the conformational sampling, and achieves quick convergence in free energy estimates. We demonstrate this by computing the conformational free energy landscapes of solvated alanine tripeptide and Trp-cage mini protein in explicit water.

Mean force based temperature accelerated sliced sampling: Efficient reconstruction of high dimensional free energy landscapes.

Temperature accelerated sliced sampling (TASS) is an efficient method to compute high dimensional free energy landscapes. The original TASS method employs the weighted histogram analysis method (WHAM) which is an iterative post-processing to reweight and stitch high dimensional probability distributions in sliced windows that are obtained in the presence of restraining biases. The WHAM necessitates that TASS windows lie close to each other for proper overlap of distributions and span the collective variable space of interest. On the other hand, increase in number of TASS windows implies more number of simulations, and thus it affects the efficiency of the method. To overcome this problem, we propose herein a new mean-force (MF) based reweighting scheme called TASS-MF, which enables accurate computation with a fewer number of windows devoid of the WHAM post-processing. Application of the technique is demonstrated for alanine di- and tripeptides in vacuo to compute their two- and four-dimensional free energy landscapes, the latter of which is formidable in conventional umbrella sampling and metadynamics. The landscapes are computed within a kcal mol-1 accuracy, ensuring a safe usage for broad applications in computational chemistry.

Mechanism of Nucleotidyltransfer Reaction and Role of Mg2+ Ion in Sugar Nucleotidyltransferases.

Sugar nucleotidyl transferases (SNTs) catalyze nucleotidyltransfer reactions to form sugar-nucleotides and pyrophosphate in the presence of two Mg2+ ions (Mg2+A and Mg2+B). We unveil the mechanism and free energetics of nucleotidyl transfer reaction in an SNT called GlmU through hybrid quantum mechanics-molecular mechanics molecular dynamics simulations and free energy calculations. The study identifies the roles of the active site residues and the Mg2+ ions in catalyzing the reaction. Of great significance, we are able to compare the free energy barrier for the reaction with that for the Mg2+-assisted release of the product (i.e., pyrophosphate) into the solution, shedding light on the general mechanistic and kinetic aspects of catalysis by SNTs.

Efficient computation of free energy surfaces of chemical reactions using ab initio molecular dynamics with hybrid functionals and plane waves.

Ab initio molecular dynamics (AIMD) simulations employing density functional theory (DFT) and plane waves are routinely carried out using density functionals at the level of generalized gradient approximation (GGA). AIMD simulations employing hybrid density functionals are of great interest as it offers a more accurate description of structural and dynamic properties than the GGA functionals. However, the computational cost for carrying out calculations using hybrid functionals and plane wave basis set is at least two orders of magnitude higher than that using GGA functionals. Recently, we proposed a strategy that combined the adaptively compressed exchange operator formulation and the multiple time step integration scheme to reduce the computational cost by an order of magnitude [J. Chem. Phys. 151, 151102 (2019)]. In this work, we demonstrate the application of this method to study chemical reactions, in particular, formamide hydrolysis in an alkaline aqueous medium. By actuating our implementation with the well-sliced metadynamics scheme, we can compute the two-dimensional free energy surface of this reaction at the level of hybrid-DFT. This work also investigates the accuracy of the PBE0 (hybrid) and the PBE (GGA) functionals in predicting the free energetics of this chemical reaction.

Elucidating the Molecular Basis of Avibactam-Mediated Inhibition of Class†A ß-Lactamases.

Disseminating antibiotic resistance rendered by bacteria against the widely used ß-lactam antibiotics is a serious concern for public health care. The development of inhibitors for drug-resistant ß-lactamase enzymes is vital to combat this rapidly escalating problem. Recently, the U.S. Food and Drug Administration approved a non-ß-lactam inhibitor called avibactam for the treatment of complicated intra-abdominal and urinary tract infections caused by drug-resistant Gram-negative bacteria. This work sheds light on the molecular origin of the inhibitory effect of avibactam against the drug-resistant CTX-M variant of class A ß-lactamases. In particular, we probed the structural evolution, dynamics features, and energetics along the acylation and deacylation reaction pathways through enhanced sampling molecular dynamics methods and free-energy calculations. We scrutinized the roles of active site residues, the nature of the carbamoyl linkage formed in the inhibitor-enzyme covalent intermediate, and other structural features of the inhibitor molecule. By unraveling the reasons behind the inhibition of all the deacylation routes, we can explain various experimental structural and kinetics data, and propose a way to design new inhibitors based on the ß-lactam framework.

Speeding-up ab initio molecular dynamics with hybrid functionals using adaptively compressed exchange operator based multiple timestepping.

Ab initio molecular dynamics (AIMD) simulations using hybrid density functionals and plane waves are of great interest owing to the accuracy of this approach in treating condensed matter systems. On the other hand, such AIMD calculations are not routinely carried out since the computational cost involved in applying the Hartree-Fock exchange operator is very high. In this work, we make use of a strategy that combines adaptively compressed exchange operator formulation and multiple time step integration to significantly reduce the computational cost of these simulations. We demonstrate the efficiency of this approach for a realistic condensed matter system.

Molecular insights into avibactam mediated class C ß-lactamase inhibition: competition between reverse acylation and hydrolysis through desulfation.

Avibactam is one of the promising next generation ß-lactamase inhibitors due to its exceptional inhibition against wide-spectrum serine ß-lactamases. The unusual reversible acylation mechanism has particularly gained interest to explain the inhibition mechanism of avibactam. We explore the mechanism of acylation and deacylation involving avibactam in class-C ß-lactamases (CBLs) through hybrid quantum mechanical/molecular mechanical (QM/MM) enhanced sampling molecular dynamics (MD) simulations. Based on these computations, we probe the kinetic stability of the acyl-enzyme complex formed by avibactam and CBLs, thereby gaining molecular level insights into the avibactam-mediated inhibition of CBLs.

Enhanced sampling and free energy calculations with hybrid functionals and plane waves for chemical reactions.

Plane wave basis sets offer many advantages in ab initio molecular dynamics due to their efficiency and simplicity. In combination with hybrid density functionals, they become computationally expensive due to the evaluation of the Hartree-Fock exchange energy. The computational cost can be significantly reduced by screening the Kohn-Sham orbital products after localizing the orbitals in real space. However, such a procedure introduces apparent errors in the wavefunctions and nuclear forces resulting in unstable dynamics. It is shown here that a noise-stabilized dynamics approach can overcome this problem and at the same time permits using insufficiently converged wavefunctions for evaluating atomic forces. In this way, we achieve significant speed up even for a small system containing about 100 atoms. After benchmarking the accuracy and efficiency of this approach, we use it in combination with well-sliced metadynamics to compute the free energy barrier of formamide hydrolysis in alkaline aqueous medium. These results provide insight into the error of the Perdew-Burke-Ernzerhof functional in predicting the free energy barrier for hydrolysis reactions in water.

Hydrolysis of cephalexin and meropenem by New Delhi metallo-ß-lactamase: the substrate protonation mechanism is drug dependent.

Emergence of antibiotic resistance due to New Delhi metallo-ß-lactamase (NDM-1) bacterial enzymes is of great concern due to their ability to hydrolyze a wide range of antibiotics. There are ongoing efforts to obtain the atomistic details of the hydrolysis mechanism in order to develop inhibitors for NDM-1. In particular, it remains elusive how drug molecules of different families of antibiotics are hydrolyzed by NDM-1 in an efficient manner. Here we report the detailed molecular mechanism of NDM-1 catalyzed hydrolysis of cephalexin, a cephalosporin family drug, and meropenem, a carbapenem family drug. This study employs molecular dynamics (MD) simulations using hybrid quantum mechanical/molecular mechanical (QM/MM) methods at the density functional theory (DFT) level, based on which reaction pathways and the associated free energies are obtained. We find that the mechanism and the free energy barrier for the ring-opening step are the same for both the drug molecules, while the subsequent protonation step differs. In particular, we observe that the mechanism of the protonation step depends on the R2 group of the drug molecule. Our simulations show that allylic carbon protonation occurs in the case of the cephalexin drug molecule where Lys211 is the proton donor, and the proton transfer occurs via a water chain formed (only) at the ring-opened intermediate structure. Based on the free energy profiles, the overall kinetics of drug hydrolysis is discussed. Finally, we show that the proposed mechanisms and free energy profiles could explain various experimental observations.

CPMD/GULP QM/MM interface for modeling periodic solids: Implementation and its application in the study of Y-zeolite supported Rhn clusters.

We report here the development of hybrid quantum mechanics/molecular mechanics (QM/MM) interface between the plane-wave density functional theory based CPMD code and the empirical force-field based GULP code for modeling periodic solids and surfaces. The hybrid QM/MM interface is based on the electrostatic coupling between QM and MM regions. The interface is designed for carrying out full relaxation of all the QM and MM atoms during geometry optimizations and molecular dynamics simulations, including the boundary atoms. Both Born-Oppenheimer and Car-Parrinello molecular dynamics schemes are enabled for the QM part during the QM/MM calculations. This interface has the advantage of parallelization of both the programs such that the QM and MM force evaluations can be carried out in parallel to model large systems. The interface program is first validated for total energy conservation and parallel scaling performance is benchmarked. Oxygen vacancy in α-cristobalite is then studied in detail and the results are compared with a fully QM calculation and experimental data. Subsequently, we use our implementation to investigate the structure of rhodium cluster (Rhn ; n = 2 to 6) formed from Rh(C2 H4 )2 complex adsorbed within a cavity of Y-zeolite in a reducible atmosphere of H2 gas. © 2016 Wiley Periodicals, Inc.

Sampling free energy surfaces as slices by combining umbrella sampling and metadynamics.

Metadynamics (MTD) is a very powerful technique to sample high-dimensional free energy landscapes, and due to its self-guiding property, the method has been successful in studying complex reactions and conformational changes. MTD sampling is based on filling the free energy basins by biasing potentials and thus for cases with flat, broad, and unbound free energy wells, the computational time to sample them becomes very large. To alleviate this problem, we combine the standard Umbrella Sampling (US) technique with MTD to sample orthogonal collective variables (CVs) in a simultaneous way. Within this scheme, we construct the equilibrium distribution of CVs from biased distributions obtained from independent MTD simulations with umbrella potentials. Reweighting is carried out by a procedure that combines US reweighting and Tiwary-Parrinello MTD reweighting within the Weighted Histogram Analysis Method (WHAM). The approach is ideal for a controlled sampling of a CV in a MTD simulation, making it computationally efficient in sampling flat, broad, and unbound free energy surfaces. This technique also allows for a distributed sampling of a high-dimensional free energy surface, further increasing the computational efficiency in sampling. We demonstrate the application of this technique in sampling high-dimensional surface for various chemical reactions using ab initio and QM/MM hybrid molecular dynamics simulations. Further, to carry out MTD bias reweighting for computing forward reaction barriers in ab initio or QM/MM simulations, we propose a computationally affordable approach that does not require recrossing trajectories. © 2016 Wiley Periodicals, Inc.

A potential with low point charges for pure siliceous zeolites.

A modified CHARMM force-field (ZHB potential) with low point charges for silica was previously proposed by Zimmerman et al. (J. Chem. Theory Comput. 2011, 7, 1695). The ZHB potential is advantageous for quantum mechanics/molecular mechanics simulations as it minimizes the electron spill-out problems. In the same spirit, here we propose a modified ZHB potential (MZHB) by reformulating its bonding potential, while retaining the nonbonding potential as in the ZHB force-field. We show that several structural and dynamic properties of silica, like the IR spectrum, distribution functions, mechanical properties, and negative thermal expansion computed using the MZHB potential agree well with experimental data. Further, transferability of MZHB is also tested for reproducing the crystallographic structures of several polymorphs of silica. © 2015 Wiley Periodicals, Inc.

Interaction potential models for bulk ZnS, ZnS nanoparticle, and ZnS nanoparticle-PMMA from first-principles.

An ab initio derived transferable polarizable force-field has been developed for Zinc sulphide (ZnS) nanoparticle (NP) and ZnS NP-PMMA nanocomposite. The structure and elastic constants of bulk ZnS using the new force-field are within a few percent of experimental observables. The new force-field show remarkable ability to reproduce structures and nucleation energies of nanoclusters (Zn1S1-Zn12S12) as validated with that of the density functional theory calculations. A qualitative agreement of the radial distribution functions of Zn-O, in a ZnS nanocluster-PMMA system, obtained using molecular mechanics molecular dynamics (MD) and ab initio MD (AIMD) simulations indicates that the ZnS-PMMA interaction through Zn-O bonding is explained satisfactorily by our force-field.

Challenges in modelling homogeneous catalysis: new answers from ab initio molecular dynamics to the controversy over the Wacker process.

We present here a review of the mechanistic studies of the Wacker process stressing the long controversy about the key reaction steps. We give an overview of the previous experimental and theoretical studies on the topic. Then we describe the importance of the most recent Ab Initio Molecular Dynamics (AIMD) calculations in modelling organometallic reactivity in water. As a prototypical example of homogeneous catalytic reactions, the Wacker process poses serious challenges to modelling. The adequate description of the multiple role of the water solvent is very difficult by using static quantum chemical approaches including cluster and continuum solvent models. In contrast, such reaction systems are suitable for AIMD, and by combining with rare event sampling techniques, the method provides reaction mechanisms and the corresponding free energy profiles. The review also highlights how AIMD has helped to obtain a novel understanding of the mechanism and kinetics of the Wacker process.