Hydrogen and muonium in diamond: A path-integral molecular dynamics simulation

Hydrogen bonded electrolytes that exhibit accelerated proton transport via sequential reactive hops have drawn interest for their promise in clean energy applications. Molecular dynamics simulations of these electrolytes offer the opportunity to uncover microscopic mechanistic details that could be used to design and tune the properties of candidate electrolyte technologies. However, accurately modeling the proton transfer reactions and transport properties that give rise to high charge conductivites in these electrolytes proves computationally challenging because of the need to perform lengthy condensed phase simulations, treating both the electronic and nuclear degrees of freedom quantum mechanically. In this paper, we demonstrate that such a modeling task can be efficiently achieved with the use of density functional theory (DFT)-trained machine learning potentials (MLP) to accelerate path integral molecular dynamics (PIMD) simulations. We highlight the practical utility of this approach by using it to benchmark how closely PIMD simulations employing different DFT exchange-correlation functionals reproduce the composition-dependent densities, diffusion coefficients, and electrical conductivities of mixtures consisting of imidazole and levulinic acid. Even with the speedup afforded by our MLPs, PIMD simulations remain quite expensive. In order to render PIMD more computationally tractable, we introduce and benchmark the accuracy of a ring polymer contraction approach that leverages a computationally efficient short-range MLP to accelerate our PIMD simulations by an additional factor of four.

Close encounters: The geometric isotope effect (GIE) for the hydrogen bonds in positively and negatively charged ammonia dimers is studied by ab initio path integral MD simulation. The hydrogen-bonded H atom in the N2H7+ cation is located around the center of two N atoms, while the hydrogen-bonded D atom in the N2D7+ cation is closer to one of the N atoms (see figure).

Starting with separated atoms on the surface of helium clusters 4HeN, and as a first step to assess the formation of rubidium dimers Rb2 in the triplet state, we perform Path Integral Molecular Dynamics simulations in the NVT canonical ensemble. Based on an accurate potential energy surface (PES) for the He–Rb2(3Σ u +) interaction [Guillon et al., J. Chem. Phys. 136, 174307 (2012)], the total PES is analytically described as the addition of pair interactions. The i-PI code [Ceriotti et al., Comput. Phys. Commun. 185, 1019 (2014)] was used to perform the simulations. At a temperature of 2 K, clusters containing up to N = 70 helium atoms, with a number up to 200 beads per particle to describe quantum effects, were considered.

The QM/MM BEST method presented in the first article of this series [M. Shiga and M. Masia, J. Chem. Phys. 139, 044120 (2013)] has been applied herein to simulate the whole series of hydrated alkali ions. In this article we show how to overcome the sampling bottleneck for QM/MM simulations by using our method with multiple time scale algorithm (MTS-BEST). We extend the use of MTS-BEST to ab initio QM/MM path integral molecular dynamics simulations, thus demonstrating that one could obtain a complete quantum description of the primary subsystem based on first principles. We highlight that the MTS-BEST approach could be generally applied to hybrid multiscale simulation of diffusive systems, thus extending its relevance to a broad class of simulation techniques beyond QM/MM. We show that it is important to account for electron correlation to better reproduce the hydration structural properties such as the ion-water radial distribution functions, and the anisotropic angular distributions around the ion.

Ab initio path integral simulation study on 16O/ 18O isotope effect in water and hydronium ion

A local pseudopotential for the interaction of an electron with a water molecule in the electronic ground state is developed. The potential contains Coulomb, adiabatic polarization, exclusion, and exchange contributions. The potential is suitable for a description of excess electron states in water clusers, and for studies of electron solvation in water. Quantum path integral molecular dynamics simulations of electron localization in water clusters employing this potential yield results in agreement with available experimental data and all-electron quantum chemical calculations.

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.

We describe how the fast quasi-centroid molecular dynamics (f-QCMD) method can be applied to condensed-phase systems by approximating the quasi-centroid potential of mean force as a sum of inter- and intramolecular corrections to the classical interaction potential. The corrections are found by using a regularized iterative Boltzmann inversion procedure to recover the inter- and intramolecular quasi-centroid distribution functions obtained from a path integral molecular dynamics simulation. The resulting methodology is found to give good agreement with a previously published QCMD dipole absorption spectrum for liquid water and satisfactory agreement for ice. It also gives good agreement with spectra from a recent implementation of CMD that uses a precomputed elevated temperature potential of mean force. Modern centroid molecular dynamics methods, therefore, appear to be reaching a consensus regarding the impact of nuclear quantum effects on the vibrational spectra of water and ice.

The dynamics of the K atom attached to the 300-atom helium cluster upon photoexcitation has been studied theoretically based on the quantum molecular dynamics method. With the use of a hybrid method between time-dependent quantum dynamics and semiclassical path integral centroid molecular dynamics, the calculation includes the electronic state mixing for the K(2P) states as well as the quantum motions of helium atoms. It was found from the calculated trajectories that the K atom mostly desorbs by itself from the surface of the helium cluster without depending on the types of initial electronic excited states. It was also found that the K*He exciplex can be formed only as a minor product, following the initial excitation to the Π-type electronic state alone. Helium solvation structures of the K*Hen exciplexes were also investigated by path integral molecular dynamics simulations. The result implies that stable exciplexes can be formed up to n = 6 on the lowest excited state until the first solvation shell around the 4p orbital of K* is completed. Meanwhile, it was found that only two helium atoms can be bound to the K(2P) atom on the second lowest excited state.

Nuclear quantum effect on intramolecular hydrogen bond of hydrogen maleate anion: An ab initio path integral molecular dynamics study

We carried out Path Integral Molecular Dynamics simulations that describe the microscopic properties of two isomerization processes taking place in the MeOH(H2O)2H+ trimer, at T = 50 K. In particular, we focused attention on the free energies associated with the exchange of the connective/dangling characteristics of a pair of protons located at key subunits of the trimer. In one of the processes, the isomerization produces a modification in the topology of the cluster's connectivity pattern, from branched-like to chain-like motifs. In contrast, along the other transformation, reactant and product states are both chain-like and equivalent. Changes in associated free energies were computed following reversible paths described in terms of order parameters involving angular degrees of freedom. As a common feature, along both isomerizations, we registered concomitant migrations of the excess proton. The strongly confining, single-well characteristics of the potential energy surface along the asymmetric stretch coordinate promote compact, ∼0.35 Å long structures for the isomorphic polymer associated with the itinerant proton, which persist along the complete isomerization path. These observations suggest the absence of tunneling contributions to the resulting mechanisms that control the proton transfer process. Estimates for the corresponding isomerization rates are also computed.

The dynamic, quantum structure of [Mg(H2)(n=1-6)](2+)complexes is investigated via ab initio path integral molecular dynamics simulations. These complexes represent the strong, σ-complex regime of metal-H2 interactions and are representative of bonding motifs found in metal-organic frameworks. Significant nuclear motion within the coordination sphere is observed, even though the ligands remain largely intact. Quantum effects are found to be important in the H-H and metal-H2 stretch coordinates, but the remaining motion in the molecule is well represented by classical simulations. Nearly free rotation of the dihydrogen moiety is observed in all complexes. Statistical averages and distributions of structural parameters are found to deviate nontrivially from the same parameters in static, equilibrium structures.

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

The structural properties and thermal transport of superionic ice are critical to understanding the ice giant planets. However, their nuclear quantum effects (NQEs) remain unknown. Here, we systematically investigated the NQEs of body-centered cubic superionic ice based on a deep potential model combined with path-integral molecular dynamics simulations. The NQEs are found to significantly modify the structural properties of superionic ice. The quantum delocalization of protons leads to a contraction of the oxygen lattice. The radial distribution function analysis reveals that NQEs cause a leftward shift of the O-O peak, elongation of the O-H bond, merger of the H-H peaks, and anomalous temperature independence of the first O-H and H-H peaks. The total thermal conductivity obtained from quantum simulations is significantly higher than that from classical simulations, with the primary contribution arising from the enhanced heat convection due to proton diffusion rather than from heat conduction dominated by lattice vibrations. This study elucidates the microscopic mechanism of NQE-regulated thermal transport in superionic ice through the enhancement of proton convection, providing an important quantum-correction basis for modeling the thermal physical properties of ice under extreme conditions.

Thermal properties of graphene monolayers are studied by path-integral molecular dynamics simulations, which take into account the quantization of vibrational modes in the crystalline membrane and allow one to consider anharmonic effects in these properties. This system was studied at temperatures in the range from 12 to 2000 K and zero external stress, by describing the interatomic interactions through the LCBOPII effective potential. We analyze the internal energy and specific heat and compare the results derived from the simulations with those yielded by a harmonic approximation for the vibrational modes. This approximation turns out to be rather precise up to temperatures of about 400 K. At higher temperatures, we observe an influence of the elastic energy due to the thermal expansion of the graphene sheet. Zero-point and thermal effects on the in-plane and "real" surface of graphene are discussed. The thermal expansion coefficient α of the real area is found to be positive at all temperatures, in contrast to the expansion coefficient αp of the in-plane area, which is negative at low temperatures and becomes positive for T ≳ 1000 K.