On the short-time limit of ring polymer molecular dynamics

We examine the short-time accuracy of a class of approximate quantum dynamical techniques that includes the centroid molecular dynamics (CMD) and ring polymer molecular dynamics (RPMD) methods. Both of these methods are based on the path integral molecular dynamics (PIMD) technique for calculating the exact static equilibrium properties of quantum mechanical systems. For Kubo-transformed real-time correlation functions involving operators that are linear functions of positions or momenta, the RPMD and (adiabatic) CMD approximations differ only in the choice of the artificial mass matrix of the system of ring polymer beads that is employed in PIMD. The obvious ansatz for a general method of this type is therefore to regard the elements of the PIMD (or Parrinello-Rahman) mass matrix as an adjustable set of parameters that can be chosen to improve the accuracy of the resulting approximation. We show here that this ansatz leads uniquely to the RPMD approximation when the criterion that is used to select the mass matrix is the short-time accuracy of the Kubo-transformed correlation function. In particular, we show that the leading error in the RPMD position autocorrelation function is O(t(8)) and the error in the velocity autocorrelation function is O(t(6)), for a general anharmonic potential. The corresponding errors in the CMD approximation are O(t(6)) and O(t(4)), respectively.     

A comparative study of imaginary time path integral based methods for quantum dynamics

The recently introduced approximate many-body quantum simulation method, ring polymer molecular dynamics (RPMD), is compared to the centroid molecular dynamics method (CMD). Comparisons of simulation results for liquid para-hydrogen at two state points and liquid ortho-deuterium at one state point are presented. The calculated quantum correlation functions for the two methods are shown to be in good agreement with one another for a large portion of the time spectrum. However, as the quantum mechanical nature of the system increases, RPMD is less accurate in predicting the kinetic energy of the system than is CMD. A simplified and highly efficient algorithm is proposed which largely corrects this deficiency.     

Comparison of approximate quantum simulation methods applied to normal liquid helium at 4 K

The Feynman-Kleinert linearized path integral molecular dynamics (FK-LPI), ring polymer molecular dynamics (RPMD), and centroid molecular dynamics (CMD) methods are applied to the simulation of normal liquid helium. Comparisons of the simulation results at the T = 4 K and rho = 0.01873 A-3 state point are presented. The calculated quantum correlation functions for the three methods show significant differences, both in the short time and in the intermediate regions of the spectrum. Our simulation results are also compared to the recent results of other approximate quantum simulation methods. We find that FK-LPI qualitatively agrees with other approximate quantum simulation results while CMD and RPMD predict a qualitatively different impulsive rebound in the velocity autocorrelation function. Frequency space analysis reveals that RPMD exhibits a broad high-frequency tail similar to that from quantum mode coupling theory and numerical analytic continuation approaches, while FK-LPI provides a somewhat more rapid decay at high frequency than any of these three methods. CMD manifests a high-frequency component that is greatly reduced compared with the other methods.     

Quantum mechanical correlation functions, maximum entropy analytic continuation, and ring polymer molecular dynamics

The maximum entropy analytic continuation (MEAC) and ring polymer molecular dynamics (RPMD) methods provide complementary approaches to the calculation of real time quantum correlation functions. RPMD becomes exact in the high temperature limit, where the thermal time betavariant Planck's over 2pi tends to zero and the ring polymer collapses to a single classical bead. MEAC becomes most reliable at low temperatures, where betavariant Planck's over 2pi exceeds the correlation time of interest and the numerical imaginary time correlation function contains essentially all of the information that is needed to recover the real time dynamics. We show here that this situation can be exploited by combining the two methods to give an improved approximation that is better than either of its parts. In particular, the MEAC method provides an ideal way to impose exact moment (or sum rule) constraints on a prior RPMD spectrum. The resulting scheme is shown to provide a practical solution to the "nonlinear operator problem" of RPMD, and to give good agreement with recent exact results for the short-time velocity autocorrelation function of liquid parahydrogen. Moreover these improvements are obtained with little extra effort, because the imaginary time correlation function that is used in the MEAC procedure can be computed at the same time as the RPMD approximation to the real time correlation function. However, there are still some problems involving long-time dynamics for which the RPMD+MEAC combination is inadequate, as we illustrate with an example application to the collective density fluctuations in liquid orthodeuterium.     

Thermal Gaussian molecular dynamics for quantum dynamics simulations of many-body systems: application to liquid para-hydrogen

A new method, here called thermal Gaussian molecular dynamics (TGMD), for simulating the dynamics of quantum many-body systems has recently been introduced [I. Georgescu and V. A. Mandelshtam, Phys. Rev. B 82, 094305 (2010)]. As in the centroid molecular dynamics (CMD), in TGMD the N-body quantum system is mapped to an N-body classical system. The associated both effective Hamiltonian and effective force are computed within the variational Gaussian wave-packet approximation. The TGMD is exact for the high-temperature limit, accurate for short times, and preserves the quantum canonical distribution. For a harmonic potential and any form of operator A?, it provides exact time correlation functions C(AB)(t) at least for the case of B, a linear combination of the position, x, and momentum, p, operators. While conceptually similar to CMD and other quantum molecular dynamics approaches, the great advantage of TGMD is its computational efficiency. We introduce the many-body implementation and demonstrate it on the benchmark problem of calculating the velocity time auto-correlation function for liquid para-hydrogen, using a system of up to N = 2592 particles.     

A refined ring polymer molecular dynamics theory of chemical reaction rates

We further develop the ring polymer molecular dynamics (RPMD) method for calculating chemical reaction rates [I. R. Craig and D. E. Manolopoulos, J. Chem. Phys. 122, 084106 (2005)]. We begin by showing how the rate coefficient we obtained before can be calculated in a more efficient way by considering the side functions of the ring-polymer centroids, rather than averaging over the side functions of the individual ring-polymer beads. This has two distinct advantages. First, the statistics of the phase-space average over the ring-polymer coordinates and momenta are greatly improved. Second, the resulting flux-side correlation function converges to its long-time limit much more rapidly. Indeed the short-time limit of this flux-side correlation function already provides a "quantum transition state theory" approximation to the final rate coefficient. In cases where transition state recrossing effects are negligible, and the transition state dividing surface is put in the right place, the RPMD rate is therefore obtained almost instantly. We then go on to show that the long-time limit of the new flux-side correlation function, and hence the fully converged RPMD reaction rate, is rigorously independent of the choice of the transition state dividing surface. This is especially significant because the optimum dividing surface can often be very difficult to determine for reactions in complex chemical systems.     

Ab initio centroid path integral molecular dynamics: application to vibrational dynamics of diatomic molecular systems

An ab initio centroid molecular dynamics (CMD) method is developed by combining the CMD method with the ab initio molecular orbital method. The ab initio CMD method is applied to vibrational dynamics of diatomic molecules, H2 and HF. For the H2 molecule, the temperature dependence of the peak frequency of the vibrational spectral density is investigated. The results are compared with those obtained by the ab initio classical molecular dynamics method and exact quantum mechanical treatment. It is shown that the vibrational frequency obtained from the ab initio CMD approaches the exact first excitation frequency as the temperature lowers. For the HF molecule, the position autocorrelation function is also analyzed in detail. The present CMD method is shown to well reproduce the exact quantum result for the information on the vibrational properties of the system.     

A centroid molecular dynamics study of liquid para-hydrogen and ortho-deuterium

Centroid molecular dynamics (CMD) is applied to the study of collective and single-particle dynamics in liquid para-hydrogen at two state points and liquid ortho-deuterium at one state point. The CMD results are compared with the results of classical molecular dynamics, quantum mode coupling theory, a maximum entropy analytic continuation approach, pair-product forward- backward semiclassical dynamics, and available experimental results. The self-diffusion constants are in excellent agreement with the experimental measurements for all systems studied. Furthermore, it is shown that the method is able to adequately describe both the single-particle and collective dynamics of quantum liquids.     

Quantum statistics and classical mechanics: real time correlation functions from ring polymer molecular dynamics

We propose an approximate method for calculating Kubo-transformed real-time correlation functions involving position-dependent operators, based on path integral (Parrinello-Rahman) molecular dynamics. The method gives the exact quantum mechanical correlation function at time zero, exactly satisfies the quantum mechanical detailed balance condition, and for correlation functions of the form C(Ax)(t) and C(xB)(t) it gives the exact result for a harmonic potential. It also works reasonably well at short times for more general potentials and correlation functions, as we illustrate with some example calculations. The method provides a consistent improvement over purely classical molecular dynamics that is most apparent in the low-temperature regime.     

Inclusion of inversion symmetry in centroid molecular dynamics: a possible avenue to recover quantum coherence

Inversion symmetry is included in the operator formulation of the centroid molecular dynamics (CMD). This work involves the development of a symmetry-adapted CMD (SA-CMD), here particularly for symmetrization and antisymmetrization projections. A symmetry-adapted quasidensity operator, as defined by Blinov and Roy [J. Chem. Phys. 115, 7822 (2001)], is employed to obtain the centroid representation of quantum mechanical operators. Numerical examples are given for a single particle confined to one-dimensional symmetric quartic and symmetric double-well potentials. Two SA-CMD simulations are performed separately for both projections, and centroid position autocorrelation functions are obtained. For each projection, the quality of the approximation as well as the accuracy are similar to those of regular CMD. It is shown that individual trajectories from two separate SA-CMD simulations can be properly combined to recover trajectories for Boltzmann statistics. Position autocorrelation functions are compared to the exact quantum mechanical ones. This explicit account of inversion symmetry provides a qualitative improvement on the conventional CMD approach and allows the recovery of some quantum coherence.     

Model for dynamics of inhomogeneous and bulk fluids

An accurate model for the density of states (DOS) for strongly inhomogeneous and bulk fluids has been proposed based on gamma distributions. The contribution to the density of states from the collective dynamics is modeled as an incomplete gamma distribution and the high frequency region is obtained from the solution of the memory equation using a sech memory kernel. Using only the frequency moments as input, the model parameters for the collective dynamics are obtained by matching moments of the resulting distribution. The model results in an analytical expression for the self-diffusivity of the fluid. We present results for soft sphere fluids confined in slit-shaped pores as well as bulk soft sphere liquids. Comparisons of the DOS, velocity autocorrelation functions, and memory kernels with molecular dynamics simulations reveal that the model predicts features in the DOS over the entire frequency range and is able to capture changes in the DOS as a function of fluid density and temperature. As a result the predicted VACFs, memory kernels, and self-diffusivities are accurately predicted over a wide range of conditions. Since the frequency moments for bulk liquids can be obtained from pair correlation functions, our method provides a direct route from fluid structure to dynamics. For fluids confined in slit-shaped pores, where the frequency moments are obtained from molecular dynamics simulations, the predicted self-diffusivities capture the resulting oscillations due to variations in the solvation pressure, and in the case of smooth walled pores, the predictions are superior to those obtained using kinetic theory.     

Fast centroid molecular dynamics: a force-matching approach for the predetermination of the effective centroid forces

A fast centroid molecular dynamics (CMD) methodology is proposed in which the effective centroid forces are predetermined through a force-matching algorithm applied to a standard path integral molecular dynamics simulation. The resulting method greatly reduces the computational cost of generating centroid trajectories, thus extending the applicability of CMD. The method is applied to the study of liquid para-hydrogen at two state points and liquid ortho-deuterium at one state point. The static and dynamical results are compared to those obtained from full adiabatic CMD simulations and found to be in excellent agreement for all three systems; the transport properties are also compared to experiment and found to have a similar level of agreement.     

Quantum effects in light and heavy liquid water: A rigid-body centroid molecular dynamics study

The centroid molecular dynamics (CMD) method is applied to the study of liquid water in the context of the rigid-body approximation. This rigid-body CMD technique, which is significantly more efficient than the standard CMD method, is implemented on the TIP4P model for water and used to examine isotopic effects in the equilibrium and dynamical properties of liquid H(2)O and D(2)O. The results obtained with this approach compare remarkably well with those determined previously with path integrals simulations as well as those obtained from the standard CMD method employing flexible models. In addition, an examination of the impact of quantization on the rotational and librational motion of the water molecule is also reported.     

An efficient and accurate implementation of centroid molecular dynamics using a Gaussian approximation

An approximate method for Centroid Molecular Dynamics (CMD) is presented which uses a Gaussian approximation. The resulting method, called Gaussian CMD (GCMD), is 100-1000 times faster than CMD because it replaces explicit path-integral sampling, which is the most time-consuming part of CMD, with a Gaussian averaging, which can be done analytically. Several methods for computing the Gaussian width parameter in the GCMD approach are also presented. This new method is shown to give satisfactory results for the position correlation function in simple one-dimensional systems when CMD itself is consistent with the exact result. The GCMD and CMD results are also compared for the case of 1-dimensional systems coupled to harmonic baths, with good success. GCMD is further compared to CMD with good success for liquid para-hydrogen at two different temperatures, 14 K and 25 K, and for ortho-deuterium at 20.7 K.     

Nonadiabatic Ehrenfest molecular dynamics within the projector augmented-wave method

We derive equations for nonadiabatic Ehrenfest molecular dynamics within the projector augmented-wave (PAW) formalism. The discretization of the electrons is time-dependent as the augmentation functions depend on the positions of the nuclei. We describe the implementation of the Ehrenfest molecular dynamics equations within the real-space PAW method. We demonstrate the applicability of our method by studying the vibration of NaCl, the torsional rotation of H(2)C=NH(2)(+) in both the adiabatic and the nonadiabatic regimes, and the hydrogen bombardment of C(40)H(16).     

Dynamic mechanism of fatty acid transport across cellular membranes through FadL: molecular dynamics simulations

FadL is an important member of the family of fatty acid transport proteins within membranes. In this study, 11 conventional molecular dynamics (CMD) and 25 steered molecular dynamics (SMD) simulations were performed to investigate the dynamic mechanism of transport of long-chain fatty acids (LCFAs) across FadL. The CMD simulations addressed the intrinsically dynamic behavior of FadL. Both the CMD and SMD simulations revealed that a fatty acid molecule can move diffusively to a high-affinity site (HAS) from a low-affinity site (LAS). During this process, the swing motion of the L3 segment and the hydrophobic interaction between the fatty acid and FadL could play important roles. Furthermore, 22 of the SMD simulations revealed that fatty acids can pass through the gap between the hatch domain and the transmembrane domain (TMD) by different pathways. SMD simulations identified nine possible pathways for dodecanoic acid (DA) threading the barrel of FadL. The binding free energy profiles between DA and FadL along the MD trajectories indicate that all of the possible pathways are energetically favorable for the transport of fatty acids; however, one pathway (path VI) might be the most probable pathway for DA transport. The reasonability and reliability of this study were further demonstrated by correlating the MD simulation results with the available mutagenesis results. On the basis of the simulations, a mechanism for the full-length transport process of DA from the extracellular side to the periplasmic space mediated by FadL is proposed.     

Electron-nuclear correlations for photo-induced dynamics in molecular dimers

Ultrafast photoinduced dynamics of electronic excitation in molecular dimers is drastically affected by the dynamic reorganization of inter- and intra- molecular nuclear configuration modeled by a quantized nuclear degree of freedom. The dynamics of the electronic population and nuclear coherence is analyzed by solving the chain of coupled differential equations for population inversion, electron-vibrational correlation, etc. Intriguing results are obtained in the approximation of a small change of the nuclear equilibrium upon photoexcitation. In the limiting case of resonance between the electronic energy gap and the frequency of the nuclear mode these results are justified by comparison to the exactly solvable Jaynes-Cummings model. It is found that the photoinduced processes in the model dimer are arranged according to their time scales: (i) Fast scale of nuclear motion, (ii) intermediate scale of dynamical redistribution of electronic population between excited states as well as growth and dynamics of electron-nuclear correlation, (iii) slow scale of electronic population approach to the quasi-equilibrium distribution, decay of electron-nuclear correlation, and decrease of the amplitude of mean coordinate oscillation. The latter processes are accompanied by a noticeable growth of the nuclear coordinate dispersion associated with the overall nuclear wave packet width. The demonstrated quantum relaxation features of the photoinduced vibronic dynamics in molecular dimers are obtained by a simple method, applicable to systems with many degrees of freedom.     

Local density augmentation and dynamic properties of hydrogen-and non-hydrogen-bonded supercritical fluids: a molecular dynamics study

The local density inhomogeneities in neat supercritical fluids were investigated via canonical molecular dynamics simulations. The selected systems under investigation were the polar and hydrogen-bonded fluid methanol as well as the quadrupolar non-hydrogen-bonded carbon dioxide one. Effective local densities, local density augmentation, and enhancement factors were calculated at state points along an isotherm close to the critical temperature of each system (T(r)=1.03). The results obtained reveal strong influence of the polarity and hydrogen bonding upon the intensity of the local density augmentation. It is found that this effect is sufficiently larger in the case of the polar and associated methanol in comparison to those predicted for carbon dioxide. For both fluids the local density augmentation values are maximized in the bulk density region near 0.7rho(c), a result that is in agreement with experiment. In addition, the local density dynamics of each fluid were investigated in terms of the appropriate time correlation functions. The behavior of these functions reveals that the bulk density dependence of the local density reorganization times is very sensitive to the specific intermolecular interactions and to the size of the local region. Also, the estimated local density reorganization time as a function of bulk density of each fluid was further analyzed and successfully related to two different time-scale relaxation mechanisms. Finally, the results obtained indicate a possible relationship between the single-molecule reorientational dynamics and the local density reorganization ones.