Ab initio molecular dynamics with nuclear quantum effects at classical cost: Ring polymer contraction for density functional theory

Water-air interface revisited by means of path-integral ab initio molecular dynamics

Although nuclear quantum effects (NQEs) have been considered on bulk liquid water, the impact of these latter on the air-water interface has not yet been reported. Herein, by performing and comparing ab initio molecular dynamics (AIMD) and path integral AIMD (PI-AIMD) simulations, we reveal the impact of NQEs on structural, dynamical and electronic properties as well as IR spectra of the air-water interface at room temperature. NQEs, being able to describe a more accurate proton delocalization in H-bonded system than AIMD, reveal a different structural arrangement and dynamical behaviour of both bulk and interfacial water molecules in comparison to AIMD results. A more de-structured and de-bound water arrangement and coordination are identified when the quantum nature of nuclei are considered for both bulk and interfacial water molecules. Structural properties, such as inter-/intra-molecular bond lengths, coordination numbers and H-bonding angles of bulk and interfacial water molecules here calculated, are affected by NQEs mitigating the overstructured description given by AIMD. Further evidences of an AIMD overstructured description of bulk water are in the computed IR spectra, where an increased absorption peak intensity and an increased strength of the hydrogen-bond network are alleviated by NQEs. In addition, NQEs show a valuable impact on the electronic structure of the air-water interface, reducing the total valence bandwidth and the electronic energy band-gap when passing from bulk to interfacial water. This work proves how NQEs significantly affect properties and features of the air-water interface, that are essential to accurately describe H-bonded interfacial systems.

Path Integral Simulations of Condensed-Phase Vibrational Spectroscopy

Recent theoretical and algorithmic developments have improved the accuracy with which path integral dynamics methods can include nuclear quantum effects in simulations of condensed-phase vibrational spectra. Such methods are now understood to be approximations to the delocalized classical Matsubara dynamics of smooth Feynman paths, which dominate the dynamics of systems such as liquid water at room temperature. Focusing mainly on simulations of liquid water and hexagonal ice, we explain how the recently developed quasicentroid molecular dynamics (QCMD), fast-QCMD, and temperature-elevated path integral coarse-graining simulations (Te PIGS) methods generate classical dynamics on potentials of mean force obtained by averaging over quantum thermal fluctuations. These new methods give very close agreement with one another, and the Te PIGS method has recently yielded excellent agreement with experimentally measured vibrational spectra for liquid water, ice, and the liquid-air interface. We also discuss the limitations of such methods.

Multiple-time scale integration method based on an interpolated potential energy surface for ab initio path integral molecular dynamics

A new multiple-time scale integration method is presented that propagates ab initio path integral molecular dynamics (PIMD). This method uses a large time step to generate an approximate geometrical configuration whose energy and gradient are evaluated at the level of an ab initio method, and then, a more precise integration scheme, e.g., the Bulirsch-Stoer method or velocity Verlet integration with a smaller time step, is used to integrate from the previous step using the computationally efficient interpolated potential energy surface constructed from two consecutive points. This method makes the integration of PIMD more efficient and accurate compared with the velocity Verlet integration. A Nosé-Hoover chain thermostat combined with this new multiple-time scale method has good energy conservation even with a large time step, which is usually challenging in velocity Verlet integration for PIMD due to the very small chain mass when a large number of beads are used. The new method is used to calculate infrared spectra and free energy profiles to demonstrate its accuracy and capabilities.

Intra-cluster Charge Migration upon Hydration of Protonated Formic Acid Revealed by Anharmonic Analysis of Cold Ion Vibrational Spectra

The rates of many chemical reactions are accelerated when carried out in micron-sized droplets, but the molecular origin of the rate acceleration remains unclear. One example is the condensation reaction of 1,2-diaminobenzene with formic acid to yield benzimidazole. The observed rate enhancements have been rationalized by invoking enhanced acidity at the surface of methanol solvent droplets with low water content to enable protonation of formic acid to generate a cationic species (protonated formic acid or PFA) formed by attachment of a proton to the neutral acid. Because PFA is the key feature in this reaction mechanism, vibrational spectra of cryogenically cooled, microhydrated PFA·(H2O)n=1-6 were acquired to determine how the extent of charge localization depends on the degree of hydration. Analysis of these highly anharmonic spectra with path integral ab initio molecular dynamics simulations reveals the gradual displacement of the excess proton onto the water network in the microhydration regime at low temperatures with n = 3 as the tipping point for intra-cluster proton transfer.

Accurate SCC-DFTB Parametrization of Liquid Water with Improved Atomic Charges and Iterative Boltzmann Inversion

This work presents improvements of the description of liquid water within the self-consistent-charge density-functional based tight-binding scheme combining the use of Weighted Mulliken (WMull) charges and optimized O-H repulsive potential through the iterative Boltzmann inversion (IBI) process. The quality of the newly developed models is validated considering pair radial distribution functions (RDFs), as well as other structural, energetic, thermodynamic, and dynamic properties. The use of WMull charges certainly improves the agreement with experimental data, however leading to over-structured RDFs at short distance, that can be further improved by considering an optimized O-H repulsive potential obtained by the IBI process. Three different schemes were used to optimize this potential: (i) optimization including short O-H distances. This led to accurate RDFs as well as improved self-diffusion coefficient and heat of vaporization, while the proton transfer energy barrier is severely deteriorated; (ii) optimization starting at long distance. The proton transfer energy barrier is recovered while the heat of vaporization is deteriorated and the O-H RDF is less accurate at short distance; (iii) optimization within the path-integral molecular dynamics scheme which allows us to exclude nuclear quantum effects from the repulsive potential. The latter potential, in conjunction with the WMull improved atomic charges, provides similar results as (i) for structural, dynamic, and thermodynamic properties while recovering a large part of the proton transfer energy barrier. It therefore offers a good compromise to study both dynamic properties and chemistry within liquid water at a quantum chemical level.

Developing machine-learned potentials to simultaneously capture the dynamics of excess protons and hydroxide ions in classical and path integral simulations

The transport of excess protons and hydroxide ions in water underlies numerous important chemical and biological processes. Accurately simulating the associated transport mechanisms ideally requires utilizing ab initio molecular dynamics simulations to model the bond breaking and formation involved in proton transfer and path-integral simulations to model the nuclear quantum effects relevant to light hydrogen atoms. These requirements result in a prohibitive computational cost, especially at the time and length scales needed to converge proton transport properties. Here, we present machine-learned potentials (MLPs) that can model both excess protons and hydroxide ions at the generalized gradient approximation and hybrid density functional theory levels of accuracy and use them to perform multiple nanoseconds of both classical and path-integral proton defect simulations at a fraction of the cost of the corresponding ab initio simulations. We show that the MLPs are able to reproduce ab initio trends and converge properties such as the diffusion coefficients of both excess protons and hydroxide ions. We use our multi-nanosecond simulations, which allow us to monitor large numbers of proton transfer events, to analyze the role of hypercoordination in the transport mechanism of the hydroxide ion and provide further evidence for the asymmetry in diffusion between excess protons and hydroxide ions.

Implementation of the self-consistent phonons method with ab initio potentials (AI-SCP)

The self-consistent phonon (SCP) method allows one to include anharmonic effects when treating a many-body quantum system at thermal equilibrium. The system is then described by an effective temperature-dependent harmonic Hamiltonian, which can be used to estimate its various dynamic and static properties. In this paper, we combine SCP with ab initio (AI) potential energy evaluation in which case the numerical bottleneck of AI-SCP is the evaluation of Gaussian averages of the AI potential energy and its derivatives. These averages are computed efficiently by the quasi-Monte Carlo method utilizing low-discrepancy sequences leading to a fast convergence with respect to the number, S, of the AI energy evaluations. Moreover, a further substantial (an-order-of-magnitude) improvement in efficiency is achieved once a numerically cheap approximation of the AI potential is available. This is based on using a perturbation theory-like (the two-grid) approach in which it is the average of the difference between the AI and the approximate potential that is computed. The corresponding codes and scripts are provided.

Routine Molecular Dynamics Simulations Including Nuclear Quantum Effects: From Force Fields to Machine Learning Potentials

We report the implementation of a multi-CPU and multi-GPU massively parallel platform dedicated to the explicit inclusion of nuclear quantum effects (NQEs) in the Tinker-HP molecular dynamics (MD) package. The platform, denoted Quantum-HP, exploits two simulation strategies: the Ring-Polymer Molecular Dynamics (RPMD) that provides exact structural properties at the cost of a MD simulation in an extended space of multiple replicas and the adaptive Quantum Thermal Bath (adQTB) that imposes the quantum distribution of energy on a classical system via a generalized Langevin thermostat and provides computationally affordable and accurate (though approximate) NQEs. We discuss some implementation details, efficient numerical schemes, and parallelization strategies and quickly review the GPU acceleration of our code. Our implementation allows an efficient inclusion of NQEs in MD simulations for very large systems, as demonstrated by scaling tests on water boxes with more than 200,000 atoms (simulated using the AMOEBA polarizable force field). We test the compatibility of the approach with Tinker-HP's recently introduced Deep-HP machine learning potentials module by computing water properties using the DeePMD potential with adQTB thermostatting. Finally, we show that the platform is also compatible with the alchemical free energy estimation capabilities of Tinker-HP and fast enough to perform simulations. Therefore, we study how NQEs affect the hydration free energy of small molecules solvated with the recently developed Q-AMOEBA water force field. Overall, the Quantum-HP platform allows users to perform routine quantum MD simulations of large condensed-phase systems and will help to shed new light on the quantum nature of important interactions in biological matter.

Incorporating Nuclear Quantum Effects in Molecular Dynamics with a Constrained Minimized Energy Surface

The accurate incorporation of nuclear quantum effects in large-scale molecular dynamics (MD) simulations remains a significant challenge. Recently, we combined constrained nuclear-electronic orbital (CNEO) theory with classical MD and obtained a new approach (CNEO-MD) that can accurately and efficiently incorporate nuclear quantum effects into classical simulations. In this Letter, we provide the theoretical foundation for CNEO-MD by developing an alternative formulation of the equations of motion for MD. In this new formulation, the expectation values of quantum nuclear positions evolve classically on an effective energy surface that is obtained from a constrained energy minimization procedure when solving for the quantum nuclear wave function, thus enabling the incorporation of nuclear quantum effects in classical MD simulations. For comparison with other existing approaches, we examined a series of model systems and found that this new MD approach is significantly more accurate than the conventional way of performing classical MD and generally outperforms centroid MD and ring-polymer MD in describing vibrations in these model systems.

2D spectroscopies from condensed phase dynamics: Accessing third-order response properties from equilibrium multi-time correlation functions

The third-order response lies at the heart of simulating and interpreting nonlinear spectroscopies ranging from two dimensional infrared (2D-IR) to 2D electronic (2D-ES), and 2D sum frequency generation (2D-SFG). The extra time and frequency dimensions in these spectroscopies provides access to rich information on the electronic and vibrational states present, the coupling between them, and the resulting rates at which they exchange energy that are obscured in linear spectroscopy, particularly for condensed phase systems that usually contain many overlapping features. While the exact quantum expression for the third-order response is well established it is incompatible with the methods that are practical for calculating the atomistic dynamics of large condensed phase systems. These methods, which include both classical mechanics and quantum dynamics methods that retain quantum statistical properties while obeying the symmetries of classical dynamics, such as LSC-IVR, Centroid Molecular Dynamics (CMD) and Ring Polymer Molecular Dynamics (RPMD) naturally provide short-time approximations to the multi-time symmetrized Kubo transformed correlation function. Here, we show how the third-order response can be formulated in terms of equilibrium symmetrized Kubo transformed correlation functions. We demonstrate the utility and accuracy of our approach by showing how it can be used to obtain the third-order response of a series of model systems using both classical dynamics and RPMD. In particular, we show that this approach captures features such as anharmonically induced vertical splittings and peak shifts while providing a physically transparent framework for understanding multidimensional spectroscopies.

Structural and Dynamic Properties of Solvated Hydroxide and Hydronium Ions in Water from Ab Initio Modeling

Predicting the asymmetric structure and dynamics of solvated hydroxide and hydronium in water has been a challenging task from ab initio molecular dynamics (AIMD). The difficulty mainly comes from a lack of accurate and efficient exchange-correlation functional in elucidating the amphiphilic nature and the ubiquitous proton transfer behaviors of the two ions. By adopting the strongly-constrained and appropriately normed (SCAN) meta-GGA functional in AIMD simulations, we systematically examine the amphiphilic properties, the solvation structures, the electronic structures, and the dynamic properties of the two water ions. In particular, we compare these results to those predicted by the PBE0-TS functional, which is an accurate yet computationally more expensive exchange-correlation functional. We demonstrate that the general-purpose SCAN functional provides a reliable choice in describing the two water ions. Specifically, in the SCAN picture of water ions, the appearance of the fourth and fifth hydrogen bonds near hydroxide stabilizes the pot-like shape solvation structure and suppresses the structural diffusion, while the hydronium stably donates three hydrogen bonds to its neighbors. We apply a detailed analysis of the proton transfer mechanism of the two ions and find the two ions exhibit substantially different proton transfer patterns. Our AIMD simulations indicate hydroxide diffuses slower than hydronium in water, which is consistent with the experiments.

Using Machine Learning to Greatly Accelerate Path Integral Ab Initio Molecular Dynamics

Ab initio molecular dynamics (AIMD) has become one of the most popular and robust approaches for modeling complicated chemical, liquid, and material systems. However, the formidable computational cost often limits its widespread application in simulations of the largest-scale systems. The situation becomes even more severe in cases where the hydrogen nuclei may be better described as quantized particles using a path integral representation. Here, we present a computational approach that combines machine learning with recent advances in path integral contraction schemes, and we achieve a 2 orders of magnitude acceleration over direct path integral AIMD simulation while at the same time maintaining its accuracy.

Affordable Ab Initio Path Integral for Thermodynamic Properties via Molecular Dynamics Simulations Using Semiempirical Reference Potential

Path integral molecular dynamics (PIMD) is becoming a routinely applied method for incorporating the nuclear quantum effect in computer simulations. However, direct PIMD simulations at an ab initio level of theory are formidably expensive. Using the protonated 1,8-bis(dimethylamino)naphthalene molecule as an example, we show in this work that the computational expense for the intramolecular proton transfer between the two nitrogen atoms can be remarkably reduced by implementing the idea of reference-potential methods. The simulation time can be easily extended to a scale of nanoseconds while maintaining the accuracy on an ab initio level of theory for thermodynamic properties. In addition, postprocessing can be carried out in parallel on massive computer nodes. A 545-fold reduction in the total CPU time can be achieved in this way as compared to a direct PIMD simulation at the same ab initio level of theory.

Nuclear Quantum Effects in Liquid Water at Near Classical Computational Cost Using the Adaptive Quantum Thermal Bath

We demonstrate the accuracy and efficiency of a recently introduced approach to account for nuclear quantum effects (NQEs) in molecular simulations: the adaptive quantum thermal bath (adQTB). In this method, zero-point energy is introduced through a generalized Langevin thermostat designed to precisely enforce the quantum fluctuation-dissipation theorem. We propose a refined adQTB algorithm with improved accuracy and report adQTB simulations of liquid water. Through extensive comparison with reference path integral calculations, we demonstrate that it provides excellent accuracy for a broad range of structural and thermodynamic observables as well as infrared vibrational spectra. The adQTB has a computational cost comparable to that of classical molecular dynamics, enabling simulations of up to millions of degrees of freedom.

Direct observation of ultrafast hydrogen bond strengthening in liquid water

Water is one of the most important, yet least understood, liquids in nature. Many anomalous properties of liquid water originate from its well-connected hydrogen bond network¹, including unusually efficient vibrational energy redistribution and relaxation². An accurate description of the ultrafast vibrational motion of water molecules is essential for understanding the nature of hydrogen bonds and many solution-phase chemical reactions. Most existing knowledge of vibrational relaxation in water is built upon ultrafast spectroscopy experiments^2-7. However, these experiments cannot directly resolve the motion of the atomic positions and require difficult translation of spectral dynamics into hydrogen bond dynamics. Here, we measure the ultrafast structural response to the excitation of the OH stretching vibration in liquid water with femtosecond temporal and atomic spatial resolution using liquid ultrafast electron scattering. We observed a transient hydrogen bond contraction of roughly 0.04 Å on a timescale of 80 femtoseconds, followed by a thermalization on a timescale of approximately 1 picosecond. Molecular dynamics simulations reveal the need to treat the distribution of the shared proton in the hydrogen bond quantum mechanically to capture the structural dynamics on femtosecond timescales. Our experiment and simulations unveil the intermolecular character of the water vibration preceding the relaxation of the OH stretch.

Investigating Tunneling-Controlled Chemical Reactions through Ab Initio Ring Polymer Molecular Dynamics

We use the ab initio ring polymer molecular dynamics (RPMD) approach to investigate tunneling-controlled reactions in methylhydroxycarbene. Nuclear tunneling effects enable molecules to overcome the barriers which cannot be overcome classically. Under low-temperature conditions, intrinsic quantum tunneling effects can facilitate the chemical reaction in a pathway that is favored neither thermodynamically nor kinetically. This behavior is referred to as the tunneling-controlled chemical reaction and is regarded as the third paradigm of chemical reaction controls. In this work, we use the ab initio RPMD approach to incorporate the tunneling effects in our quantum dynamics simulations and investigate the reaction kinetics of two competitive reaction pathways at various temperatures. The reaction rate constants obtained here agree extremely well with the experimentally measured rates. We demonstrate the feasibility of using ab initio RPMD rate calculations in a realistic molecular system and provide an interesting and important example for future investigations of reaction mechanisms dominated by quantum tunneling effects.

Nuclear Quantum Effect and Its Temperature Dependence in Liquid Water from Random Phase Approximation via Artificial Neural Network

We report structural and dynamical properties of liquid water described by the random phase approximation (RPA) correlation together with the exact exchange energy (EXX) within density functional theory. By utilizing thermostated ring polymer molecular dynamics, we examine the nuclear quantum effects and their temperature dependence. We circumvent the computational limitation of performing direct first-principles molecular dynamics simulation at this high level of electronic structure theory by adapting an artificial neural network model. We show that the EXX+RPA level of theory accurately describes liquid water in terms of both dynamical and structural properties.

Dynamics & Spectroscopy with Neutrons-Recent Developments & Emerging Opportunities

This work provides an up-to-date overview of recent developments in neutron spectroscopic techniques and associated computational tools to interrogate the structural properties and dynamical behavior of complex and disordered materials, with a focus on those of a soft and polymeric nature. These have and continue to pave the way for new scientific opportunities simply thought unthinkable not so long ago, and have particularly benefited from advances in high-resolution, broadband techniques spanning energy transfers from the meV to the eV. Topical areas include the identification and robust assignment of low-energy modes underpinning functionality in soft solids and supramolecular frameworks, or the quantification in the laboratory of hitherto unexplored nuclear quantum effects dictating thermodynamic properties. In addition to novel classes of materials, we also discuss recent discoveries around water and its phase diagram, which continue to surprise us. All throughout, emphasis is placed on linking these ongoing and exciting experimental and computational developments to specific scientific questions in the context of the discovery of new materials for sustainable technologies.

Pure non-local machine-learned density functional theory for electron correlation

Density-functional theory (DFT) is a rigorous and (in principle) exact framework for the description of the ground state properties of atoms, molecules and solids based on their electron density. While computationally efficient density-functional approximations (DFAs) have become essential tools in computational chemistry, their (semi-)local treatment of electron correlation has a number of well-known pathologies, e.g. related to electron self-interaction. Here, we present a type of machine-learning (ML) based DFA (termed Kernel Density Functional Approximation, KDFA) that is pure, non-local and transferable, and can be efficiently trained with fully quantitative reference methods. The functionals retain the mean-field computational cost of common DFAs and are shown to be applicable to non-covalent, ionic and covalent interactions, as well as across different system sizes. We demonstrate their remarkable possibilities by computing the free energy surface for the protonated water dimer at hitherto unfeasible gold-standard coupled cluster quality on a single commodity workstation.

Semilocal density functionals, while computationally efficient, do not account for non-local correlation. Here, the authors propose a machine-learning approach to DFT that leads to non-local and transferable functionals applicable to non-covalent, ionic and covalent interactions across system of different sizes.

Diabat method for polymorph free energies: Extension to molecular crystals

Lattice-switch Monte Carlo and the related diabat methods have emerged as efficient and accurate ways to compute free energy differences between polymorphs. In this work, we introduce a one-to-one mapping from the reference positions and displacements in one molecular crystal to the positions and displacements in another. Two features of the mapping facilitate lattice-switch Monte Carlo and related diabat methods for computing polymorph free energy differences. First, the mapping is unitary so that its Jacobian does not complicate the free energy calculations. Second, the mapping is easily implemented for molecular crystals of arbitrary complexity. We demonstrate the mapping by computing free energy differences between polymorphs of benzene and carbamazepine. Free energy calculations for thermodynamic cycles, each involving three independently computed polymorph free energy differences, all return to the starting free energy with a high degree of precision. The calculations thus provide a force field independent validation of the method and allow us to estimate the precision of the individual free energy differences.

Nuclear Quantum Effects Largely Influence Molecular Dissociation and Proton Transfer in Liquid Water under an Electric Field

Proton transfer in liquid water controls acid-base chemistry, crucial enzyme reactions, and the functioning of fuel cells. Externally applied static electric fields in water are capable of dissociating molecules and transferring protons across the H-bond network. However, the impact of nuclear quantum effects (NQEs) on these fundamental field-induced phenomena has not yet been reported. By comparing state-of-the-art ab initio molecular dynamics (AIMD) and path integral AIMD simulations of water under electric fields, I show that quantum delocalization of the proton lowers the molecular ionization threshold to approximately one-third. Moreover, also the water behavior as a protonic semiconductor is considerably modified by the inclusion of NQEs. In fact, when the quantum nature of the nuclei is taken into account, the proton conductivity is ∼50% larger. This work proves that NQEs sizably affect the protolysis phenomenon and proton transfer in room-temperature liquid water.

Simulation of Liquids with the Tight-Binding Density-Functional Approach and Improved Atomic Charges

Theoretical description of liquids, especially liquid water, is an ongoing subject with important implications in various domains such as homogeneous catalysis; solvation of molecular, ionic, and biomolecular species; and reactivity. Various formalisms exist to describe liquids, each one displaying its own balance between accuracy and computational cost that defines its range of applications. The present article revisits the ability of the density-functional-based tight-binding (SCC-DFTB) approach to model liquids by focusing on liquid water and liquid benzene under ambient conditions. To do so, we benchmark a recent correction for the SCC-DFTB atomic charges that allows for a drastic improvement of the pair radial distribution functions of liquid water as compared to both experimental data and density-functional theory results performed in the generalized-gradient approximation. We also report the coupling of the deMonNano and i-PI codes to perform path-integral molecular dynamics. This allows us to rationalize the impact of nuclear quantum effects on the SCC-DFTB description of liquid water. This study evidences the rather good ability of SCC-DFTB to describe liquid water and liquid benzene. As the first example of application, we also present results for a benzene molecule solvated in water with the perspectives of further studies devoted to solvent/water interfaces.

Elucidating the Proton Transport Pathways in Liquid Imidazole with First-Principles Molecular Dynamics

Imidazole is a promising anhydrous proton conductor with a high conductivity comparable to that of water at a similar temperature relative to its melting point. Previous theoretical studies of the mechanism of proton transport in imidazole have relied either on empirical models or on ab initio trajectories that have been too short to draw significant conclusions. Here, we present the results of multiple time-step ab initio molecular dynamics simulations of an excess proton in liquid imidazole reaching 1 ns in total simulation time. We find that the proton transport is dominated by structural diffusion, with the diffusion constant of the proton defect being ∼8 times higher than that of self-diffusion of the imidazole molecules. By using correlation function analysis, we decompose the mechanism for proton transport into a series of first-order processes and show that the proton transport mechanism occurs over three distinct time and length scales. Although the mechanism at intermediate times is dominated by hopping along pseudo-one-dimensional chains, at longer times the overall rate of diffusion is limited by the re-formation of these chains. These results provide a more complete picture of the traditional idealized Grotthuss structural diffusion mechanism.

NMR crystallography of molecular organics

Developments of NMR methodology to characterise the structures of molecular organic structures are reviewed, concentrating on the previous decade of research in which density functional theory-based calculations of NMR parameters in periodic solids have become widespread. With a focus on demonstrating the new structural insights provided, it is shown how "NMR crystallography" has been used in a spectrum of applications from resolving ambiguities in diffraction-derived structures (such as hydrogen atom positioning) to deriving complete structures in the absence of diffraction data. As well as comprehensively reviewing applications, the different aspects of the experimental and computational techniques used in NMR crystallography are surveyed. NMR crystallography is seen to be a rapidly maturing subject area that is increasingly appreciated by the wider crystallographic community.

Accurate Description of Nuclear Quantum Effects with High-Order Perturbed Path Integrals (HOPPI)

Imaginary time path-integral (PI) simulations that account for nuclear quantum effects (NQE) beyond the harmonic approximation are increasingly employed together with modern electronic-structure calculations. Existing PI methods are applicable to molecules, liquids, and solids; however, the computational cost of such simulations increases dramatically with decreasing temperature. To address this challenge, here, we propose to combine high-order PI factorization with perturbation theory (PT). Already for conventional second-order PI simulations, the PT ansatz increases the accuracy 2-fold compared to fourth-order schemes with the same settings. In turn, applying PT to high-order path integrals (HOPI) further improves the efficiency of simulations for molecular and condensed matter systems especially at low temperatures. We present results for bulk liquid water, the aspirin molecule, and the CH₅⁺ molecule. Perturbed HOPI simulations remain both efficient and accurate down to 20 K and provide a convenient method to estimate the convergence of quantum-mechanical observables.

A new approach to estimate atomic energies

A new approach to estimate atomic energies is introduced. The method is based on utilization of experimental ionization energies as well as conversion of "n-electron atomic systems" to n "one-electron systems." Sample detail calculations are presented with typical graphs to show the distribution of different types of energy within an atom. The breakdown of atomic energies into kinetic, electron-nucleus attraction, and electron-electron repulsion is shown within an atom as well the total of energies of each type for elements. Then in a following step, the variations in kinetic, electron-electron, and electron-nucleus interaction energies of electrons as evidence for atomic shell changes are presented. Furthermore, the article overviews the spatial gaps between orbitals as an added evidence for existence of electronic shells. The findings in this article have significant implications for the structure of atoms and the layout of periodic table. Graphical abstractElectronic energy variations of barium (Ba) for kinetic, electron-electron repulsion, and electron-nucleus energies versus electron number.

Atomic shells according to ionization energies

This article relies only on experimental data rather than getting involved with theories, calculations, approximations, or interpolations. Experimental ionization energies of all atoms in the periodic table are collected and utilized to discover the order of electronic shells. The assumption in this paper is mainly the energy difference between atomic shells. In other words, one should observe an abrupt change in the energy moving from one atomic shell to another. Electronic energies within an atom are either kinetic or potential. Potential energy can be further broken into "electron-nucleus attraction" and "electron-electron repulsion" energies. The ionization energy that holds an electron onto an atom is equal to the balance of electronic energies. Electronic energies should show jumps and drops when moving from one atomic shell to another and ionization energy is the only part of total energy that can be obtained experimentally. Hence, the variations of experimental ionization energies between consecutive electrons were utilized to investigate the order of atomic shells. The variations of ionization energy were drawn versus the electron numbers to find abrupt changes in the energy, which are so-called "peaks". Then the observed peaks in the graphs were recorded as evidence for the order of atomic shells. The observed order of peaks did not completely match and support the order of atomic shells given by the well-known aufbau (or Madelung) rule. Thus, the observation is reported and a consistent view of the periodic table according to the new order is presented.

Coarse-graining of many-body path integrals: Theory and numerical approximations

Feynman's imaginary time path integral approach to quantum statistical mechanics provides a theoretical formalism for including nuclear quantum effects (NQEs) in simulation of condensed matter systems. Sinitskiy and Voth [J. Chem. Phys. 143, 094104 (2015)] have presented the coarse-grained path integral (CG-PI) theory, which provides a reductionist coarse-grained representation of the imaginary time path integral based on the quantum-classical isomorphism. In this paper, the many-body generalization of the CG-PI theory is presented. It is shown that the N interacting particles obeying quantum Boltzmann statistics can be represented as a system of N pairs of classical-like pseudoparticles coupled to each other analogous to the pseudoparticle pair of the one-body theory. Moreover, we present a numerical CG-PI (n-CG-PI) method applying a simple approximation to the coupling scheme between the pseudoparticles due to numerical challenges of directly implementing the full many-body CG-PI theory. Structural correlations of two liquid systems are investigated to demonstrate the performance of the n-CG-PI method. Both the many-body CG-PI theory and the n-CG-PI method not only present reductionist views of the many-body quantum Boltzmann statistics but also provide theoretical and numerical insight into how to explicitly incorporate NQEs in the representation of condensed matter systems with minimal additional degrees of freedom.

Structure and Dynamics of Water at the Water-Air Interface Using First-Principles Molecular Dynamics Simulations. II. NonLocal vs Empirical van der Waals Corrections

van der Waals (vdW) correction schemes have been recognized to be essential for an accurate description of liquid water in first-principles molecular dynamics simulation. The description of the structure and dynamics of water is governed by the type of the vdW corrections. So far, two vdW correction schemes have been often used: empirical vdW corrections and nonlocal vdW corrections. In this paper, we assess the influence of the empirical vs nonlocal vdW correction schemes on the structure and dynamics of water at the water–air interface. Since the structure of water at the water–air interface is established by a delicate balance of hydrogen bond formation and breaking, the simulation at the water–air interface provides a unique platform to testify as to the heterogeneous interaction of water. We used the metrics [Ohto et al. J. Chem. Theory Comput., 2019, 15, 595−60230468702] which are directly connected with the sum-frequency generation spectroscopic measurement. We find that the overall performance of nonlocal vdW methods is either similar or worse compared to the empirical vdW methods. We also investigated the performance of the optB88-DRSLL functional, which showed slightly less accuracy than the revPBE-D3 method. We conclude that the revPBE-D3 method shows the best performance for describing the interfacial water.

Understanding the effect of nanoconfinement on the structure of water hydrogen bond networks

Dynamic hydrogen bond trails in water confined between two phospholipid membranes traced by the information flow model.

Using a monomer potential energy surface to perform approximate path integral molecular dynamics simulation of ab initio water at near-zero added cost

It is now established that nuclear quantum motion plays an important role in determining water's hydrogen bonding, structure, and dynamics. Such effects are important to include in density functional theory (DFT) based molecular dynamics simulation of water. The standard way of treating nuclear quantum effects, path integral molecular dynamics (PIMD), multiplies the number of energy/force calculations by the number of beads required. In this work we introduce a method whereby PIMD can be incorporated into a DFT simulation with little extra cost and little loss in accuracy. The method is based on the many body expansion of the energy and has the benefit of including a monomer level correction to the DFT energy. Our method calculates intramolecular forces using the highly accurate monomer potential energy surface developed by Partridge-Schwenke, which is cheap to evaluate. Intermolecular forces and energies are calculated with DFT only once per timestep using the centroid positions. We show how our method may be used in conjunction with a multiple time step algorithm for an additional speedup and how it relates to ring polymer contraction and other schemes that have been introduced recently to speed up PIMD simulations. We show that our method, which we call "monomer PIMD", correctly captures changes in the structure of water found in a full PIMD simulation but at much lower computational cost.

Converged Colored Noise Path Integral Molecular Dynamics Study of the Zundel Cation Down to Ultralow Temperatures at Coupled Cluster Accuracy

For a long time, performing converged path integral simulations at ultralow but finite temperatures of a few Kelvin has been a nearly impossible task. However, recent developments in advanced colored noise thermostatting schemes for path integral simulations, namely, the Path Integral Generalized Langevin Equation Thermostat (PIGLET) and the Path Integral Quantum Thermal Bath (PIQTB), have been able to greatly reduce the computational cost of these simulations, thus making the ultralow temperature regime accessible in practice. In this work, we investigate the influence of these two thermostatting schemes on the description of hydrogen-bonded systems at temperatures down to a few Kelvin as encountered, for example, in helium nanodroplet isolation or tagging photodissociation spectroscopy experiments. For this purpose, we analyze the prototypical hydrogen bond in the Zundel cation (H₅O₂⁺) as a function of both oxygen-oxygen distance and temperature in order to elucidate how the anisotropic quantum delocalization and, thus, the shape of the shared proton adapts depending on the donor-acceptor distance. The underlying electronic structure of the Zundel cation is described in terms of Behler's Neural Network Potentials of essentially converged Coupled Cluster accuracy, CCSD(T*)-F12a/AVTZ. In addition, the performances of the PIQTB and PIGLET methods for energetic, structural, and quantum delocalization properties are assessed and directly compared. Overall, our results emphasize the validity and practical usefulness of these two modern thermostatting approaches for path integral simulations of hydrogen-bonded systems even at ultralow temperatures.

The Quest for Accurate Liquid Water Properties from First Principles

Developing accurate ab initio molecular dynamics (AIMD) models that capture both electronic reorganization and nuclear quantum effects associated with hydrogen bonding is key to quantitative understanding of bulk water and its anomalies as well as its role as a universal solvent. For condensed phase simulations, AIMD has typically relied on the generalized gradient approximation (GGA) of density functional theory (DFT) as the underlying model chemistry for the potential energy surface, with nuclear quantum effects (NQEs) sometimes modeled by performing classical molecular dynamics simulations at elevated temperatures. Here we show that the properties of liquid water obtained from the meta-GGA B97M-rV functional, when evaluated using accelerated path integral molecular dynamics simulations, display accuracy comparable to a computationally expensive dispersion-corrected hybrid functional, revPBE0-D3. We show that the meta-GGA DFT functional reproduces bulk water properties including radial distribution functions, self-diffusion coefficients, and infrared spectra with comparable accuracy of a much more expensive hybrid functional. This work demonstrates that the underlying quality of a good DFT functional requires evaluation with quantum nuclei and that high-temperature simulations are a poor proxy for properly treating NQEs.

On the importance of accounting for nuclear quantum effects in ab initio calibrated force fields in biological simulations

In many important processes in chemistry, physics, and biology the nuclear degrees of freedom cannot be described using the laws of classical mechanics. At the same time, the vast majority of molecular simulations that employ wide-coverage force fields treat atomic motion classically. In light of the increasing desire for and accelerated development of quantum mechanics (QM)-parameterized interaction models, we reexamine whether the classical treatment is sufficient for a simple but crucial chemical species: alkanes. We show that when using an interaction model or force field in excellent agreement with the "gold standard" QM data, even very basic simulated properties of liquid alkanes, such as densities and heats of vaporization, deviate significantly from experimental values. Inclusion of nuclear quantum effects via techniques that treat nuclear degrees of freedom using the laws of classical mechanics brings the simulated properties much closer to reality.

Unraveling electronic absorption spectra using nuclear quantum effects: Photoactive yellow protein and green fluorescent protein chromophores in water

Many physical phenomena must be accounted for to accurately model solution-phase optical spectral line shapes, from the sampling of chromophore-solvent configurations to the electronic-vibrational transitions leading to vibronic fine structure. Here we thoroughly explore the role of nuclear quantum effects, direct and indirect solvent effects, and vibronic effects in the computation of the optical spectrum of the aqueously solvated anionic chromophores of green fluorescent protein and photoactive yellow protein. By analyzing the chromophore and solvent configurations, the distributions of vertical excitation energies, the absorption spectra computed within the ensemble approach, and the absorption spectra computed within the ensemble plus zero-temperature Franck-Condon approach, we show how solvent, nuclear quantum effects, and vibronic transitions alter the optical absorption spectra. We find that including nuclear quantum effects in the sampling of chromophore-solvent configurations using ab initio path integral molecular dynamics simulations leads to improved spectral shapes through three mechanisms. The three mechanisms that lead to line shape broadening and a better description of the high-energy tail are softening of heavy atom bonds in the chromophore that couple to the optically bright state, widening the distribution of vertical excitation energies from more diverse solvation environments, and redistributing spectral weight from the 0-0 vibronic transition to higher energy vibronic transitions when computing the Franck-Condon spectrum in a frozen solvent pocket. The absorption spectra computed using the combined ensemble plus zero-temperature Franck-Condon approach yield significant improvements in spectral shape and width compared to the spectra computed with the ensemble approach. Using the combined approach with configurations sampled from path integral molecular dynamics trajectories presents a significant step forward in accurately modeling the absorption spectra of aqueously solvated chromophores.

Decoding the spectroscopic features and time scales of aqueous proton defects

Acid solutions exhibit a variety of complex structural and dynamical features arising from the presence of multiple interacting reactive proton defects and counterions. However, disentangling the transient structural motifs of proton defects in the water hydrogen bond network and the mechanisms for their interconversion remains a formidable challenge. Here, we use simulations treating the quantum nature of both the electrons and nuclei to show how the experimentally observed spectroscopic features and relaxation time scales can be elucidated using a physically transparent coordinate that encodes the overall asymmetry of the solvation environment of the proton defect. We demonstrate that this coordinate can be used both to discriminate the extremities of the features observed in the linear vibrational spectrum and to explain the molecular motions that give rise to the interconversion time scales observed in recent nonlinear experiments. This analysis provides a unified condensed-phase picture of the proton structure and dynamics that, at its extrema, encompasses proton sharing and spectroscopic features resembling the limiting Eigen [H₃O(H₂O)₃]⁺ and Zundel [H(H₂O)₂]⁺ gas-phase structures, while also describing the rich variety of interconverting environments in the liquid phase.

From classical to quantum and back: Hamiltonian adaptive resolution path integral, ring polymer, and centroid molecular dynamics

Path integral-based methodologies play a crucial role for the investigation of nuclear quantum effects by means of computer simulations. However, these techniques are significantly more demanding than corresponding classical simulations. To reduce this numerical effort, we recently proposed a method, based on a rigorous Hamiltonian formulation, which restricts the quantum modeling to a small but relevant spatial region within a larger reservoir where particles are treated classically. In this work, we extend this idea and show how it can be implemented along with state-of-the-art path integral simulation techniques, including path-integral molecular dynamics, which allows for the calculation of quantum statistical properties, and ring-polymer and centroid molecular dynamics, which allow the calculation of approximate quantum dynamical properties. To this end, we derive a new integration algorithm that also makes use of multiple time-stepping. The scheme is validated via adaptive classical-path-integral simulations of liquid water. Potential applications of the proposed multiresolution method are diverse and include efficient quantum simulations of interfaces as well as complex biomolecular systems such as membranes and proteins.

First-Principles Simulations of Liquid Water Using a Dielectric-Dependent Hybrid Functional

We carried out first-principles simulations of liquid water under ambient conditions using a dielectric-dependent hybrid functional, where the fraction of exact exchange is set equal to the inverse of the high-frequency dielectric constant of the liquid. We found excellent agreement with experiment for the oxygen-oxygen partial correlation function at the experimental equilibrium density and 311 ± 3 K. Other structural and dynamical properties, such as the diffusion coefficient, molecular dipole moments, and vibrational spectra, are also in good agreement with experiment. Our results, together with previous findings on electronic properties of the liquid with the same functional, show that the dielectric-dependent hybrid functional accurately describes both the structural and electronic properties of liquid water.

Computational Methods for Ab Initio Molecular Dynamics

Ab initio molecular dynamics is an irreplaceable technique for the realistic simulation of complex molecular systems and processes from first principles. This paper proposes a comprehensive and self-contained review of ab initio molecular dynamics from a computational perspective and from first principles. Quantum mechanics is presented from a molecular dynamics perspective. Various approximations and formulations are proposed, including the Ehrenfest, Born–Oppenheimer, and Hartree–Fock molecular dynamics. Subsequently, the Kohn–Sham formulation of molecular dynamics is introduced as well as the afferent concept of density functional. As a result, Car–Parrinello molecular dynamics is discussed, together with its extension to isothermal and isobaric processes. Car–Parrinello molecular dynamics is then reformulated in terms of path integrals. Finally, some implementation issues are analysed, namely, the pseudopotential, the orbital functional basis, and hybrid molecular dynamics.

Hydrogen-bond structure dynamics in bulk water: insights from ab initio simulations with coupled cluster theory

An accurate and efficient ab initio molecular dynamics (AIMD) simulation of liquid water was made possible using the fragment-based approach (J. F. Liu, X. He and J. Z. H. Zhang, Phys. Chem. Chem. Phys., 2017, 19, 11931-11936). In this study, we advance the AIMD simulations using the fragment-based coupled cluster (CC) theory, more accurately revealing the structural and dynamical properties of liquid water under ambient conditions. The results show that the double-donor hydrogen-bond configurations in liquid water are nearly in balance with the single-donor configurations, with a slight bias towards the former. Our observation is in contrast to the traditional tetrahedral water structure. The hydrogen-bond switching dynamics in liquid water are very fast, with a hydrogen-bond life time of around 0.78 picoseconds, determined using AIMD simulation at the CCD/aug-cc-pVDZ level. This time scale is remarkably shorter than the ∼3.0 picoseconds that is commonly obtained from traditional nonpolarized force fields and density functional theory (DFT) based first-principles simulations. Additionally, the obtained radial distribution functions, triplet oxygen angular distribution, diffusion coefficient, and the dipole moment of the water molecule are uniformly in good agreement with the experimental observations. The current high-level AIMD simulation sheds light on the understanding of the structural and dynamical properties of liquid water.