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.     

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.     

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.     

Symmetry aspects of diatomics-in-molecules (DIM) calculations. II. Construction of symmetry-adapted diatomic fragment subspaces

A new general procedure is presented for the construction of symmetry adapted diatomic fragment sub-spaces of spin-adapted polyatomic DIM bases obtained by direct diagonalization of the polyatomic spin-square operator. The spin part of the construction is accomplished by the diagonalization in the DIM basis of the diatomic spin-square operators. The spatial part of relevant transformations is given in terms of quantum numbers of the constituent atomic functions. Simple expressions for the spatial symmetry and spin designation of the resulting symmetry-adapted diatomic fragment states are obtained. Special attention is paid to consistent phases of diatomic fragment functions corresponding to equivalent fragments. The approach eliminates such mental operations from the construction of symmetry-adapted diatomic fragment subspaces which are difficult to formulate in a way suitable for straightforward computer implementation.     

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.     

New conditions for validity of the centroid molecular dynamics and ring polymer molecular dynamics

Validity of the centroid molecular dynamics (CMD) and ring polymer molecular dynamics (RPMD) in quantum liquids is studied on an assumption that momenta of liquid particles relax fast. The projection operator method allows one to derive the generalized Langevin equation including a memory effect for the full-quantum canonical (Kubo-transformed) correlation function. Similar equations for the CMD and RPMD correlation functions can be derived too. The comparison of these equations leads to conditions under which the RPMD and CMD correlation functions agree approximately with the full-quantum canonical correlation function. The condition for the RPMD is that the memory effects of the full-quantum and RPMD equations vanish quickly with the same time constants. The CMD correlation function requires additional conditions concerning static correlation.     

Implementation of molecular symmetry in valence bond calculation

A novel strategy for the construction of many-electron symmetry-adapted wave function is proposed for ab initio valence bond (VB) calculations and is implemented for valence bond self-consistent filed (VBSCF) and breathing orbital valence bond (BOVB) methods with various orbital optimization algorithms. Symmetry-adapted VB functions are constructed by the projection operator of symmetry group. The many-electron symmetry-adapted wave function is expressed in terms of symmetry-adapted VB functions, and thus the VB calculations can be performed with the molecular symmetry restriction. Test results show that molecular symmetry reduces the computational cost of both the iteration numbers and CPU time. Furthermore, excited states with specific symmetry can be conveniently obtained in VB calculations by using symmetry-adapted VB functions.     

Path-integral centroid dynamics for general initial conditions: a nonequilibrium projection operator formulation

The formulation of path-integral centroid dynamics is extended to the quantum dynamics of density operators evolving from general initial states by means of the nonequilibrium projection operator technique. It is shown that the new formulation provides a basis for applying the method of centroid dynamics to nonequilibrium situations and that it allows the derivation of new formal relations, which can be useful in improving current equilibrium centroid dynamics methods. A simple approximation of uniform relaxation for the unprojected portion of the Liouville space propagator leads to a class of practically solvable equations of motion for the centroid variables, but with an undetermined parameter of relaxation. This new class of equations encompasses the centroid molecular-dynamics (CMD) method as a limiting case, and can be applied to both equilibrium and nonequilibrium situations. Tests for the equilibrium dynamics of one-dimensional model systems demonstrate that the new equations with appropriate choice of the relaxation parameter are comparable to the CMD method.     

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.     

Optimized Monte Carlo sampling in forward-backward semiclassical dynamics

Forward-backward semiclassical dynamics (FBSD) provides a rigorous and powerful methodology for calculating time correlation functions in condensed phase systems characterized by substantial quantum mechanical effects associated with zero-point motion, quantum dispersion, or identical particle exchange symmetries. The efficiency of these simulations arises from the use of classical trajectories to capture all dynamical information. However, full quantization of the density operator makes these calculations rather expensive compared to fully classical molecular dynamics simulations. This article discusses the convergence properties of various correlation functions and introduces an optimal Monte Carlo sampling scheme that leads to a significant reduction of statistical error. A simple and efficient procedure for normalizing the FBSD results is also discussed. Illustrative examples on model systems are presented.     

Symmetrized correlation function for liquid para-hydrogen using complex-time pair-product propagators

We present a simple and efficient method for calculating symmetrized time correlation functions of neat quantum fluids. Using the pair-product approximation to each complex-time quantum mechanical propagator, symmetrized correlation functions are written in terms of a double integral for each degree of freedom with a purely positive integrand. At moderate temperatures and densities, where the pair-product approximation to the Boltzmann operator is sufficiently accurate, the method leads to quantitative results for the early time part of the correlation function. The method is tested extensively on liquid para-hydrogen at 25 K and used to obtain accurate quantum mechanical results for the initial 0.2 ps segment of the symmetrized velocity autocorrelation function of this system, as well as the incoherent dynamic structure factor at certain momentum transfer values.     

On the calculation of single-particle time correlation functions from Bose-Einstein centroid dynamics

The calculation of single-particle time correlation functions using the Bose-Einstein centroid dynamics formalism is discussed. A new definition of the quasidensity operator is used to calculate the centroid force on a given particle for an anharmonic system. The force includes correlation effects due to quantum statistics and is used for the calculation of the classical-like dynamics of phase-space centroid variables within the centroid molecular dynamics approximation. Time correlation functions are then obtained for single-particle quantities. These correspond to the double-Kubo transform of exact quantum-mechanical correlation functions. The centroid dynamics results are compared to those of exact basis-set calculations and a good agreement is found. The level of accuracy is in fact the same as what was observed earlier for the calculation of center-of-mass correlation functions for Fermi-Dirac and Bose-Einstein statistics, and for any correlation function for Boltzmann statistics. These results show that it is now possible to use Bose-Einstein centroid molecular dynamics to calculate single-particle correlation functions for systems where quantum exchange effects are present.     

Quantum effects in ice Ih

Quantum and classical simulations are carried out on ice Ih over a range of temperatures utilizing the TIP4P water model. The rigid-body centroid molecular dynamics method employed allows for the investigation of equilibrium and dynamical properties of the quantum system. The impact of quantization on the local structure, as measured by the radial and spatial distribution functions, as well as the energy is presented. The effects of quantization on the lattice vibrations, associated with the molecular translations and librations, are also reported. Comparison of quantum and classical simulation results indicates that shifts in the average potential energy are equivalent to rising the temperature about 80 K and are therefore non-negligible. The energy shifts due to quantization and the quantum mechanical uncertainties observed in ice are smaller than the values previously reported for liquid water. Additionally, we carry out a comparative study of melting in our classical and quantum simulations and show that there are significant differences between classical and quantum ice.     

Symmetry adaption reduced to tabulated quantities

The aim of this paper is to give an algebraic formula for symmetry-adapted linear combinations avoiding intuitive or laborious projection operator techniques. By utilization of the tabulated Clebsch-Gordan coefficients and surface harmonics of the point-groups the symmetry-adapted linear combinations are given in formula (4). A five-step algorithm is proposed and the example of a tetrahedron worked out. The relation to methods using site symmetry is discussed.     

Structure-breaking effects of solvated Rb(I) in dilute aqueous solution--an ab initio QM/MM MD approach

Structural properties of the hydrated Rb(I) ion have been investigated by ab initio quantum mechanical/molecular mechanical (QM/MM) molecular dynamics (MD) simulations at the double-zeta HF quantum mechanical level. The first shell coordination number was found to be 7.1, and several other structural parameters such as angular distribution functions, radial distribution functions and tilt- and theta-angle distributions allowed the full characterization of the hydration structure of the Rb(I) ion in dilute aqueous solution. Velocity autocorrelation functions were used to calculate librational and vibrational motions, ion-ligand motions, as well as reorientation times. Different dynamical parameters such as water reorientation, mean ligand residence time, the number of ligand exchange processes, and rate constants were also analyzed. The mean ligand residence time for the first shell was determined as tau = 2.0 ps.     

Dynamical properties of the water molecules in the hydration shells of Fe(II) and Fe(III) ions: ab initio QM/MM molecular dynamics simulations

Dynamical properties, librational and vibrational motions of water molecules in the first and second hydration shells of the Fe(II) and Fe(III) ion were evaluated by means of velocity autocorrelation functions obtained by combined quantum mechanical/molecular mechanical molecular dynamics (QM/MM-MD) simulations. The frequencies of rotation around three principal axes and the frequencies of intramolecular vibration of the water molecules in the first hydration shells obtained from the simulations are blue-shifted for both ions compared to those observed experimentally for liquid water. The intramolecular geometry of water molecules in the quantum mechanically treated region (ion plus first hydration shell) shows shorter O–H bonds and wider H–O–H angles than the bulk solvent.     

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.     

Quantum corrections to classical time-correlation functions: hydrogen bonding and anharmonic floppy modes

Several simple quantum correction factors for classical line shapes, connecting dipole autocorrelation functions to infrared spectra, are compared to exact quantum data in both the frequency and time domain. In addition, the performance of the centroid molecular dynamics approach to line shapes and time-correlation functions is compared to that of these a posteriori correction schemes. The focus is on a tunable model that is able to describe typical hydrogen bonding scenarios covering continuously phenomena from tunneling via low-barrier hydrogen bonds to centered hydrogen bonds with an emphasis on floppy modes and anharmonicities. For these classes of problems, the so-called "harmonic approximation" is found to perform best in most cases, being, however, outperformed by explicit centroid molecular dynamics calculations. In addition, a theoretical analysis of quantum correction factors is carried out within the framework of the fluctuation-dissipation theorem. It can be shown that the harmonic approximation not only restores the detailed balance condition like all other correction factors, but that it is the only one that also satisfies the fluctuation-dissipation theorem. Based on this analysis, it is proposed that quantum corrections of response functions in general should be based on the underlying Kubo-transformed correlation functions.