Path Integral Molecular Dynamics Simulations Research Articles

Relationship between Hydrogen-Bonding Motifs and the 1b1 Splitting in the X-ray Emission Spectrum of Liquid Water.

The split of the 1b1 peak observed in the X-ray emission (XE) spectrum of liquid water has been the focus of intense research. Although several hypotheses have been proposed to explain the origin of this split, a consensus has not yet been reached. Here, we introduce a novel theoretical/computation approach which, combining path-integral molecular dynamics simulations with the MB-pol model and time-dependent density functional theory calculations, predicts the 1b1 splitting in liquid water and not in crystalline ice, in agreement with the experimental observations. A systematic analysis of the underlying local structure of liquid water at ambient conditions indicates that several different hydrogen-bonding motifs contribute to the overall XE line shape in the energy range corresponding to emissions from the 1b1 orbitals. This suggests that it is not possible to unambiguously attribute the split of the 1b1 peak to only two specific structural arrangements of the underlying hydrogen-bonding network.

Nuclear quantum effects on autoionization of water isotopologs studied by ab initio path integral molecular dynamics.

In this study, we investigate the nuclear quantum effects (NQEs) on the acidity constant (pKA) of liquid water isotopologs under the ambient condition by path integral molecular dynamics (PIMD) simulations. We compared simulations using a fully explicit solvent model with a classical polarizable force field, density functional tight binding, and ab initio density functional theory, which correspond to empirical, semiempirical, and ab initio PIMD simulations, respectively. The centroid variable with respect to the proton coordination number of a water molecule was restrained to compute the gradient of the free energy, which measures the reversible work of the proton abstraction for the quantum mechanical system. The free energy curve obtained by thermodynamic integration was used to compute the pKA value based on probabilistic determination. This technique not only reproduces the pKA value of liquid D2O experimentally measured (14.86) but also allows for a theoretical prediction of the pKA values of liquid T2O and aqueous HDO and HTO, which are unknown due to their scarcity. It is also shown that the NQEs on the free energy curve can result in a downshift of 4.5 ± 0.9 pKA units in the case of liquid water, which indicates that the NQEs plays an indispensable role in the absolute determination of pKA. The results of this study can help inform further extensions into the calculation of the acidity constants of isotope substituted species with high accuracy.

The metastable structures of anthracene-argon clusters inside helium nanodroplets

Recent experiments have revealed that anthracene forms rather diverse complexes with multiple argon atoms upon embedding in helium nanodroplets. In this work we have modeled such complexes and investigated the relative stability of various isomeric forms using classical force fields reparametrized on electronic structure data at the level of spin-component scaled MP2, nuclear quantum effects being accounted for by means of path-integral molecular dynamics simulations. We generally find that the structures obtained are very sensitive to the details of the potential, but also to the inclusion of zero-point motion, for which the harmonic approximation is found to be generally accurate. Inside He $$_{1000}$$ clusters, the complexes also exhibit different relative stabilities, the quantum solvent usually favoring more asymmetric complexes. No evidence for fluxionality at 1 K in a given cluster is found. Using ring-polymer molecular dynamics trajectories, the stability of an argon dimer initially distant from the anthracene molecule was also investigated. Such metastable configurations are found to have a lifetime of the order of nanoseconds, confirming the possibility of trapping them in experiments.

Quantum Nuclei at Weakly Bonded Interfaces: The Case of Cyclohexane on Rh(111)

AbstractThe electronic properties of interfaces can depend on their isotopic constitution. One known case is that of cyclohexane physisorbed on Rh(111), in which isotope effects have been measured on the work function change and desorption energies. These effects can only be captured by calculations including nuclear quantum effects (NQE). In this paper, this interface is addressed employing dispersion‐inclusive density‐functional theory coupled to a quasi‐harmonic (QH) approximation for NQE, as well as to fully anharmonic ab initio path integral molecular dynamics (PIMD). The QH approximation is able to capture that deuterated cyclohexane has a smaller adsorption energy and lies about 0.01 Å farther from the Rh(111) surface than its isotopologue, which can be correlated to the isotope effect in the work function change. An investigation of the validity of the QH approximation relying on PIMD simulations, leads to the conclusion that although this interface is highly impacted by anharmonic quantum fluctuations in the molecular layer and at bonding sites, these anharmonic contributions play a minor role when analyzing isotope effects at low temperatures. Nevertheless, anharmonic quantum fluctuations cause an increase in the distance between the molecular layer and Rh(111), a consequent smaller overall work function change, and intricate changes in orbital hybridization.

First-principles equation of state database for warm dense matter computation.

We put together a first-principles equation of state (FPEOS) database for matter at extreme conditions by combining results from path integral Monte Carlo and density functional molecular dynamics simulations of the elements H, He, B, C, N, O, Ne, Na, Mg, Al, and Si as well as the compounds LiF,B_{4}C,BN,CH_{4},CH_{2},C_{2}H_{3},CH,C_{2}H,MgO,andMgSiO_{3}. For all these materials, we provide the pressure and internal energy over a density-temperature range from ∼0.5 to 50 gcm^{-3} and from ∼10^{4} to 10^{9}K, which are based on ∼5000 different first-principles simulations. We compute isobars, adiabats, and shock Hugoniot curves in the regime of L- and K-shell ionization. Invoking the linear mixing approximation, we study the properties of mixtures at high density and temperature. We derive the Hugoniot curves for water and alumina as well as for carbon-oxygen, helium-neon, and CH-silicon mixtures. We predict the maximal shock compression ratios of H_{2}O, H_{2}O_{2}, Al_{2}O_{3},CO,andCO_{2} to be 4.61, 4.64, 4.64, 4.89, and 4.83, respectively. Finally we use the FPEOS database to determine the points of maximum shock compression for all available binary mixtures. We identify mixtures that reach higher shock compression ratios than their end members. We discuss trends common to all mixtures in pressure-temperature and particle-shock velocity spaces. In the Supplemental Material, we provide all FPEOS tables as well as computer codes for interpolation, Hugoniot calculations, and plots of various thermodynamic functions.

Nuclear quantum effects: Their relevance in neutron diffraction studies of liquid water

Corrections for nuclear quantum effects (NQE) have been calculated for classical molecular dynamics (MD) simulation models of light (H2O), heavy (D2O) and ‘null’ [(H2O)0.64(D2O)0.36] water. New path integral molecular dynamics (PIMD) simulations have also been conducted for the same systems. NQEs have somewhat smaller influence on the O-D and D-D partial radial distribution functions of heavy water than on the O-H and H-H ones of light water. After correcting for NQEs the O-‘H’ bondlengths in light and heavy water have become different: the O-D ones are about 0.5% shorter than the O-H ones. Following NQE corrections, the total RDF of ‘null’ water does show hydrogen related features at the position of the intramolecular O-‘H’ peak. Based on a cross-check procedure involving the NQE-corrected total and partial radial distribution functions, it can be stated that concerning the structure of liquid water, the assumption that H and D are equal is valid to a very good approximation for intermolecular correlations. These findings are also supported by our path integral molecular dynamics simulations.

Two-temperature warm dense hydrogen as a test of quantum protons driven by orbital-free density functional theory electronic forces

We consider a steady-state (but transient) situation in which a warm dense aggregate is a two-temperature system with equilibrium electrons at temperature Te, ions at Ti, and Te ≠ Ti. Such states are achievable by pump–probe experiments. For warm dense hydrogen in such a two-temperature situation, we investigate nuclear quantum effects (NQEs) on structure and thermodynamic properties, thereby delineating the limitations of ordinary ab initio molecular dynamics. We use path integral molecular dynamics (PIMD) simulations driven by orbital-free density functional theory (OFDFT) calculations with state-of-the-art noninteracting free-energy and exchange-correlation functionals for the explicit temperature dependence. We calibrate the OFDFT calculations against conventional (explicit orbitals) Kohn–Sham DFT. We find that when the ratio of the ionic thermal de Broglie wavelength to the mean interionic distance is larger than about 0.30, the ionic radial distribution function is meaningfully affected by the inclusion of NQEs. Moreover, NQEs induce a substantial increase in both the ionic and electronic pressures. This confirms the importance of NQEs for highly accurate equation-of-state data on highly driven hydrogen. For Te > 20 kK, increasing Te in the warm dense hydrogen has slight effects on the ionic radial distribution function and equation of state in the range of densities considered. In addition, we confirm that compared with thermostatted ring-polymer molecular dynamics, the primitive PIMD algorithm overestimates electronic pressures, a consequence of the overly localized ionic description from the primitive scheme.

Anomalous Freezing of Low-Dimensional Water Confined in Graphene Nanowrinkles.

Various properties of water are affected by confinement as the space-filling of the water molecules is very different from bulk water. In our study, we challenged the creation of a stable system in which water molecules are permanently locked in nanodimensional graphene traps. For that purpose, we developed a technique, nitrocellulose-assisted transfer of graphene grown by chemical vapor deposition, which enables capturing of the water molecules below an atomically thin graphene membrane structured into a net of regular wrinkles with a lateral dimension of about 4 nm. After successfully confining water molecules below a graphene monolayer, we employed cryogenic Raman spectroscopy to monitor the phase changes of the confined water as a function of the temperature. In our experiment system, the graphene monolayer structured into a net of fine wrinkles plays a dual role: (i) it enables water confinement and (ii) serves as an extremely sensitive probe for phase transitions involving water via graphene-based spectroscopic monitoring of the underlying water structure. Experimental findings were supported with classical and path integral molecular dynamics simulations carried out on our experimental system. Results of simulations show that surface premelting of the ice confined within the wrinkles starts at ∼200 K and the melting process is complete at ∼240 K, which is far below the melting temperature of bulk water ice. The processes correspond to changes in the doping and strain in the graphene tracked by Raman spectroscopy. We conclude that water can be confined between graphene structured into nanowrinkles and silica substrate and its phase transitions can be tracked via Raman spectral feature of the encapsulating graphene. Our study also demonstrated that peculiar behavior of liquids under spatial confinement can be inspected via the optical response of atomically thin graphene sensors.

Quantum effects and 1H NMR chemical shifts of a bifurcated short hydrogen bond.

The monoprotonated compound N,N',N''-tris(p-tolyl)azacalix[3](2,6)pyridine (TAPH) contains an intramolecular hydrogen bond that is formed from three N atoms in its cavity. Constrained by the macrocyclic molecular structure, the separations between the N atoms in this bifurcated hydrogen bond are about 2.6 Å, considerably shorter than those typically observed for hydrogen bonded systems in the condensed phases. As such, TAPH exhibits significantly elongated N-H lengths in its hydrogen bond and a downfield 1H NMR chemical shift of 22.1 ppm. In this work, we carry out ab initio molecular dynamics and ab initio path integral molecular dynamics simulations of TAPH in the acetonitrile solution to reveal the geometry and proton sharing conditions of the bifurcated short hydrogen bond and uncover how the interplay of electronic and nuclear quantum effects gives rise to its far downfield 1H chemical shift. Taking a linear short hydrogen bond as a reference, we demonstrate the distinct features of competing quantum effects and electronic shielding effects in the bifurcated hydrogen bond of TAPH. We further use the degree of deshielding on the proton as a measure of the hydrogen bonding interactions and evaluate the strength of the bifurcated short hydrogen bond as compared to its linear counterpart.

Towards Accurate Predictions of Proton NMR Spectroscopic Parameters in Molecular Solids.

The factors contributing to the accuracy of quantum-chemical calculations for the prediction of proton NMR chemical shifts in molecular solids are systematically investigated. Proton chemical shifts of six solid amino acids with hydrogen atoms in various bonding environments (CH, CH2 , CH3 , OH, SH and NH3 ) were determined experimentally using ultra-fast magic-angle spinning and proton-detected 2D NMR experiments. The standard DFT method commonly used for the calculations of NMR parameters of solids is shown to provide chemical shifts that deviate from experiment by up to 1.5 ppm. The effects of the computational level (hybrid DFT functional, coupled-cluster calculation, inclusion of relativistic spin-orbit coupling) are thoroughly discussed. The effect of molecular dynamics and nuclear quantum effects are investigated using path-integral molecular dynamics (PIMD) simulations. It is demonstrated that the accuracy of the calculated proton chemical shifts is significantly better when these effects are included in the calculations.

Nuclear Quantum Effect on the Geometry of NH4+(H2O)

Abstract The nuclear quantum effect of a singly hydrated ammonium ion cluster, NH4+(H2O), is investigated using an ab initio on-the-fly path integral molecular dynamics (PIMD) simulation with a BHandHLYP/6-31++G(d,p) level of calculations. Owing to the nuclear quantum effect with a PIMD simulation, the bond length of N–O (RNO) is shortened and the distribution of the angle of O–N–Hb ($\theta _{\text{ONH}_{\text{b}}}$) is greater. These results indicate that the nuclear quantum effect has two apparent contradictory contributions. One is the strengthening of the hydrogen bond arising from the zero-point energy (ZPE), whereas the other is facilitating the rotation of NH4+ in water owing to a quantum fluctuation.

Nuclear quantum effects on the hydrogen bond donor-acceptor exchange in water-water and water-methanol dimers.

We present results from path integral molecular dynamics simulations that describe effects from the explicit incorporation of nuclear quantum fluctuations on the topology of the free energy associated with the geared exchange of hydrogen bonds in the water-water dimer. Compared to the classical treatment, our results reveal important reductions in the free energy barriers and changes at a qualitative level in the overall profile. Most notable are those manifested by a plateau behavior, ascribed to nuclear tunneling, which bridges reactant and product states, contrasting with the usual symmetric double-well profile. The characteristics of the proton localizations along the pathway are examined. An imaginary time analysis of the rotational degrees of freedom of the partners in the dimer at the vicinities of transition states shows a clear "anticorrelation" between intermolecular interactions coupling beads localized in connective and dangling basins of attractions. As such, the transfer is operated by gradual concerted inter-basin migrations in opposite directions, at practically no energy costs. Modifications operated by partial deuteration and by the asymmetries in the hydrogen bonding characteristics prevailing in water-methanol heterodimers are also examined.

Isotopic effects in structural properties of graphene

Isotopic effects are relevant to understand several properties of solids, and have been thoroughly analyzed along the years. These effects may depend on the dimensionality of the considered solid. Here we assess their magnitude for structural properties of graphene, a paradigmatic two-dimensional material. We use path-integral molecular dynamics simulations, a well-suited technique to quantify the influence of nuclear quantum effects on equilibrium variables, especially in cases where anharmonic effects are important. Emphasis is put on interatomic distances and mean-square displacements, as well as on the in-plane area of the graphene layer. At low temperature, the relative difference in C–C distance for 13C and 14C, with respect to 12C, is found to be − 2.5 and − 4.6 × 10−4, respectively, larger than in three-dimensional carbon-based materials such as diamond. For the in-plane area, the relative changes amount to − 3.9 and − 6.9 × 10−4. The magnitude of anharmonicity in the lattice vibrations is estimated by comparing the internal energy and atomic vibrational amplitudes with those derived from a harmonic approximation.

Quantum driven proton diffusion in brucite-like minerals under high pressure

Transport of hydrogen in hydrous minerals under high pressure is a key step for the water cycle within the Earth interior. Brucite Mg(OH)2 is one of the simplest minerals containing hydroxyl groups and is believed to decompose under the geological condition of the deep Earth’s mantle. In the present study, we investigate the proton diffusion in brucite under high pressure, which results from a complex interplay between two processes: the O–H reorientations motion around the c axis and O–H covalent bond dissociations. First-principle path-integral molecular dynamics simulations reveal that the increasing pressure tends to lock the former motion, while, in contrast, it activates the latter which is mainly triggered by nuclear quantum effects. These two competing effects therefore give rise to a pressure sweet spot for proton diffusion within the mineral. In brucite Mg(OH)2, proton diffusion reaches a maximum for pressures close to 70GPa, while the structurally similar portlandite Ca(OH)2 never shows proton diffusion within the pressure range and time scale that we explored. We analyze the different behavior of brucite and portlandite, which might constitute two prototypes for other minerals with same structure.

Path integral molecular dynamics for fermions: Alleviating the sign problem with the Bogoliubov inequality.

We present a method for performing path integral molecular dynamics (PIMD) simulations for fermions and address its sign problem. PIMD simulations are widely used for studying many-body quantum systems at thermal equilibrium. However, they assume that the particles are distinguishable and neglect bosonic and fermionic exchange effects. Interacting fermions play a key role in many chemical and physical systems, such as electrons in quantum dots and ultracold trapped atoms. A direct sampling of the fermionic partition function is impossible using PIMD since its integrand is not positive definite. We show that PIMD simulations for fermions are feasible by employing our recently developed method for bosonic PIMD and reweighting the results to obtain fermionic expectation values. The approach is tested against path integral Monte Carlo (PIMC) simulations for up to seven electrons in a two-dimensional quantum dot for a range of interaction strengths. However, like PIMC, the method suffers from the sign problem at low temperatures. We propose a simple approach for alleviating it by simulating an auxiliary system with a larger average sign and obtaining an upper bound to the energy of the original system using the Bogoliubov inequality. This allows fermions to be studied at temperatures lower than would otherwise have been feasible using PIMD, as demonstrated in the case of a three-electron quantum dot. Our results extend the boundaries of PIMD simulations of fermions and will hopefully stimulate the development of new approaches for tackling the sign problem.

Nuclear quantum effect for hydrogen adsorption on Pt(111)

Nuclear quantum effect and many-body interaction importantly interplay in the hydrogen on the Pt(111) system under the high coverage conditions of electrochemical interest, as revealed by our ab initio path integral and ring polymer molecular dynamics simulations done at room temperature. At the full monolayer coverage, hydrogen atoms are close-packed either at the atop sites or the fcc sites owing to their strong repulsion and the nearly degenerate nature of the adsorption sites. While at the 2/3 monolayer, they are delocalized over the fcc and hcp sites via the bridge sites because of the hopping. The quantum many-body effect is thus crucially important in determining the coverage dependence and provides a clue for reconciling the long-standing controversy on this system.

Superconducting transition temperatures of metallic liquids

We develop a nonperturbative approach for calculating the superconducting transition temperatures (Tc's) of liquids. The electron-electron scattering amplitude induced by electron-phonon coupling (EPC), from which an effective pairing interaction can be inferred, is related to the fluctuation of the T matrix of electron scattering induced by ions. By applying the relation, EPC parameters can be extracted from a path-integral molecular dynamics simulation. For determining Tc, the linearized Eliashberg equations are reestablished nonperturbatively. We apply the approach to estimate Tc's of metallic hydrogen liquids. It indicates that metallic hydrogen liquids in the pressure regime from 0.5 to 1.5TPa have Tc's well above their melting temperatures and therefore are superconducting liquids.1 MoreReceived 29 January 2019Revised 28 January 2020Accepted 30 January 2020DOI:https://doi.org/10.1103/PhysRevResearch.2.013340Published by the American Physical Society under the terms of the Creative Commons Attribution 4.0 International license. Further distribution of this work must maintain attribution to the author(s) and the published article's title, journal citation, and DOI.Published by the American Physical SocietyPhysics Subject Headings (PhySH)Research AreasElectron-phonon couplingMethods in superconductivityPressure effectsSuperconducting phase transitionTransition temperaturePhysical SystemsHigh-temperature superconductorsLiquidsUnconventional superconductorsTechniquesDensity functional theoryEliashberg theoryGreen's function methodsCondensed Matter, Materials & Applied Physics

Vibrationally Mediated Stabilization of Electrons in Nonpolar Matter.

We explore solvation of electrons in nonpolar matter, here represented by butadiene clusters. Isolated butadiene supports only the existence of transient anions (resonances). Two-dimensional electron energy loss spectroscopy shows that the resonances lead to an efficient vibrational excitation of butadiene, which can result into the almost complete loss of energy of the interacting electron. Cluster-beam experiments show that molecular clusters of butadiene form stable anions, however only at sizes of more than 9 molecular units. We have calculated the distribution of electron affinities of clusters using classical and path integral molecular dynamics simulations. There is almost a continuous transition from the resonant to the bound anions with an increase in cluster size. The comparison of the classical and quantum dynamics reveals that the electron binding is strongly supported by molecular vibrations, brought about by nuclear zero-point motion and thermal agitation. We also inspected the structure of the solvated electron, finding it well localized.

Nuclear quantum effects in graphane

Graphane is a quasi-two-dimensional material consisting of a single layer of fully hydrogenated graphene, with a C:H ratio of 1. We study nuclear quantum effects in the so-called chair-graphane by using path-integral molecular dynamics (PIMD) simulations. The interatomic interactions are modeled by a tight-binding potential model fitted to density-functional calculations. Finite-temperature properties are studied in the range from 50 to 1500 K. To assess the magnitude of nuclear quantum effects in the properties of graphane, classical molecular dynamics simulations have been also performed. These quantum effects are significant in structural properties such as interatomic distances and layer area at finite temperatures. The in-plane compressibility of graphane is found to be about twice larger than that of graphene, and at low temperature it is 9% higher than the classical calculation. The thermal expansion coefficient resulting from PIMD simulations vanishes in the zero-temperature limit, in agreement with the third law of Thermodynamics.

Combined Effects from Solvation and Nuclear Quantum Fluctuations on Autoionization Mechanisms in Aqueous Clusters.

Using path-integral molecular dynamics simulations, we examine isomerization paths involving collective proton transfers in [H2O]5 and [H2O]8 clusters under cryogenic conditions. We focused attention on combined effects derived from solvation and nuclear quantum fluctuations on the characteristics of free energy barriers and relative stabilities of reactants and products. In particular, we analyzed two different processes: the first one involves the exchange of donor-acceptor hydrogen bond roles along cyclic moieties, whereas the second one corresponds to charge separation leading to stable [H3O]+[OH]- ion pairs. In the first case, the explicit incorporation of quantum tunneling introduces important modifications in the classical free energy profile. The resulting quantum profile presents two main contributions, one corresponding to compressions of O-O distances and a second one ascribed to nuclear tunneling of the light protons. Solvation effects promote a moderate polarization of the cyclic structures and a partial loss of concertedness in the collective modes, most notably, at the onset of tunneling. Still, the latter effects are also sufficiently strong to promote the stabilization of ion pairs along the classical trajectories. In contrast, the explicit incorporation of nuclear quantum fluctuations leads to charge separated configurations that are marginally stable. As such, the latter states could also be regarded as short-lived intermediate states along the reactive exchange path.