Path Integral Molecular Dynamics Research Articles

Zero-cost corrections to influence functional coefficients from bath response functions.

Recent work has shown that it is possible to circumvent the calculation of the spectral density and directly generate the coefficients of the discretized influence functionals using data from classical trajectory simulations. However, the accuracy of this procedure depends on the validity of the high temperature approximation. In this work, an alternative derivation based on the Kubo formalism is provided. This enables the calculation of additional correction terms that increases the range of applicability of the procedure to lower temperatures. Because it is based on the Kubo-transformed correlation function, this approach allows the direct use of correlation functions obtained from methods such as ring-polymer molecular dynamics and centroid molecular dynamics in determining the influence functional coefficients for subsequent system-solvent simulations. The accuracy of the original procedure and the corrected procedure is investigated across a range of parameters. It is interesting that the correction term comes at zero additional cost. Furthermore, it is possible to improve upon the correction using zero-cost physical intuition and heuristics making the method even more accurate.

Equilibrium and Dynamical Characteristics of Hydrogen Bond Bifurcations in Water-Water and Water-Ammonia Dimers: A Path Integral Molecular Dynamics Study.

We present path integral molecular dynamics results that describe the effects of nuclear quantum fluctuations on equilibrium and dynamical characteristics pertaining to bifurcation pathways in hydrogen bonded dimers combining water and ammonia, at cryogenic temperatures of the order of 20 K. Along these isomerizations, the hydrogen atoms in the molecules acting as hydrogen-bond donors interchange their original dangling/connective characters. Our results reveal that the resulting quantum transition paths comprise three stages: the initial and final ones involve overall rotations during which the two protons retain their classical-like characteristics. Effects from quantum fluctuation are clearly manifested in the changes operated at the intermediate passages over transition states, as the spatial extents of the protons stretch over typical lengths comparable to the distances between connective and dangling basins of attractions. Consequently, the classical over-the-hill path is replaced by a tunneling controlled mechanism which, within the path integral perspective, can be cast in terms of concerted inter-basin migrations of polymer beads from dangling-to-connective and from connective-to-dangling, at practically no energy costs. We also estimated the characteristic timescales describing such interconversions within the approximate ring polymer rate theory. Effects derived from full and partial deuteration are also discussed.

Clustering, collision, and relaxation dynamics in pure and doped helium nanoclusters: Density- vs particle-based approaches.

The clustering, collision, and relaxation dynamics of pristine and doped helium nanodroplets is theoretically investigated in cases of pickup and clustering of heliophilic argon, collision of heliophobic cesium atoms, and coalescence of two droplets brought into contact by their mutual long-range van der Waals interaction. Three approaches are used and compared with each other. The He time-dependent density functional theory method considers the droplet as a continuous medium and accounts for its superfluid character. The ring-polymer molecular dynamics method uses a path-integral description of nuclear motion and incorporates zero-point delocalization while bosonic exchange effects are ignored. Finally, the zero-point averaged dynamics approach is a mixed quantum-classical method in which quantum delocalization is described by attaching a frozen wavefunction to each He atom, equivalent to classical dynamics with effective interaction potentials. All three methods predict that the growth of argon clusters is significantly hindered by the helium host droplet due to the impeding shell structure around the dopants and kinematic effects freezing the growing cluster in metastable configurations. The effects of superfluidity are qualitatively manifested by different collision dynamics of the heliophilic atom at high velocities, as well as quadrupole oscillations that are not seen with particle-based methods, for droplets experiencing a collision with cesium atoms or merging with each other.

Artificial neural network-based path integral simulations of hydrogen isotope diffusion in palladium

The contribution of nuclear quantum effects (NQEs) to the kinetics and dynamics of interstitial H isotopes in face-centered cubic Pd was intensively investigated using several path-integral techniques, along with a newly developed machine-learning interatomic potential based on artificial neural networks for Pd–H alloys. The diffusion coefficients (D) of protium, deuterium, and tritium in Pd were predicted over a wide temperature range (50–1500 K) based on quantum transition-state theory (QTST) combined with path-integral molecular-dynamics simulations. The importance of NQEs even at high temperatures was illustrated in terms of the characteristic temperature dependence of the activation free energies for H-isotope migration in Pd. This illuminates the overall picture of anomalous D crossovers among the three H isotopes in Pd. In addition, the D of protium in Pd was directly computed using two approximate quantum-dynamics methods based on Feynman’s path-integral theory, i.e. centroid molecular dynamics (CMD) and ring-polymer molecular dynamics (RPMD), in the temperature range 370–1500 K. The D values obtained from the CMD and RPMD simulations were very similar and agreed better with the reported experimental values than the QTST results in this temperature range. Our machine learning-based path-integral calculations elucidate the underlying quantum nature of the ‘reversed S’-type nonlinear behavior of D for the three H isotopes in Pd on the Arrhenius plots.

BIGDML—Towards accurate quantum machine learning force fields for materials

Machine-learning force fields (MLFF) should be accurate, computationally and data efficient, and applicable to molecules, materials, and interfaces thereof. Currently, MLFFs often introduce tradeoffs that restrict their practical applicability to small subsets of chemical space or require exhaustive datasets for training. Here, we introduce the Bravais-Inspired Gradient-Domain Machine Learning (BIGDML) approach and demonstrate its ability to construct reliable force fields using a training set with just 10–200 geometries for materials including pristine and defect-containing 2D and 3D semiconductors and metals, as well as chemisorbed and physisorbed atomic and molecular adsorbates on surfaces. The BIGDML model employs the full relevant symmetry group for a given material, does not assume artificial atom types or localization of atomic interactions and exhibits high data efficiency and state-of-the-art energy accuracies (errors substantially below 1 meV per atom) for an extended set of materials. Extensive path-integral molecular dynamics carried out with BIGDML models demonstrate the counterintuitive localization of benzene–graphene dynamics induced by nuclear quantum effects and their strong contributions to the hydrogen diffusion coefficient in a Pd crystal for a wide range of temperatures.

Evidence for Dissociation and Ionization in Shock Compressed Nitrogen to 800GPa.

Triple bonding in the nitrogen molecule (N_{2}) is among the strongest chemical bonds with a dissociation enthalpy of 9.8 eV/molecule. Nitrogen is therefore an excellent test bed for theoretical and numerical methods aimed at understanding how bonding evolves under the influence of the extreme pressures and temperatures of the warm dense matter regime. Here, we report laser-driven shock experiments on fluid molecular nitrogen up to 800GPa and 4.0 g/cm^{3}. Line-imaging velocimetry measurements and impedance matching method with a quartz reference yield shock equation of state data of initially precompressed nitrogen. Comparison with numerical simulations using path integral MonteCarlo and density functional theory molecular dynamics reveals clear signatures of chemical dissociation and the onset of L-shell ionization. Combining data along multiple shock Hugoniot curves starting from densities between 0.76 and 1.29 g/cm^{3}, our study documents how pressure and density affect these changes in chemical bonding and provides benchmarks for future theoretical developments in this regime, with applications for planetary interior modeling, high energy density science, and inertial confinement fusion research.

Kinetics study of the OH + SiH4 hydrogen abstraction reaction: A theoretical analysis

AbstractThe temperature dependence of the rate constants and the kinetics isotope effects of the OH + SiH4 hydrogen abstraction reaction were theoretically studied. Three different kinetics approaches were used in this study: variational transition state theory with multidimensional tunneling (VTST/MT) corrections, quasi‐classical trajectory (QCT) calculations, and ring polymer molecular dynamics (RPMD). These studies were performed on a global analytical potential energy surface based on explicitly correlated high level ab initio calculations developed in the present work. While the experimental rate constants are practically independent with temperature, the VTST/MT results show a pronounced dependence, the QCT partially corrects this behavior, and RPMD rates are in better agreement with experiment. In this case, in the experimental common temperature range, 300–600 K, the T‐expression for the RPMD rate constants is: k(T) = 6.40 × 10−13.T0.52exp(−0.10/T). The kinetics isotope effects, H/D, are small in this temperature range, revealing the small importance of tunneling effects in this very exothermic and low barrier reaction. Different factors, such as anharmonicity, the zero‐point violation problem or limitations of the surface, without ruling out experimental uncertainties, were analyzed to explain the discrepancies.

Numerical calculation of Green's function and momentum distribution for spin-polarized fermions by path integral molecular dynamics.

Most recently, path integral molecular dynamics (PIMD) has been successfully applied to perform simulations of identical bosons and fermions by Hirshberg et al. In this work, we demonstrate that PIMD can be developed to calculate Green's function and extract momentum distributions for spin-polarized fermions. In particular, we show that the momentum distribution calculated by PIMD has potential applications to numerous quantum systems, e.g., ultracold fermionic atoms in optical lattices.

Nuclear quantum effects on the dynamics and glass behavior of a monatomic liquid with two liquid states.

We perform path integral molecular dynamics (PIMD) simulations of a monatomic liquid that exhibits a liquid-liquid phase transition and liquid-liquid critical point. PIMD simulations are performed using different values of Planck's constant h, allowing us to study the behavior of the liquid as nuclear quantum effects (NQE, i.e., atoms delocalization) are introduced, from the classical liquid (h = 0) to increasingly quantum liquids (h > 0). By combining the PIMD simulations with the ring-polymer molecular dynamics method, we also explore the dynamics of the classical and quantum liquids. We find that (i) the glass transition temperature of the low-density liquid (LDL) is anomalous, i.e., Tg LDL(P) decreases upon compression. Instead, (ii) the glass transition temperature of the high-density liquid (HDL) is normal, i.e., Tg HDL(P) increases upon compression. (iii) NQE shift both Tg LDL(P) and Tg HDL(P) toward lower temperatures, but NQE are more pronounced on HDL. We also study the glass behavior of the ring-polymer systems associated with the quantum liquids studied (via the path-integral formulation of statistical mechanics). There are two glass states in all the systems studied, low-density amorphous ice (LDA) and high-density amorphous ice (HDA), which are the glass counterparts of LDL and HDL. In all cases, the pressure-induced LDA-HDA transformation is sharp, reminiscent of a first-order phase transition. In the low-quantum regime, the LDA-HDA transformation is reversible, with identical LDA forms before compression and after decompression. However, in the high-quantum regime, the atoms become more delocalized in the final LDA than in the initial LDA, raising questions on the reversibility of the LDA-HDA transformation.

Quantum and temperature effects on the crystal structure of superhydride LaH10 : A path integral molecular dynamics study

By classical and path-integral molecular dynamics simulations, we study the pressure-temperature ($P$-$T$) phase diagram of LaH$_{10}$ to clarify the impact of temperature and atomic zero-point motions. We calculate the XRD pattern and analyze the space group of the crystal structures. For 125 GPa $\leq P\leq$ 150 GPa and $T=300$ K, we show that a highly symmetric $Fm\bar{3}m$ structure, for which superconductivity is particularly favored, is stabilized only by the temperature effect. On the other hand, for $T=200$ K, the interplay between the temperature and quantum effects is crucial to realize the $Fm\bar{3}m$ structure. For $P=$100 GPa and $T=$300 K, we find that the system is close to the critical point of the structural phase transition between the $Fm\bar{3}m$ structure and those with lower symmetries.

Learning the Quantum Centroid Force Correction in Molecular Systems: A Localized Approach.

Molecular mechanics (MM) is a powerful tool to study the properties of molecular systems in the fields of biology and materials science. With the development of ab initio force field and the application of ab initio potential energy surface, the nuclear quantum effect (NQE) is becoming increasingly important for the robustness of the simulation. However, the state-of-the-art path-integral molecular dynamics simulation, which incorporates NQE in MM, is still too expensive to conduct for most biological and material systems. In this work, we analyze the locality of NQE, using both analytical and numerical approaches, and conclude that NQE is an extremely localized phenomenon in nonreactive molecular systems. Therefore, we can use localized machine learning (ML) models to predict quantum force corrections both accurately and efficiently. Using liquid water as example, we show that the ML facilitated centroid MD can reproduce the NQEs in both the thermodynamical and the dynamical properties, with a minimal increase in computational time compared to classical molecular dynamics. This simple approach thus largely decreases the computational cost of quantum simulations, making it really accessible to the studies of large-scale molecular systems.

Microscopic collective dynamics in liquid neon-deuterium mixtures: Inelastic neutron scattering and quantum simulations.

In this paper a combined neutron scattering and quantum simulation study of the collective dynamics in liquid Ne-D_{2} mixtures, at a temperature of T=30K and in the wave-vector transfer range 4 nm^{-1}<q<51nm^{-1}, is presented. Two D_{2} concentrations are investigated, one close to 25% molar and the other close to 50% molar, together with pure Ne. The dynamic structure factor for the centers of mass of the two molecular species is extracted from the neutron scattering data and subsequently compared with that obtained from three different quantum simulation methods, such as ring polymer molecular dynamics and two slightly different versions of the Feynman-Kleinert approach. A general agreement is found, even though some discrepancies both among simulations, and between simulations and experimental data, can be observed. In order to clarify the physical meaning of the present spectroscopic results, an analysis of the longitudinal current spectral maxima is carried out showing the peculiarities of the D_{2} center-of-mass dynamics in these mixtures. A comparison with the centroid molecular dynamics results obtained for the D_{2} center-of-mass self-dynamics in the same liquid mixtures is finally proposed.

Importance of vibrational anharmonicity for electron-phonon coupling in molecular crystals

Accurate predictions of electronic band gaps are key to the computational design of functional materials with tailored optical and thermoelectric properties. While an accurate description of the underlying electronic structure is indispensable, electron-phonon interactions also play a prominent role in determining the band gap. The harmonic vibrational approximation is used almost universally to describe electron-phonon coupling and provides a satisfactory description in most systems. Here, we demonstrate that this is not the case for molecular crystals due to the presence of strongly anharmonic motions. We demonstrate the breakdown of the harmonic approximation in acene molecular crystals, and show that converged band gaps can be obtained by rigorously accounting for vibrational anharmonicity using path-integral molecular dynamics. Finally, we characterize the most anharmonic vibrational modes and their contributions to the band-gap renormalization.

Nuclear Quantum Effects in H2 Adsorption Dynamics on a Small Water Cluster Studied with Ring-Polymer Molecular Dynamics Simulations

Nuclear Quantum Effects in H₂ Adsorption Dynamics on a Small Water Cluster Studied with Ring-Polymer Molecular Dynamics Simulations

Quantum Simulations of Vibrational Strong Coupling via Path Integrals.

A quantum simulation of vibrational strong coupling (VSC) in the collective regime via thermostated ring-polymer molecular dynamics (TRPMD) is reported. For a collection of liquid-phase water molecules resonantly coupled to a single lossless cavity mode, the simulation shows that as compared with a fully classical calculation, the inclusion of nuclear and photonic quantum effects does not lead to a change in the Rabi splitting but does broaden polaritonic line widths roughly by a factor of 2. Moreover, under thermal equilibrium, both quantum and classical simulations predict that the static dielectric constant of liquid water is largely unchanged inside vs outside the cavity. This result disagrees with a recent experiment demonstrating that the static dielectric constant of liquid water can be resonantly enhanced under VSC, suggesting either limitations of our approach or perhaps other experimental factors that have not yet been explored.

Low-Barrier Hydrogen Bond in Fujikurin A-D: A Computational Study.

The compounds Fujikurin A, B, and D, recently isolated from Fusarium fujikuroi, possess intramolecular low-barrier hydrogen bonds (LBHBs), which are hydrogen bonds with a very low-energy barrier for proton transfer. The isolated compounds have a hydrogen-bonded proton that appears to rapidly switch between two equilibrium states via a transition state (TS). To understand the characteristics of these intramolecular LBHBs in detail, we performed path integral molecular dynamics (PIMD) simulations, which can consider nuclear quantum effects (NQEs) under a finite temperature. The PIMD simulations predicted that the NQE completely washed out the energy barrier for the proton transfer reaction. Consequently, a single-well shape emerged in the results, along with the effective free-energy potential surface for the hydrogen-bonded proton distribution. Thus, we conclude that the hydrogen-bonded proton in Fujikurin does not in fact transfer between two equilibrium structures but widely delocalizes around the global minimum structure involving the TS region.

Evidence of a liquid\u2013liquid phase transition in H_2O and D_2O from path-integral molecular dynamics simulations

We perform path-integral molecular dynamics (PIMD), ring-polymer MD (RPMD), and classical MD simulations of H_2O and D_2O using the q-TIP4P/F water model over a wide range of temperatures and pressures. The density rho (T), isothermal compressibility kappa _T(T), and self-diffusion coefficients D(T) of H_2O and D_2O are in excellent agreement with available experimental data; the isobaric heat capacity C_P(T) obtained from PIMD and MD simulations agree qualitatively well with the experiments. Some of these thermodynamic properties exhibit anomalous maxima upon isobaric cooling, consistent with recent experiments and with the possibility that H_2O and D_2O exhibit a liquid-liquid critical point (LLCP) at low temperatures and positive pressures. The data from PIMD/MD for H_2O and D_2O can be fitted remarkably well using the Two-State-Equation-of-State (TSEOS). Using the TSEOS, we estimate that the LLCP for q-TIP4P/F H_2O, from PIMD simulations, is located at P_c = 167 pm 9 MPa, T_c = 159 pm 6 K, and rho _c = 1.02 pm 0.01 g/cm^3. Isotope substitution effects are important; the LLCP location in q-TIP4P/F D_2O is estimated to be P_c = 176 pm 4 MPa, T_c = 177 pm 2 K, and rho _c = 1.13 pm 0.01 g/cm^3. Interestingly, for the water model studied, differences in the LLCP location from PIMD and MD simulations suggest that nuclear quantum effects (i.e., atoms delocalization) play an important role in the thermodynamics of water around the LLCP (from the MD simulations of q-TIP4P/F water, P_c = 203 pm 4 MPa, T_c = 175 pm 2 K, and rho _c = 1.03 pm 0.01 g/cm^3). Overall, our results strongly support the LLPT scenario to explain water anomalous behavior, independently of the fundamental differences between classical MD and PIMD techniques. The reported values of T_c for D_2O and, particularly, H_2O suggest that improved water models are needed for the study of supercooled water.

Equilibrium-nonequilibrium ring-polymer molecular dynamics for nonlinear spectroscopy.

Two-dimensional Raman and hybrid terahertz-Raman spectroscopic techniques provide invaluable insight into molecular structures and dynamics of condensed-phase systems. However, corroborating experimental results with theory is difficult due to the high computational cost of incorporating quantum-mechanical effects in the simulations. Here, we present the equilibrium-nonequilibrium ring-polymer molecular dynamics (RPMD), a practical computational method that can account for nuclear quantum effects on the two-time response function of nonlinear optical spectroscopy. Unlike a recently developed approach based on the double Kubo transformed (DKT) correlation function, our method is exact in the classical limit, where it reduces to the established equilibrium-nonequilibrium classical molecular dynamics method. Using benchmark model calculations, we demonstrate the advantages of the equilibrium-nonequilibrium RPMD over classical and DKT-based approaches. Importantly, its derivation, which is based on the nonequilibrium RPMD, obviates the need for identifying an appropriate Kubo transformed correlation function and paves the way for applying real-time path-integral techniques to multidimensional spectroscopy.

Path integral molecular dynamics simulations for Green's function in a system of identical bosons.

Path integral molecular dynamics (PIMD) has been successfully applied to perform simulations of large bosonic systems in a recent study [Hirshberg et al., Proc. Natl. Acad. Sci. U. S. A. 116, 21445 (2019)]. In this work, we extend PIMD techniques to study Green's function for bosonic systems. We demonstrate that the development of the original PIMD method enables us to calculate Green's function and extract momentum distribution from our simulations. We also apply our method to systems of identical interacting bosons to study Berezinskii-Kosterlitz-Thouless transition around its critical temperature.

Revealing Intra- and Intermolecular Interactions Determining Physico-Chemical Features of Selected Quinolone Carboxylic Acid Derivatives.

The intra- and intermolecular interactions of selected quinolone carboxylic acid derivatives were studied in monomers, dimers and crystals. The investigated compounds are well-recognized as medicines or as bases for further studies in drug design. We employed density functional theory (DFT) in its classical formulation to develop gas-phase and solvent reaction field (PCM) models describing geometric, energetic and electronic structure parameters for monomers and dimers. The electronic structure was investigated based on the atoms in molecules (AIM) and natural bond orbital (NBO) theories. Special attention was devoted to the intramolecular hydrogen bonds (HB) present in the investigated compounds. The characterization of energy components was performed using symmetry-adapted perturbation theory (SAPT). Finally, the time-evolution methods of Car–Parrinello molecular dynamics (CPMD) and path integral molecular dynamics (PIMD) were employed to describe the hydrogen bond dynamics as well as the spectroscopic signatures. The vibrational features of the O-H stretching were studied using Fourier transformation of the autocorrelation function of atomic velocity. The inclusion of quantum nuclear effects provided an accurate depiction of the bridged proton delocalization. The CPMD and PIMD simulations were carried out in the gas and crystalline phases. It was found that the polar environment enhances the strength of the intramolecular hydrogen bonds. The SAPT analysis revealed that the dispersive forces are decisive factors in the intermolecular interactions. In the electronic ground state, the proton-transfer phenomena are not favourable. The CPMD results showed generally that the bridged proton is localized at the donor side, with possible proton-sharing events in the solid-phase simulation of stronger hydrogen bridges. However, the PIMD enabled the quantitative estimation of the quantum effects inclusion—the proton position was moved towards the bridge midpoint, but no qualitative changes were detected. It was found that the interatomic distance between the donor and acceptor atoms was shortened and that the bridged proton was strongly delocalized.