Path Integral Molecular Dynamics Simulations Research Articles

Coarse-Graining of Imaginary Time Feynman Path Integrals: Inclusion of Intramolecular Interactions and Bottom-up Force-Matching.

Feynman's imaginary time path integral formalism of quantum statistical mechanics and the corresponding quantum-classical isomorphism provide a tangible way of incorporating nuclear quantum effect (NQE) in the simulation of condensed matter systems using well-developed classical simulation techniques. Our previous work has presented the many-body coarse-graining of path integral (CG-PI) theory that builds an isomorphism between the quantum partition function of N distinguishable particles and the classical partition function of 2N pseudoparticles. In this present work, we develop a generalized version of the many-body CG-PI theory that incorporates many-body interactions in the force field. Based on the new derivation, we provide a numerical CG-PI (n-CG-PI) modeling strategy parametrized from the underlying path integral molecular dynamics (PIMD) trajectories using force matching and Boltzmann inversion. The n-CG-PI models for two liquid systems are shown to capture well both the intramolecular and intermolecular structural correlations of the reference PIMD simulations. The generalized derivation of the many-body CG-PI theory and the n-CG-PI model presented in this work extend the scope of the CG-PI formalism by generalizing the previously limited theory to incorporate force fields of realistic molecular systems.

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.

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.

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.

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.

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.

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.

Nonlocal Screening Dictates the Radiative Lifetimes of Excitations in Lead Halide Perovskites.

We use path integral molecular dynamics simulations and theory to elucidate the interactions between charge carriers, as mediated by a lead halide perovskite lattice. We find that the charge-lattice coupling of MAPbI3 results in a repulsive interaction between electrons and holes at intermediate distances. The effective interaction is understood using a Gaussian field theory, whereby the underlying soft, polar lattice contributes a nonlocal screening between quasiparticles. Path integral calculations of this nonlocal screening model are used to rationalize the small exciton binding energy and low radiative recombination rate observed experimentally and are compared to traditional Wannier-Mott and Fröhlich models, which fail to do so. These results clarify the origin of the high power conversion efficiencies in lead halide perovskites. Emergent repulsive electron-hole interactions provide a design principle for optimizing soft, polar semiconductors.

Nuclear Quantum Effects on the Charge-Density Wave Transition in NbX2 (X = S, Se).

Understanding the origin of charge-density wave (CDW) instability is important for manipulating novel collective electronic states. Many layered transition metal dichalcogenides (TMDs) share similarity in the structural and electronic instability, giving rise to diverse CDW phases and superconductivity. It is still puzzling that even isostructural and isoelectronic TMDs show distinct CDW features. For instance, bulk NbSe2 exhibits CDW order at low temperature, while bulk NbS2 displays no CDW instability. The CDW transitions in single-layer NbS2 and NbSe2 are also different. In the classic limit, we investigate the electron correlation effects on the dimensionality dependence of the CDW ordering. By performing ab initio path integral molecular dynamics simulations and comparative analyses, we further revealed significant nuclear quantum effects in these systems. Specifically, the quantum motion of sulfur anions significantly reduces the CDW transition temperature in both bulk and single-layer NbS2, resulting in distinct CDW features in the NbS2 and NbSe2 systems.

Simulating Raman spectra of hydrogen hydrates using first-principles path-integral ring-polymer molecular dynamics

A variant of first-principles (FP) path-integral ring-polymer molecular dynamics (RPMD) simulations was used to investigate the C0,C1, and C2 H2 hydrates at 280 K aiming at reproducing their experimental Raman spectra. We find that our implementation based on white-noise Langevin thermostat applied to the H2O sub-system only in the formulation of FP-RPMD is able to reproduce with high accuracy the characteristic features of the measured Raman spectra of H2 hydrates in all C0,C1, and C2 phases, demonstrating potential of our computational approach for investigating H-containing materials under pressure.

Structures of liquid and aqueous water isotopologues at ambient temperature from ab initio path integral simulations

The heavy hydrogen isotopes D and T are found in trace amounts in water, and when their concentration increases they can play an intricate role in modulating the physical properties of the liquid. We present an analysis of the microscopic structures of ambient light water (H2O(l)), heavy water (D2O(l)), T2O(l), HDO(aq) and HTO(aq) studied by ab initio path integral molecular dynamics (PIMD). Unlike previous ab initio PIMD investigations of H2O(l) and D2O(l) [Chen et al., Phys. Rev. Lett., 2003, 91, 215503] [Machida et al., J. Chem. Phys., 2017, 148, 102324] we find that D2O(l) is more structured than H2O(l), as is predicted by the experiment. The agreement between the experiment and our simulation for H2O(l) and D2O(l) allows us to accurately predict the intra- and intermolecular structures of T2O(l) HDO(aq) and HTO(aq). T2O(l) is found to have a similar intermolecular structure to that of D2O(l), while the intramolecular structure is more compact, giving rise to a smaller dipole moment than those of H2O(l) and D2O(l). For the mixed isotope species, HDO(aq) and HTO(aq), we find smaller dipole moments and fewer hydrogen bonds when compared with the pure species H2O and D2O. We can attribute this effect to the relative compactness of the mixed isotope species, which results in a lower dipole moment than that of the pure species.

Atomic mass dependence of the nuclear quantum effect in NH4+(H2O)

• Path integral molecular dynamics (PIMD) was used to treat the nuclear quantum effect. • The ammonium cluster has two competing quantum effects. • The quantum fluctuation along the covalent bond strengthens the intermolecular interaction. • The quantum fluctuation along of bond angle decreases intermolecular interaction. The interaction between ammonium ions and water has been studied previously via molecular dynamics. However, the nuclear quantum effect phenomenon did not receive sufficient consideration therein. In this study, we report the ab initio path integral molecular dynamics (PIMD) simulations of the isotopologues of NH 4 + ( H 2 O ) , that is, NH 4 + ( H 2 O ) , ND 4 + ( D 2 O ) , and NT 4 + ( T 2 O ) , as well as the classical MD simulations. We found that the intermolecular distance does not change monotonically. Even though the nuclear quantum effect shortens the overall intermolecular distance (N–O atoms), the distance for NH 4 + ( H 2 O ) is longer than that for NT 4 + ( T 2 O ) . Two dimensional distributions of the binding energies with the contour lines were constructed to assess this unusual quantum behavior. We find that the elongation of the covalent bond increases the binding energy, while the quantum fluctuation of the O–N–H b bond angle decreases it. The bond length is determined by the delicate balance of such competing quantum effects.

Fast quasi-centroid molecular dynamics.

We describe a fast implementation of the quasi-centroid molecular dynamics (QCMD) method in which the quasi-centroid potential of mean force is approximated as a separable correction to the classical interaction potential. This correction is obtained by first calculating quasi-centroid radial and angular distribution functions in a short path integral molecular dynamics simulation and then using iterative Boltzmann inversion to obtain an effective classical potential that reproduces these distribution functions in a classical NVT simulation. We illustrate this approach with example applications to the vibrational spectra of gas phase molecules, obtaining excellent agreement with QCMD reference calculations for water and ammonia and good agreement with the quantum mechanical vibrational spectrum of methane.

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.

Competitive nuclear quantum effect and H/D isotope effect on torsional motion of H2O2: An ab initio path integral molecular dynamics study

Ab initio path integral molecular dynamics (PIMD) simulations were performed to understand the nuclear quantum effect and thermal effect on the torsional motion of H2O2 (D2O2). A broad distribution of geometrical parameters (OH and OO bond lengths [ROH and ROO] and HOO bond angle [∠HOO]) was obtained in quantum ab initio PIMD simulations, rather than in classical MD simulations. However, a modest distribution around the trans-configuration TS (∠HOOH = 180°) was found in the quantum ab initio PIMD simulations of H2O2 at low temperatures. Thus, the nuclear quantum effect slightly enhanced the HOOH torsional motion at low temperatures. A series of simulations also revealed that the thermal effect enhanced the HOOH torsional motion. However, it was restrained by the primary nuclear quantum effect on the ROH in quantum simulations at high temperatures. Therefore, an interesting competitive nuclear quantum effect on the torsional motion of H2O2 was observed.

Ab initio study of nuclear quantum effects on sub- and supercritical water.

The structures of water in the ambient, subcritical, and supercritical conditions at various densities were studied systematically by ab initio path integral molecular dynamics simulations. It was found that the nuclear quantum effects (NQEs) have a significant impact on the structure of hydrogen bonds in close contact, not only in the ambient condition but also in the sub- and supercritical conditions. The NQEs on the structure beyond the hydrogen bond contact are important in ambient water, but not much for water in the sub- and supercritical conditions. The NQEs are furthermore important for determining the number of hydrogen bonds in the ambient conditions, and this role is, however, diminished in the sub- and supercritical conditions. The NQEs do, nevertheless, show their importance in determining the intramolecular structure of water and the close contact structures of the hydrogen bonds, even at sub- and supercritical conditions. Using the RPBE-D3 functional, the computed radial distribution functions for ambient water are in excellent agreement with experimental data, upgrading our previous results using the BLYP-D2 functional [Machida et al., J. Chem. Phys. 148, 102324 (2018)]. The computed radial distribution functions for water in the sub- and supercritical conditions were carefully compared with experiment. In particular, we found that the first peak in hydrogen pair distribution functions matches only when the NQEs are taken into account.

Adsorption of helium on a charged propeller molecule: hexaphenylbenzene

Physisorption on planar or curved graphitic surfaces or aromatic rings has been investigated by various research groups, but in these studies, the substrate was usually strictly rigid. Here, we report a combined experimental and theoretical study of helium adsorption on cationic hexaphenylbenzene (HPB), a propeller-shaped molecule. The orientation of its propeller blades is known to be sensitive to the environment, with substantial differences between the molecule in the gas phase and in the crystalline solid. Mass spectra of He_{n}HPB^{+}, synthesized in helium nanodroplets, indicate enhanced stability for ions containing n = 2, 4, 14, 28, 42, 44, or 46 helium atoms. Path-integral molecular dynamics simulations reveal a significant dependence of the dissociation energy on the details of the HPB geometry. Good agreement between the experimental data and calculated dissociation energies is obtained, provided that the symmetry of HPB^{+} is reduced from D_{6} to D_{2}, such a lower symmetry being suggested from quantum chemical calculations as arising upon electron removal.Graphic

Adsorption of Helium on Small Cationic PAHs: Influence of Hydrocarbon Structure on the Microsolvation Pattern.

The adsorption of up to ∼100 helium atoms on cations of the planar polycyclic aromatic hydrocarbons (PAHs) anthracene, phenanthrene, fluoranthene, and pyrene was studied by combining helium nanodroplet mass spectrometry with classical and quantum computational methods. Recorded time-of-flight mass spectra reveal a unique set of structural features in the ion abundance as a function of the number of attached helium atoms for each of the investigated PAHs. Path-integral molecular dynamics simulations were used with a polarizable potential to determine the underlying adsorption patterns of helium around the studied PAH cations and in good general agreement with the experimental data. The calculated structures of the helium–PAH complexes indicate that the arrangement of adsorbed helium atoms is highly sensitive toward the structure of the solvated PAH cation. Closures of the first solvation shell around the studied PAH cations are suggested to lie between 29 and 37 adsorbed helium atoms depending on the specific PAH cation. Helium atoms are found to preferentially adsorb on these PAHs following the commensurate pattern common for graphitic surfaces, in contrast to larger carbonaceous molecules like corannulene, coronene, and fullerenes that exhibit a 1 × 1 commensurate phase.

Quantum effects in the structural and elastic properties of graphite: Path-integral simulations

Graphite, as a well-known carbon-based solid, is a paradigmatic example of the so-called van der Waals layered materials, which display a large anisotropy in their physical properties. Here we study quantum effects in structural and elastic properties of graphite by using path-integral molecular dynamics simulations in the temperature range from 50 to 1500~K. This method takes into account quantization and anharmonicity of vibrational modes in the material. Our results are compared with those found by using classical molecular dynamics simulations. We analyze the volume and in-plane area as functions of temperature and external stress. The quantum motion is essential to correctly describe the in-plane and out-of-plane thermal expansion. Quantum effects cause also changes in the elastic properties of graphite with respect to a classical model. At low temperature we find an appreciable decrease in the linear elastic constants, mainly in $C_{12}$ and $C_{44}$. Quantum corrections in stiffness constants can be in some cases even larger than 20\%. The bulk modulus and Poisson's ratio are reduced in a 4\% and 19\%, respectively, due to zero-point motion of the C atoms. These quantum effects in structural and elastic properties of graphite are nonnegligible up to temperatures higher than 300~K.