Multi-Level Ewald: A hybrid multigrid / Fast Fourier Transform approach to the electrostatic particle-mesh problem

STORMM: Structure and topology replica molecular mechanics for chemical simulations

The Structure and TOpology Replica Molecular Mechanics (STORMM) code is a next-generation molecular simulation engine and associated libraries optimized for performance on fast, vectorized central processor units and graphics processing units (GPUs) with independent memory and tens of thousands of threads. STORMM is built to run thousands of independent molecular mechanical calculations on a single GPU with novel implementations that tune numerical precision, mathematical operations, and scarce on-chip memory resources to optimize throughput. The libraries are built around accessible classes with detailed documentation, supporting fine-grained parallelism and algorithm development as well as copying or swapping groups of systems on and off of the GPU. A primary intention of the STORMM libraries is to provide developers of atomic simulation methods with access to a high-performance molecular mechanics engine with extensive facilities to prototype and develop bespoke tools aimed toward drug discovery applications. In its present state, STORMM delivers molecular dynamics simulations of small molecules and small proteins in implicit solvent with tens to hundreds of times the throughput of conventional codes. The engineering paradigm transforms two of the most memory bandwidth-intensive aspects of condensed-phase dynamics, particle-mesh mapping, and valence interactions, into compute-bound problems for several times the scalability of existing programs. Numerical methods for compressing and streamlining the information present in stored coordinates and lookup tables are also presented, delivering improved accuracy over methods implemented in other molecular dynamics engines. The open-source code is released under the MIT license.

How to train a neural network potential

The introduction of modern Machine Learning Potentials (MLPs) has led to a paradigm change in the development of potential energy surfaces for atomistic simulations. By providing efficient access to energies and forces, they allow us to perform large-scale simulations of extended systems, which are not directly accessible by demanding first-principles methods. In these simulations, MLPs can reach the accuracy of electronic structure calculations, provided that they have been properly trained and validated using a suitable set of reference data. Due to their highly flexible functional form, the construction of MLPs has to be done with great care. In this Tutorial, we describe the necessary key steps for training reliable MLPs, from data generation via training to final validation. The procedure, which is illustrated for the example of a high-dimensional neural network potential, is general and applicable to many types of MLPs.

Acute phase α1-acid glycoprotein as a siderophore-capturing component of the human plasma: A molecular modeling study

Siderophores are ferric ion-specific organic compounds that are used by bacteria and fungi to secure their iron supply when infecting target organisms. There are a few proteins in the human body, named siderocalins, which bind these important virulence factors and so starve microorganisms of iron. In this study, we analyzed in silico if serum α₁-acid glycoprotein (AAG), the major acute phase lipocalin component of the human plasma, could functionally belong to this group. The real biological function of AAG is elusive and its concentration substantially increases in response to pathological stimuli, including bacterial infections. We computationally evaluated the potential binding of nine microbial siderophores into the β-barrel cavity of AAG and compared the results with the corresponding experimental data reported for siderophore-neutrophil gelatinase-associated lipocalin complexes. According to the results, petrobactin and Fe-BisHaCam are putative candidates to be recognized by this protein. It is proposed that AAG may function as a siderophore capturing component of the innate immune system being able to neutralize bacterial iron chelators not recognized by other siderocalins.

SAMPL7: Host-guest binding prediction by molecular dynamics and quantum mechanics

Statistical Assessment of Modeling of Proteins and Ligands (SAMPL) challenges provide routes to compare chemical quantities determined using computational chemistry approaches to experimental measurements that are shared after the competition. For this effort, several computational methods have been used to calculate the binding energies of Octa Acid (OA) and exo-Octa Acid (exoOA) host-guest systems for SAMPL7. The initial poses for molecular dynamics (MD) were generated by molecular docking. Binding free energy calculations were performed using molecular mechanics combined with Poisson-Boltzmann or generalized Born surface area solvation (MMPBSA/MMGBSA) approaches. The factors that affect the utility of the MMPBSA/MMGBSA approaches including solvation, partial charge, and solute entropy models were also analyzed. In addition to MD calculations, quantum mechanics (QM) calculations were performed using several different density functional theory (DFT) approaches. From SAMPL6 results, B3PW91-D3 was found to overestimate binding energies though it was effective for geometry optimizations, so it was considered for the DFT geometry optimizations in the current study, with single-point energy calculations carried out with B2PLYP-D3 with double-, triple-, and quadruple-ζ level basis sets. Accounting for dispersion effects, and solvation models was deemed essential for the predictions. MMGBSA and MMPBSA correlated better to experiment when used in conjunction with an empirical/linear correction.

Role of hydrogen bond capacity of solvents in reactions of amines with CO2: A computational study

Various computational methods were employed to investigate the zwitterion formation, a critical step for the reaction of monoethanolamine with CO₂, in five solvents (water, monoethanolamine, propylamine, methanol and chloroform) to probe the effect of hydrogen bond capacity of solvents on the reaction of amine with CO₂ occurring in the amine-based CO₂ capture process. The results indicate that the zwitterion can be formed in all considered solvents except chloroform. For two pairs of solvents (methanol and monoethanolamine, propylamine and chloroform) with similar dielectric constant but different hydrogen bond capacity, the solvents with higher hydrogen bond capacity (monoethanolamine and propylamine) facilitate the zwitterion formation. More importantly, kinetics parameters such as activation free energy for the zwitterion formation are more relevant to the hydrogen bond capacity than to dielectric constant of the considered solvents, clarifying the hydrogen bond capacity could be more important than dielectric constant in determining the kinetics of monoethanolamine with CO₂.

Building better enzymes: Molecular basis of improved non-natural nucleobase incorporation by an evolved DNA polymerase

Obtaining semisynthetic microorganisms that exploit the information density of "hachimoji" DNA requires access to engineered DNA polymerases. A KlenTaq variant has been reported that incorporates the "hachimoji" P:Z nucleobase pair with a similar efficiency to that seen for Watson-Crick nucleobase incorporation by the wild type (WT) KlenTaq DNA polymerase. The variant polymerase differs from WT KlenTaq by only four amino acid substitutions, none of which are located within the active site. We now report molecular dynamics (MD) simulations on a series of binary complexes aimed at elucidating the contributions of the four amino acid substitutions to altered catalytic activity. These simulations suggest that WT KlenTaq is insufficiently flexible to be able to bind AEGIS DNA correctly, leading to the loss of key protein/DNA interactions needed to position the binary complex for efficient incorporation of the "hachimoji" Z nucleobase. In addition, we test literature hypotheses about the functional roles of each amino acid substitution and provide a molecular description of how individual residue changes contribute to the improved activity of the KlenTaq variant. We demonstrate that MD simulations have a clear role to play in systematically screening DNA polymerase variants capable of incorporating different types of nonnatural nucleobases thereby limiting the number that need to be characterized by experiment. It is now possible to build DNA molecules containing nonnatural nucleobase pairs in addition to A:T and G:C. Exploiting this development in synthetic biology requires engineered DNA polymerases that can replicate nonnatural nucleobase pairs. Computational studies on a DNA polymerase variant reveal how amino acid substitutions outside of the active site yield an enzyme that replicates nonnatural nucleobase pairs with high efficiency. This work will facilitate efforts to obtain bacteria possessing an expanded genetic alphabet.

Incremental update of electrostatic interactions in adaptively restrained particle simulations

The computation of long-range potentials is one of the demanding tasks in Molecular Dynamics. During the last decades, an inventive panoply of methods was developed to reduce the CPU time of this task. In this work, we propose a fast method dedicated to the computation of the electrostatic potential in adaptively restrained systems. We exploit the fact that, in such systems, only some particles are allowed to move at each timestep. We developed an incremental algorithm derived from a multigrid-based alternative to traditional Fourier-based methods. Our algorithm was implemented inside LAMMPS, a popular molecular dynamics simulation package. We evaluated the method on different systems. We showed that the new algorithm's computational complexity scales with the number of active particles in the simulated system, and is able to outperform the well-established Particle Particle Particle Mesh (P3M) for adaptively restrained simulations. © 2018 Wiley Periodicals, Inc.

Quantum mechanical force fields for condensed phase molecular simulations

Molecular simulations are powerful tools for providing atomic-level details into complex chemical and physical processes that occur in the condensed phase. For strongly interacting systems where quantum many-body effects are known to play an important role, density-functional methods are often used to provide the model with the potential energy used to drive dynamics. These methods, however, suffer from two major drawbacks. First, they are often too computationally intensive to practically apply to large systems over long time scales, limiting their scope of application. Second, there remain challenges for these models to obtain the necessary level of accuracy for weak non-bonded interactions to obtain quantitative accuracy for a wide range of condensed phase properties. Quantum mechanical force fields (QMFFs) provide a potential solution to both of these limitations. In this review, we address recent advances in the development of QMFFs for condensed phase simulations. In particular, we examine the development of QMFF models using both approximate and ab initio density-functional models, the treatment of short-ranged non-bonded and long-ranged electrostatic interactions, and stability issues in molecular dynamics calculations. Example calculations are provided for crystalline systems, liquid water, and ionic liquids. We conclude with a perspective for emerging challenges and future research directions.

Consecutive non-natural PZ nucleobase pairs in DNA impact helical structure as seen in 50 μs molecular dynamics simulations

Z

Little is known about the influence of multiple consecutive 'non-standard' ( , 6-amino-5-nitro-2(1H)-pyridone, and , 2-amino-imidazo[1,2-a]-1,3,5-triazin-4(8H)-one) nucleobase pairs on the structural parameters of duplex DNA. nucleobase pairs follow standard rules for Watson-Crick base pairing but have rearranged hydrogen bonding donor and acceptor groups. Using the X-ray crystal structure as a starting point, we have modeled the motions of a DNA duplex built from a self-complementary oligonucleotide (5΄-CTTATPPPZZZATAAG-3΄) in water over a period of 50 μs and calculated DNA local parameters, step parameters, helix parameters, and major/minor groove widths to examine how the presence of multiple, consecutive nucleobase pairs might impact helical structure. In these simulations, the -containing DNA duplex exhibits a significantly wider major groove and greater average values of stagger, slide, rise, twist and h-rise than observed for a 'control' oligonucleotide in which nucleobase pairs are replaced by . The molecular origins of these structural changes are likely associated with at least two differences between and . First, the electrostatic properties of differ from in terms of density distribution and dipole moment. Second, differences are seen in the base stacking of pairs in dinucleotide steps, arising from energetically favorable stacking of the nitro group in with π-electrons of the adjacent base.

Ambient-Potential Composite Ewald Method for ab Initio Quantum Mechanical/Molecular Mechanical Molecular Dynamics Simulation

A new approach for performing Particle Mesh Ewald in ab initio quantum mechanical/molecular mechanical (QM/MM) simulations with extended atomic orbital basis sets is presented. The new approach, the Ambient-Potential Composite Ewald (CEw) method, does not perform the QM/MM interaction with Mulliken charges nor electrostatically fit charges. Instead the nuclei and electron density interact directly with the MM environment, but in a manner that avoids the use of dense Fourier transform grids. By performing the electrostatics with the underlying QM density, the CEw method avoids self-consistent field instabilities that have been encountered with simple charge mapping procedures. Potential of mean force (PMF) profiles of the p-nitrophenyl phosphate dissociation reaction in explicit solvent are computed from PBE0/6-31G* QM/MM molecular dynamics simulations with various electrostatic protocols. The CEw profiles are shown to be stable with respect to real-space Ewald cutoff, whereas the PMFs computed from truncated and switched electrostatics produce artifacts. PBE0/6-311G**, AM1/d-PhoT, and DFTB2 QM/MM simulations are performed to generate two-dimensional PMF profiles of the phosphoryl transesterification reactions with ethoxide and phenoxide leaving groups. The semiempirical models incorrectly produce a concerted ethoxide mechanism, whereas PBE0 correctly produces a stepwise mechanism. The ab initio reaction barriers agree more closely to experiment than the semiempirical models. The failure of Mulliken-charge QM/MM-Ewald is analyzed.

Multipolar Ewald methods, 1: theory, accuracy, and performance

The Ewald, Particle Mesh Ewald (PME), and Fast Fourier–Poisson (FFP) methods are developed for systems composed of spherical multipole moment expansions. A unified set of equations is derived that takes advantage of a spherical tensor gradient operator formalism in both real space and reciprocal space to allow extension to arbitrary multipole order. The implementation of these methods into a novel linear-scaling modified “divide-and-conquer” (mDC) quantum mechanical force field is discussed. The evaluation times and relative force errors are compared between the three methods, as a function of multipole expansion order. Timings and errors are also compared within the context of the quantum mechanical force field, which encounters primary errors related to the quality of reproducing electrostatic forces for a given density matrix and secondary errors resulting from the propagation of the approximate electrostatics into the self-consistent field procedure, which yields a converged, variational, but nonetheless approximate density matrix. Condensed-phase simulations of an mDC water model are performed with the multipolar PME method and compared to an electrostatic cutoff method, which is shown to artificially increase the density of water and heat of vaporization relative to full electrostatic treatment.

Progress in developing Poisson-Boltzmann equation solvers

This review outlines the recent progress made in developing more accurate and efficient solutions to model electrostatics in systems comprised of bio-macromolecules and nano-objects, the last one referring to objects that do not have biological function themselves but nowadays are frequently used in biophysical and medical approaches in conjunction with bio-macromolecules. The problem of modeling macromolecular electrostatics is reviewed from two different angles: as a mathematical task provided the specific definition of the system to be modeled and as a physical problem aiming to better capture the phenomena occurring in the real experiments. In addition, specific attention is paid to methods to extend the capabilities of the existing solvers to model large systems toward applications of calculations of the electrostatic potential and energies in molecular motors, mitochondria complex, photosynthetic machinery and systems involving large nano-objects.

Electrostatics interactions in classical simulations

Electrostatic interactions are crucial for both the accuracy and performance of atomistic biomolecular simulations. In this chapter we review well-established methods and current developments aiming at efficiency and accuracy. Specifically, we review the classical Ewald summations, particle-particle particle-method particle-method Ewald algorithms, multigrid, fast multipole, and local methods. We also highlight some recent developments targeting more accurate, yet classical, representation of the molecular charge distribution.

Polarizable Atomic Multipole X-Ray Refinement: Particle Mesh Ewald Electrostatics for Macromolecular Crystals

Refinement of macromolecular models from X-ray crystallography experiments benefits from prior chemical knowledge at all resolutions. As the quality of the prior chemical knowledge from quantum or classical molecular physics improves, in principle so will resulting structural models. Due to limitations in computer performance and electrostatic algorithms, commonly used macromolecules X-ray crystallography refinement protocols have had limited support for rigorous molecular physics in the past. For example, electrostatics is often neglected in favor of nonbonded interactions based on a purely repulsive van der Waals potential. In this work we present advanced algorithms for desktop workstations that open the door to X-ray refinement of even the most challenging macromolecular data sets using state-of-the-art classical molecular physics. First we describe theory for particle mesh Ewald (PME) summation that consistently handles the symmetry of all 230 space groups, replicates of the unit cell such that the minimum image convention can be used with a real space cutoff of any size and the combination of space group symmetry with replicates. An implementation of symmetry accelerated PME for the polarizable atomic multipole optimized energetics for biomolecular applications (AMOEBA) force field is presented. Relative to a single CPU core performing calculations on a P1 unit cell, our AMOEBA engine called Force Field X (FFX) accelerates energy evaluations by more than a factor of 24 on an 8-core workstation with a Tesla GPU coprocessor for 30 structures that contain 240 000 atoms on average in the unit cell. The benefit of AMOEBA electrostatics evaluated with PME for macromolecular X-ray crystallography refinement is demonstrated via rerefinement of 10 crystallographic data sets that range in resolution from 1.7 to 4.5 Å. Beginning from structures obtained by local optimization without electrostatics, further optimization using AMOEBA with PME electrostatics improved agreement of the model with the data (Rfree was lowered by 0.5%), improved geometric features such as favorable (ϕ, ψ) backbone conformations, and lowered the average potential energy per residue by over 10 kcal/mol. Furthermore, the MolProbity structure validation tool indicates that the geometry of these rerefined structures is consistent with X-ray crystallographic data collected up to 2.2 Å, which is 0.9 Å better than the actual mean quality (3.1 Å). We conclude that polarizable AMOEBA-assisted X-ray refinement offers advantages to methods that neglect electrostatics and is now efficient enough for routine use.