Convergence of the Electrostatic Interaction Based on Topological Atoms

Revealing the Critical Role of Global Electron Density Transfer in the Reaction Rate of Polar Organic Reactions within Molecular Electron Density Theory

The critical role of global electron density transfer (GEDT) in increasing the reaction rate of polar organic reactions has been studied within the framework of Molecular Electron Density Theory (MEDT). To this end, the series of the polar Diels-Alder (P-DA) reactions of cyclopentadiene with cyanoethylene derivatives, for which experimental kinetic data are available, have been chosen. A complete linear correlation between the computed activation Gibbs free energies and the GEDT taking place at the polar transition state structures (TSs) is found; the higher the GEDT at the TS, the lower the activation Gibbs free energy. An interacting quantum atoms energy partitioning analysis allows for establishing a complete linear correlation between the electronic stabilization of the electrophilic ethylene frameworks and the GEDT taking place at the polar TSs. This finding supports Parr's proposal for the definition of the electrophilicity ω index. The present MEDT study establishes the critical role of the GEDT in the acceleration of polar reactions, since the electronic stabilization of the electrophilic framework with the electron density gain is greater than the destabilization of the nucleophilic one, making a net favorable electronic contribution to the decrease in the activation energy.

Combining Force Fields and Neural Networks for an Accurate Representation of Chemically Diverse Molecular Interactions

A key goal of molecular modeling is the accurate reproduction of the true quantum mechanical potential energy of arbitrary molecular ensembles with a tractable classical approximation. The challenges are that analytical expressions found in general purpose force fields struggle to faithfully represent the intermolecular quantum potential energy surface at close distances and in strong interaction regimes; that the more accurate neural network approximations do not capture crucial physics concepts, e.g., nonadditive inductive contributions and application of electric fields; and that the ultra-accurate narrowly targeted models have difficulty generalizing to the entire chemical space. We therefore designed a hybrid wide-coverage intermolecular interaction model consisting of an analytically polarizable force field combined with a short-range neural network correction for the total intermolecular interaction energy. Here, we describe the methodology and apply the model to accurately determine the properties of water, the free energy of solvation of neutral and charged molecules, and the binding free energy of ligands to proteins. The correction is subtyped for distinct chemical species to match the underlying force field, to segment and reduce the amount of quantum training data, and to increase accuracy and computational speed. For the systems considered, the hybrid ab initio parametrized Hamiltonian reproduces the two-body dimer quantum mechanics (QM) energies to within 0.03 kcal/mol and the nonadditive many-molecule contributions to within 2%. Simulations of molecular systems using this interaction model run at speeds of several nanoseconds per day.

Electron density-based protocol to recover the interacting quantum atoms components of intermolecular binding energy

A new approach for obtaining interacting quantum atoms-defined components of binding energy of intermolecular interactions, which bypasses the use of standard six-dimensional integrals and two-particle reduced density matrix (2-RDM) reconstruction, is proposed. To examine this approach, three datasets calculated within the density functional theory framework using the def2-TZVP basis have been explored. The first two, containing 53 weakly bound bimolecular associates and 13 molecular clusters taken from the crystal, were used in protocol refinement, and the third one containing other 20 bimolecular and three cluster systems served as a validation reference. In addition, to verify the performance of the proposed approach on an exact 2-RDM, calculations within the coupled cluster formalism were performed for part of the first set systems using the cc-pVTZ basis set. The process of optimization of the proposed parametric model is considered, and the role of various energy contributions in the formation of non-covalent interactions is discussed with regard to the obtained trends.

Interacting Quantum Atoms and Multipolar Electrostatic Study of XH···π Interactions

The interaction energies of nine XH···π (X = C, N, and O) benzene-containing van der Waals complexes were analyzed, at the atomic and fragment levels, using QTAIM multipolar electrostatics and the energy partitioning method interacting quantum atoms/fragment (IQA/IQF). These descriptors were paired with the relative energy gradient method, which solidifies the connection between quantum mechanical properties and chemical interpretation. This combination provides a precise understanding, both qualitative and quantitative, of the nature of these interactions, which are ubiquitous in biochemical systems. The formation of the OH···π and NH···π systems is electrostatically driven, with the Qzz component of the quadrupole moment of the benzene carbons interacting with the charges of X and H in XH. There is the unexpectedly intramonomeric role of X-H (X = O, N) where its electrostatic energy helps the formation of the complex and its covalent energy thwarts it. However, the CH···π interaction is governed by exchange-correlation energies, thereby establishing a covalent character, as opposed to the literature's designation as a noncovalent interaction. Moreover, dispersion energy is relevant, statically and in absolute terms, but less relevant compared to other energy components in terms of the formation of the complex. Multipolar electrostatics are similar across all systems.

A Molecular Electron Density Theory Study of the Domino Reaction of N-Phenyl Iminoboranes with Benzaldehyde Yielding Fused Bicyclic Compounds

The reaction of N-phenyl iminoborane with benzaldehyde yielding a fused aromatic compound, recently reported by Liu et al., has been studied within the Molecular Electron Density Theory (MEDT). Formation of the fused aromatic compound is a domino process that comprises three consecutive reactions: (i) formation of a weak molecular complex between the reagents; (ii) an intramolecular electrophilic attack of the activated carbonyl carbon of benzaldehyde on the ortho position of the N-phenyl substituent of iminoborane; and (iii) a formal 1,3-hydrogen shift yielding the final fused aromatic compound. The two last steps correspond to a Friedel-Crafts acylation reaction, the product of the second reaction being the tetrahedral intermediate of an electrophilic aromatic substitution reaction. However, the presence of the imino group adjacent to the aromatic ring strongly stabilizes the corresponding intermediate, being the reaction product when the ortho positions are occupied by t-butyl substituents. This domino reaction shows a great similitude with the Brønsted acid catalyzed Povarov reaction. Although N-phenyl iminoborane can experience a formal [2+2] cycloaddition reaction with benzaldehyde, its higher activation Gibbs free energy compared to the intramolecular electrophilic attack of the activated carbonyl carbon of benzaldehyde on the ortho position of the N-phenyl substituent, 6.6 kcal·mol-1, prevents the formation of the formal [2+2] cycloadduct. The present MEDT study provides a different vision of the molecular mechanism of these reactions based on the electron density.

How to Compute Atomistic Insight in DFT Clusters: The REG-IQA Approach

The relative energy gradient (REG) method is paired with the topological energy partitioning method interacting quantum atoms (IQA), as REG-IQA, to provide detailed and unbiased knowledge on the intra- and interatomic interactions. REG operates on a sequence of geometries representing a dynamical change of a system. Its recent application to peptide hydrolysis of the human immunodeficiency virus-1 (HIV-1) protease (PDB code: 4HVP) has demonstrated its full potential in recovering reaction mechanisms and through-space electrostatic and exchange-correlation effects, making it a compelling tool for analyzing enzymatic reactions. In this study, the computational efficiency of the REG-IQA method for the 133-atom HIV-1 protease quantum mechanical system is analyzed in every detail and substantially improved by means of three different approaches. The first approach of smaller integration grids for IQA integrations reduces the computational overhead by about a factor of 3. The second approach uses the line-simplification Ramer-Douglas-Peucker (RDP) algorithm, which outputs the minimal number of geometries necessary for the REG-IQA analysis for a predetermined root mean squared error (RMSE) tolerance. This cuts the computational time of the whole REG analysis by a factor of 2 if an RMSE of 0.5 kJ/mol is considered. The third approach consists of a "biased" or "unbiased" selection of a specific subset of atoms of the whole initial quantum mechanical model wave-function, which results in more than a 10-fold speed-up per geometry for the IQA calculation, without deterioration of the outcome of the REG-IQA analysis. Finally, to show the capability of these approaches, the findings gathered from the HIV-1 protease system are also applied to a different system named haloalcohol dehalogenase (HheC). In summary, this study takes the REG-IQA method to a computationally feasible and highly accurate level, making it viable for the analysis of a multitude of enzymatic systems.

Carbonyl-Based Redox-Active Compounds as Organic Electrodes for Batteries: Escape from Middle-High Redox Potentials and Further Improvement?

Extracting─from the vast space of organic compounds─the best electrode candidates for achieving energy material breakthrough requires the identification of the microscopic causes and origins of various macroscopic features, including notably electrochemical and conduction properties. As a first guess of their capabilities, molecular DFT calculations and quantum theory of atoms in molecules (QTAIM)-derived indicators were applied to explore the family of pyrano[3,2-b]pyran-2,6-dione (PPD, i.e., A0) compounds, expanded to A0 fused with various kinds of rings (benzene, fluorinated benzene, thiophene, and merged thiophene/benzene). A glimpse of up-to-now elusive key incidences of introducing oxygen in vicinity to the carbonyl redox center within 6MRs─as embedded in the A0 core central unit common to all A-type compounds─has been gained. Furthermore, the main driving force toward achieving modulated low redox potential/band gaps thanks to fusing the aromatic rings for the A compound series was discovered.

Non-covalent interactions from a Quantum Chemical Topology perspective

About half a century after its little-known beginnings, the quantum topological approach called QTAIM has grown into a widespread, but still not mainstream, methodology of interpretational quantum chemistry. Although often confused in textbooks with yet another population analysis, be it perhaps an elegant but somewhat esoteric one, QTAIM has been enriched with about a dozen other research areas sharing its main mathematical language, such as Interacting Quantum Atoms (IQA) or Electron Localisation Function (ELF), to form an overarching approach called Quantum Chemical Topology (QCT). Instead of reviewing the latter’s role in understanding non-covalent interactions, we propose a number of ideas emerging from the full consequences of the space-filling nature of topological atoms, and discuss how they (will) impact on interatomic interactions, including non-covalent ones. The architecture of a force field called FFLUX, which is based on these ideas, is outlined. A new method called Relative Energy Gradient (REG) is put forward, which is able, by computation, to detect which fragments of a given molecular assembly govern the energetic behaviour of this whole assembly. This method can offer insight into the typical balance of competing atomic energies both in covalent and non-covalent case studies. A brief discussion on so-called bond critical points is given, highlighting concerns about their meaning, mainly in the arena of non-covalent interactions.

Using atomic charges to model molecular polarization

We review different models for introducing electric polarization in force fields, with special focus on methods where polarization is modelled at the atomic charge level. While electric polarization has been included in several force fields, the common approach has been to focus on atomic dipole polarizability. Several approaches allow modelling electric polarization by using charge-flow between charge sites instead, but this has been less exploited, despite that atomic charges and charge-flow is expected to be more important than atomic dipoles and dipole polarizability. A number of challenges are required to be solved for charge-flow models to be incorporated into polarizable force fields, for example how to parameterize the models and how to make them computational efficient.

Cocrystal Prediction Using Machine Learning Models and Descriptors

Cocrystals are of much interest in industrial application as well as academic research, and screening of suitable coformers for active pharmaceutical ingredients is the most crucial and challenging step in cocrystal development. Recently, machine learning techniques are attracting researchers in many fields including pharmaceutical research such as quantitative structure-activity/property relationship. In this paper, we develop machine learning models to predict cocrystal formation. We extract descriptor values from simplified molecular-input line-entry system (SMILES) of compounds and compare the machine learning models by experiments with our collected data of 1476 instances. As a result, we found that artificial neural network shows great potential as it has the best accuracy, sensitivity, and F1 score. We also found that the model achieved comparable performance with about half of the descriptors chosen by feature selection algorithms. We believe that this will contribute to faster and more accurate cocrystal development.

Interacting Quantum Atoms-A Review

The aim of this review is threefold. On the one hand, we intend it to serve as a gentle introduction to the Interacting Quantum Atoms (IQA) methodology for those unfamiliar with it. Second, we expect it to act as an up-to-date reference of recent developments related to IQA. Finally, we want it to highlight a non-exhaustive, yet representative set of showcase examples about how to use IQA to shed light in different chemical problems. To accomplish this, we start by providing a brief context to justify the development of IQA as a real space alternative to other existent energy partition schemes of the non-relativistic energy of molecules. We then introduce a self-contained algebraic derivation of the methodological IQA ecosystem as well as an overview of how these formulations vary with the level of theory employed to obtain the molecular wavefunction upon which the IQA procedure relies. Finally, we review the several applications of IQA as examined by different research groups worldwide to investigate a wide variety of chemical problems.

A FFLUX Water Model: Flexible, Polarizable and with a Multipolar Description of Electrostatics

Key to progress in molecular simulation is the development of advanced models that go beyond the limitations of traditional force fields that employ a fixed, point charge‐based description of electrostatics. Taking water as an example system, the FFLUX framework is shown capable of producing models that are flexible, polarizable and have a multipolar description of the electrostatics. The kriging machine‐learning methods used in FFLUX are able to reproduce the intramolecular potential energy surface and multipole moments of a single water molecule with chemical accuracy using as few as 50 training configurations. Molecular dynamics simulations of water clusters (25–216 molecules) using the new FFLUX model reveal that incorporating charge‐quadrupole, dipole–dipole, and quadrupole–charge interactions into the description of the electrostatics results in significant changes to the intermolecular structuring of the water molecules. © 2019 The Authors. Journal of Computational Chemistry published by Wiley Periodicals, Inc.

Ewald-based methods for Gaussian integral evaluation: application to a new parameterization of GEM

The development of accurate potentials for computational simulations has been an active area of research. Our group has been involved in the development of the Gaussian electrostatic model (GEM), a force field based on molecular densities. The philosophy of GEM is based on the pioneering work of N. Gresh and co-workers of the reproduction of individual inter-molecular interaction components obtained from quantum mechanical (QM) energy decomposition analysis (EDA). The molecular densities used in GEM are represented by fitting accurate QM molecular densities using auxiliary basis sets (comprised of Hermite Gaussians). The use of these molecular densities results in the need to evaluate a large number of Gaussian integrals. We have previously shown that the particle-mesh Ewald (PME), and fast Fourier Poisson (FFP) methods can be used for efficiently evaluating these types of integrals. Here, we present the latest parameterization of GEM* and its application for an extensive study of PME and FFP for molecular dynamics (MD) simulations using a hybrid version of our potential, GEM*. The temperature dependence of various bulk properties is presented and discussed, as well as the effect of various parameters affecting the performance/accuracy of both methods.

The influence of water potential in simulation: a catabolite activator protein case study

We present a rare comparison of structures of the same protein but generated by different potentials. We used four popular water potentials (SPC, TIP3P, TIP4P, TIP5P) in conjunction with the equally popular ff99SB. However, the ff12SB protein potential was used with TI3P only. Simulations (60 ns) were run on the catabolite activator protein (CAP), which is a textbook case of allosteric interaction. Overall, all potentials generated largely similar structures but failed to reproduce a crucial structural feature determined by NMR experiment. This example shows the need to develop next-generation potentials. Graphical abstract Catabolite activator protein.

Describing Molecular Polarizability by a Bond Capacity Model

We propose a bond capacity model for describing molecular polarization in force field energy functions at the charge-only level. Atomic charges are calculated by allowing charge to flow between atom pairs according to a bond capacity and a difference in electrostatic potential. The bond capacity is closely related to the bond order and decays to zero as the bond distance is increased. The electrostatic potential is composed of an intrinsic potential, identified as the electronegativity, and a screened Coulomb potential from all other charges. The bond capacity model leads to integer fragment charges upon bond dissociation and displays linear scaling of the polarizability with system size. Bond capacity parameters can be derived from reference molecular polarizabilities, while electronegativity parameters can be derived from reference atomic charges or a reference molecular electrostatic potential. Out-of-plane polarization for planar systems is modeled by off-nuclei charge sites. The model is shown to be able to reproduce anisotropic reference molecular polarizabilities with an accuracy of ∼10% using a limited set of bond capacity parameters and can describe both inter- and intramolecular polarization.

Tailoring van der Waals dispersion interactions with external electric charges

van der Waals (vdW) dispersion interactions strongly impact the properties of molecules and materials. Often, the description of vdW interactions should account for the coupling with pervasive electric fields, stemming from membranes, ionic channels, liquids, or nearby charged functional groups. However, this quantum-mechanical effect has been omitted in atomistic simulations, even in widely employed electronic-structure methods. Here, we develop a model and study the effects of an external charge on long-range vdW correlations. We show that a positive external charge stabilizes dispersion interactions, whereas a negative charge has an opposite effect. Our analytical results are benchmarked on a series of (bio)molecular dimers and supported by calculations with high-level correlated quantum-chemical methods, which estimate the induced dispersion to reach up to 35% of intermolecular binding energy (4 kT for amino-acid dimers at room temperature). Our analysis bridges electrostatic and electrodynamic descriptions of intermolecular interactions and may have implications for non-covalent reactions, exfoliation, dissolution, and permeation through biological membranes.

The long-range convergence of the energetic properties of the water monomer in bulk water at room temperature

The Interacting Quantum Atoms (IQA) energy partitioning scheme has been applied to a set of liquid water largely spherical clusters (henceforth called spheres) of up to 9 Å radius, with a maximum cluster size of 113 molecules. This constitutes half of the commonly used 216 molecules in a typical simulation box of a liquid water box, and to our knowledge is the largest analysis of this kind ever undertaken. As well as demonstrating the topological analysis of large systems, which has only recently become computationally feasible, important long range properties of liquid water are obtained. The full topological partitioning of each sphere into atomic basins is used to consider the long-range convergence of the energetic and multipolar properties of the water molecule at the centre of each sphere. It is found that the total molecular energy converges to its 9 Å value after 7 Å, which corresponds to approximately the first three solvation shells, while the molecular dipole and quadrupole moments approximately converge after 5.5 Å, which corresponds to approximately the first two solvation shells. The effect of water molecule flexibility is also considered.

A multipolar approach to the interatomic covalent interaction energy

Interatomic exchange-correlation energies correspond to the covalent energetic contributions to an interatomic interaction in real space theories of the chemical bond, but their widespread use is severely limited due to their computationally intensive character. In the same way as the multipolar (mp) expansion is customary used in biomolecular modeling to approximate the classical Coulomb interaction between two charge densities ρA(r) and ρB(r), we examine in this work the mp approach to approximate the interatomic exchange-correlation (xc) energies of the Interacting Quantum Atoms method. We show that the full xc mp series is quickly divergent for directly bonded atoms (1-2 pairs) albeit it works reasonably well most times for 1- n (n > 2) interactions. As with conventional perturbation theory, we show numerically that the xc series is asymptotically convergent and that, a truncated xc mp approximation retaining terms up to l1+l2=2 usually gives relatively accurate results, sometimes even for directly bonded atoms. Our findings are supported by extensive numerical analyses on a variety of systems that range from several standard hydrogen bonded dimers to typically covalent or aromatic molecules. The exact algebraic relationship between the monopole-monopole xc mp term and the inter-atomic bond order, as measured by the delocalization index of the quantum theory of atoms in molecules, is also established. © 2017 Wiley Periodicals, Inc.

Toward amino acid typing for proteins in FFLUX

Continuing the development of the FFLUX, a multipolar polarizable force field driven by machine learning, we present a modern approach to atom-typing and building transferable models for predicting atomic properties in proteins. Amino acid atomic charges in a peptide chain respond to the substitution of a neighboring residue and this response can be categorized in a manner similar to atom-typing. Using a machine learning method called kriging, we are able to build predictive models for an atom that is defined, not only by its local environment, but also by its neighboring residues, for a minimal additional computational cost. We found that prediction errors were up to 11 times lower when using a model specific to the correct group of neighboring residues, with a mean prediction of ∼0.0015 au. This finding suggests that atoms in a force field should be defined by more than just their immediate atomic neighbors. When comparing an atom in a single alanine to an analogous atom in a deca-alanine helix, the mean difference in charge is 0.026 au. Meanwhile, the same difference between a trialanine and a deca-alanine helix is only 0.012 au. When compared to deca-alanine models, the transferable models are up to 20 times faster to train, and require significantly less ab initio calculation, providing a practical route to modeling large biological systems. © 2016 The Authors. Journal of Computational Chemistry Published by Wiley Periodicals, Inc.

Dimerization of quercetin, Diels-Alder vs. radical-coupling approach: a joint thermodynamics, kinetics, and topological study

Quercetin is a prototypical antioxidant and prominent member of flavonoids, a large group of natural polyphenols. The oxidation of quercetin may lead to its dimerization, which is a paradigm of the more general polyphenol oligomerization. There exist two opposing mechanisms to describe the dimerization process, namely radical-coupling or Diels-Alder reactions. This work presents a comprehensive rationalization of this dimerization process, acquired from density functional theory (DFT) calculations. It is found that the two-step radical-coupling pathway is thermodynamically and kinetically preferred over the Diels-Alder reaction. This is in agreement with the experimental results showing the formation of only one isomer, whereas the Diels-Alder mechanism would yield two isomers. The evolution in bonding, occurring during these two processes, is investigated using the atoms in molecules (AIM) and electron localization function (ELF) topological approaches. It is shown that some electron density is accumulated between the fragments in the transition state of the radical-coupling reaction, but not in the transition state of the Diels-Alder process. Graphical Abstract Quantum chemistry calculations of the dimerization process of quercetin show that a radical coupling approach is preferred to a Diels-Alder type reaction, in agreement with experimental results. Analysis of the bonding evolution highlights the reaction mechanism.

Prediction of Intramolecular Polarization of Aromatic Amino Acids Using Kriging Machine Learning

Present computing power enables novel ways of modeling polarization. Here we show that the machine learning method kriging accurately captures the way the electron density of a topological atom responds to a change in the positions of the surrounding atoms. The success of this method is demonstrated on the four aromatic amino acids histidine, phenylalanine, tryptophan, and tyrosine. A new technique of varying training set sizes to vastly reduce training times while maintaining accuracy is described and applied to each amino acid. Each amino acid has its geometry distorted via normal modes of vibration over all local energy minima in the Ramachandran map. These geometries are then used to train the kriging models. Total electrostatic energies predicted by the kriging models for previously unseen geometries are compared to the true energies, yielding mean absolute errors of 2.9, 5.1, 4.2, and 2.8 kJ mol(-1) for histidine, phenylalanine, tryptophan, and tyrosine, respectively.

Prediction of conformationally dependent atomic multipole moments in carbohydrates

The conformational flexibility of carbohydrates is challenging within the field of computational chemistry. This flexibility causes the electron density to change, which leads to fluctuating atomic multipole moments. Quantum Chemical Topology (QCT) allows for the partitioning of an "atom in a molecule," thus localizing electron density to finite atomic domains, which permits the unambiguous evaluation of atomic multipole moments. By selecting an ensemble of physically realistic conformers of a chemical system, one evaluates the various multipole moments at defined points in configuration space. The subsequent implementation of the machine learning method kriging delivers the evaluation of an analytical function, which smoothly interpolates between these points. This allows for the prediction of atomic multipole moments at new points in conformational space, not trained for but within prediction range. In this work, we demonstrate that the carbohydrates erythrose and threose are amenable to the above methodology. We investigate how kriging models respond when the training ensemble incorporating multiple energy minima and their environment in conformational space. Additionally, we evaluate the gains in predictive capacity of our models as the size of the training ensemble increases. We believe this approach to be entirely novel within the field of carbohydrates. For a modest training set size of 600, more than 90% of the external test configurations have an error in the total (predicted) electrostatic energy (relative to ab initio) of maximum 1 kJ mol(-1) for open chains and just over 90% an error of maximum 4 kJ mol(-1) for rings.

Realistic sampling of amino acid geometries for a multipolar polarizable force field

The Quantum Chemical Topological Force Field (QCTFF) uses the machine learning method kriging to map atomic multipole moments to the coordinates of all atoms in the molecular system. It is important that kriging operates on relevant and realistic training sets of molecular geometries. Therefore, we sampled single amino acid geometries directly from protein crystal structures stored in the Protein Databank (PDB). This sampling enhances the conformational realism (in terms of dihedral angles) of the training geometries. However, these geometries can be fraught with inaccurate bond lengths and valence angles due to artefacts of the refinement process of the X-ray diffraction patterns, combined with experimentally invisible hydrogen atoms. This is why we developed a hybrid PDB/nonstationary normal modes (NM) sampling approach called PDB/NM. This method is superior over standard NM sampling, which captures only geometries optimized from the stationary points of single amino acids in the gas phase. Indeed, PDB/NM combines the sampling of relevant dihedral angles with chemically correct local geometries. Geometries sampled using PDB/NM were used to build kriging models for alanine and lysine, and their prediction accuracy was compared to models built from geometries sampled from three other sampling approaches. Bond length variation, as opposed to variation in dihedral angles, puts pressure on prediction accuracy, potentially lowering it. Hence, the larger coverage of dihedral angles of the PDB/NM method does not deteriorate the predictive accuracy of kriging models, compared to the NM sampling around local energetic minima used so far in the development of QCTFF.

Accurate prediction of polarised high order electrostatic interactions for hydrogen bonded complexes using the machine learning method kriging

As intermolecular interactions such as the hydrogen bond are electrostatic in origin, rigorous treatment of this term within force field methodologies should be mandatory. We present a method able of accurately reproducing such interactions for seven van der Waals complexes. It uses atomic multipole moments up to hexadecupole moment mapped to the positions of the nuclear coordinates by the machine learning method kriging. Models were built at three levels of theory: HF/6-31G(**), B3LYP/aug-cc-pVDZ and M06-2X/aug-cc-pVDZ. The quality of the kriging models was measured by their ability to predict the electrostatic interaction energy between atoms in external test examples for which the true energies are known. At all levels of theory, >90% of test cases for small van der Waals complexes were predicted within 1 kJ mol(-1), decreasing to 60-70% of test cases for larger base pair complexes. Models built on moments obtained at B3LYP and M06-2X level generally outperformed those at HF level. For all systems the individual interactions were predicted with a mean unsigned error of less than 1 kJ mol(-1).

The hydration of serine: multipole moments versus point charges

Next-generation force fields must incorporate improved electrostatic potentials in order to increase the reliability of their predictions. A crucial decision toward this goal is to abandon point charges in favour of multipole moments centered on nuclear sites. Here we compare the geometries generated by quantum topological multipole moments with those generated by four popular point charge models (TAFF, OPLS-AA, MMFF94x and PFROSST) for a hydrated serine. A main feature of this study is the dual comparison made, both at static level (geometry optimisation via energy minimisation) and at dynamic level (via molecular dynamics and radial/spatial distribution function analysis). At static level, multipolar electrostatics best reproduces the ab initio reference geometry. At dynamic level, multipolar electrostatics produces more structure than point charge electrostatics does, over the whole range. From our previous work on liquid water [Int. J. Quantum. Chem., 2004, 99, 685], where agreement with experiment only occurs when using multipole moments, we deduce that our predictions for hydrated serine will also be closer to experiment when using multipolar electrostatics. The spatial distribution function shows that only multipolar electrostatics shows pronounced structure at long range. Even at short range there are many regions where waters appear in the system governed by multipolar electrostatics but not in that governed by point charges.

Comprehensive analysis of energy minima of the 20 natural amino acids

Energy minima of the 20 natural amino acids (capped by a peptide bond at both the N and the C termini, CH3-C(═O)-N(H)-(H)Cα(R)-C(═O)-N(H)-CH3), were obtained by ab initio geometry optimization. Starting with a large number of minima, quickly generated by MarvinView, geometry optimization at the HF/6-31G(d,p) level of theory reduced the number of minima, followed by further optimization at the B3LYP/apc-1 and MP2/cc-pVDZ levels, which caused some minima to disappear and some stable minima to migrate on the Ramachandran map. There is a relation between the number of minima and the size and the flexibility of the side chain. The energy minima of the 20 amino acids are mainly located in the regions of βL, γL, δL, and αL of the Ramachandran map. Multipole moments of atoms occurring in the fragment [-NH-Cα-C(═O)-] common to all 20 amino acids were calculated at the three levels of theory mentioned above. The near parallelism in behavior of these moments between levels of theory is beneficial toward estimating moments with the more expensive B3LYP and MP2 methods from data calculated with the cheaper HF method. Finally, we explored the transferability of properties between different amino acids: the bond length and angles of the common fragment [-NH-Cα(HαCβ)-C'(═O)-] in all amino acids except Gly and Pro. All bond lengths are highly transferable between different amino acids, and the standard deviations are small.

Revisiting H2O Nucleation around Au+ and Hg2+: The Peculiar "Pseudo-Soft" Character of the Gold Cation

In this contribution, we propose a deeper understanding of the electronic effects affecting the nucleation of water around the Au+ and Hg2+ metal cations using quantum chemistry. To do so, and in order to go beyond usual energetical studies, we make extensive use of state of the art quantum interpretative techniques combining ELF/NCI/QTAIM/EDA computations to capture all ranges of interactions stabilizing the well characterized microhydrated structures. The Electron Localization Function (ELF) topological analysis reveals the peculiar role of the Au+ outer-shell core electrons (subvalence) that appear already spatially preorganized once the addition of the first water molecule occurs. Thus, despite the addition of other water molecules, the electronic structure of Au(H2O)+ appears frozen due to relativistic effects leading to a maximal acceptation of only two waters in gold's first hydration shell. As the values of the QTAIM (Quantum Theory of Atoms in Molecules) cations's charge is discussed, the Non Covalent Interactions (NCI) analysis showed that Au+ appears still able to interact through longer range van der Waals interaction with the third or fourth hydration shell water molecules. As these types of interaction are not characteristic of either a hard or soft metal cation, we introduced the concept of a "pseudo-soft" cation to define Au+ behavior. Then, extending the study, we performed the same computations replacing Au+ with Hg2+, an isoelectronic cation. If Hg2+ behaves like Au+ for small water clusters, a topological, geometrical, and energetical transition appears when the number of water molecules increases. Regarding the HSAB theory, this transition is characteristic of a shift of Hg2+ from a pseudosoft form to a soft ion and appears to be due to a competition between the relativistic and correlation effects. Indeed, if relativistic effects are predominant, then mercury will behave like gold and have a similar subvalence/geometry; otherwise when correlation effects are predominant, Hg2+ behaves like a soft cation.

GEM*: A Molecular Electronic Density-Based Force Field for Molecular Dynamics Simulations

GEM*, a force field that combines Coulomb and Exchange terms calculated with Hermite Gaussians with the polarization, bonded, and modified van der Waals terms from AMOEBA is presented. GEM* is tested on an initial water model fitted at the same level as AMOEBA. The integrals required for the evaluation of the intermolecular Coulomb interactions are efficiently evaluated by means of reciprocal space methods. The GEM* water model is tested by comparing energies and forces for a series of water oligomers and MD simulations. Timings for GEM* compared to AMOEBA are presented and discussed.

Accuracy and tractability of a kriging model of intramolecular polarizable multipolar electrostatics and its application to histidine

We propose a generic method to model polarization in the context of high-rank multipolar electrostatics. This method involves the machine learning technique kriging, here used to capture the response of an atomic multipole moment of a given atom to a change in the positions of the atoms surrounding this atom. The atoms are malleable boxes with sharp boundaries, they do not overlap and exhaust space. The method is applied to histidine where it is able to predict atomic multipole moments (up to hexadecapole) for unseen configurations, after training on 600 geometries distorted using normal modes of each of its 24 local energy minima at B3LYP/apc-1 level. The quality of the predictions is assessed by calculating the Coulomb energy between an atom for which the moments have been predicted and the surrounding atoms (having exact moments). Only interactions between atoms separated by three or more bonds ("1, 4 and higher" interactions) are included in this energy error. This energy is compared with that of a central atom with exact multipole moments interacting with the same environment. The resulting energy discrepancies are summed for 328 atom-atom interactions, for each of the 29 atoms of histidine being a central atom in turn. For 80% of the 539 test configurations (outside the training set), this summed energy deviates by less than 1 kcal mol(-1).

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.

Biomolecular electrostatics and solvation: a computational perspective

An understanding of molecular interactions is essential for insight into biological systems at the molecular scale. Among the various components of molecular interactions, electrostatics are of special importance because of their long-range nature and their influence on polar or charged molecules, including water, aqueous ions, proteins, nucleic acids, carbohydrates, and membrane lipids. In particular, robust models of electrostatic interactions are essential for understanding the solvation properties of biomolecules and the effects of solvation upon biomolecular folding, binding, enzyme catalysis, and dynamics. Electrostatics, therefore, are of central importance to understanding biomolecular structure and modeling interactions within and among biological molecules. This review discusses the solvation of biomolecules with a computational biophysics view toward describing the phenomenon. While our main focus lies on the computational aspect of the models, we provide an overview of the basic elements of biomolecular solvation (e.g. solvent structure, polarization, ion binding, and non-polar behavior) in order to provide a background to understand the different types of solvation models.