Dynamical properties from quantum correlation functions
Real time quantum dynamics simulation of large molcular systems present a challenge to theoretical physics and chemistry. Indeed, computer simulation of time-dependent quantities for quantum systems is hindered by the exponential scaling of exact methods with the number of degrees of freedom. However, access to these quantities is crucial in many interesting problems and in particular to compute, via linear response theory and time correlation functions, the response of the system to external perturbations. The development of approximate quantum dynamical schemes is therefore an active area of research. One of the main directions currently pursued is to adapt efficient techniques for classical simulations, based on molecular dynamics (MD) and/or Monte Carlo, to the quantum case. A reasonable compromise among accuracy and computational cost is realised by the so-called quasi-classical methods for computing time correlation functions. In these methods, sampling of the thermal equilibrium density for the system is combined with generalised trajectories that provide an approximation for quantum dynamics usually valid for short times.
In this section, we will present different quasi-classical methods, containing different approximations. First, we will introduce the different forms of quantum correlation functions which is possible to compute. These expressions are equivalent in the sense that they can be related to one another analytically and therefore have the same physical content. It is however useful to discuss them all since, depending of the adopted approximation, some may be more advantageous from a computational point of view. Then, we will introduce two of the most popular quasi-classical methods, which, although not formally justified, have the advantage of converging with a small number of trajectories: centroid and ring polymer molecular dynamics (CMD and RPMD, respectively) [2, 29, 30, 31, 32, 33]. Finally, we will describe three methods, which are part of the category of the linearized quasi-classical methods. The Feynman-Kleinert linearized path integral (FK-LPI) method introduced by Poulsen and Rossky [1, 34, 35, 36, 37, 38], the linearized semi classical initial value representation (LSC-IVR) with the local Gaussian approximation (LGA) describe by Liu and Miller [4, 13, 19, 39, 40, 41] and finally our method, the phase integration method (PIM) introduced by Bonella and Ciccotti [7, 8, 9, 12, 17, 18].
Time correlation functions play a crucial role in relating macroscopic observations accessible in experiments to the microscopic dynamics of physical systems. Although the statstical mechanics required to establish this link, in particular linear response theory, is valid for classical and quantum systems, our ability to perform calculations differs dramatically in the two situations as explained before. Furthermore, for quantum systems, different definitions of correlation functions are possible due to the fact in quantum statistical we deal with operators which do not commute.
Empirical quasi-classical methods: RPMD and CMD
Now we know what we have to calculate but as explained before we cannot do it fully quantum mechanically. As a consequence we have to perform approximations to calculate these quantum correlation functions. We will present two of the most popular path integral methods to compute quantum time correlation function: RPMD and CMD. As we are going to see, these methods are not formally justified but are very efficient numerically. Then, we will describe the linearization approximation and different methods, including PIM, which are based on this approximation.
The Linearized Path Integral (LPI) representation of quantum correlation functions
I will present here the linearized path integral approximation on the standard quantum correlation function but the idea is stricly the same for the Schofield or the Kubo correlation functions. We first recall the equation 2.64, which corresponds to the definition of the standard quantum correlation function: CAB(t; ) = 1 Z Tr h e bH b AeitbH =~ b BeitbH =~ i.
The Linearized Semiclassical Initial Value Representation (LSC-IVR) method with a Local Gaussian Approximation (LGA)
The LSC-IVR [21, 57, 22] approximation is stricly equivalent to the LPI approximation detailed before. The originality of the method develop by Liu and Miller is to combine this approximation with the Local Gaussian Approximation (LGA) [13, 20] for sampling the Wigner density. As a consequence, I will insist here on the Local Gaussian Approximation and why it is a good method to compute Wigner densities.
Calculation of the Wigner function for operator b B in equation 2.127 is usually straightforward; in fact, b B is often a function only of coordinates or only of momenta, in which case its Wigner functions is simply the classical function itself. However, as shown in the previous section, calculating the Wigner transform of e bH b A is far from being trivial. Nonetheless, this is necessary in order to obtain the distribution of initial conditions of coordinates x0 and momenta p0 for the real time trajectories. Usually, local harmonic approximations of the potential (such as the Local Harmonic Approximation (LHA) of Shi and Geva , the FK-LPI of Poulsen and Rossky , …) are adapted and give good results to many quantum simulations. These methods, however, fail when the imaginary local frequency and temperature are such that ~j!j and so are not reliable for problems dominated by potential barriers at low temperature.
Edgeworth expansion for the Wigner density
The manipulations presented above give, by construction, a positive definite Wigner density. While this is usually not considered a serious problem for most realistic applications (and indeed, essentially all existing approximation schemes have the same limitation), it is interesting to explore alternative representations that can reproduce also negative features of the Wigner density if they exist. To that end, let us consider again the integral over r at the left hand side of equation 3.17. It is the characteristic function of the conditional probability c(rjr). We have previously introduced a cumulant expansion of this quantity, which, for all practical purposes, will be truncated in the following at second order thus approximating this integral by a Gaussian in p, with an r dependent variance (covariance matrix in the multidimensional case). By taking this Gaussian as a reference it is then possible to trivially rewrite the integral as: Z dr e i ~ prc(rjr) = e2 2 p2 ~2.
Kubo correlation functions for Infrared spectroscopy
We will derive in this section the PIM version of the Kubo correlation functions. The first case will be for operators linear in positions (it is possible to do it for more general position operators but we will not use it). Indeed, this particular case will show us how the computational cost associated to the Kubo expression is reduced compared to the symmetrised case on the most simple test case: the position autocorrelation function of a 1D harmonic oscillator.
Linear operators in the position, as we can see on the equations 4.16, are relevant for IR spectra and are appropriate for the isolated molecules considered here. Application of the methods to bulk systems in periodic boundary conditions is, however, simpler starting from the dipole-derivation expression in equation 4.17. I will then discuss how we can calculate Kubo correlation functions with operators linear in momentum.
Table of contents :
2 State of the art
2.1 Formalism of path integrals
2.1.1 Path integrals and the canonical density matrix
2.1.2 Polymer isomorphism and PIMD
2.1.3 Time-independent equilibrium properties in the canonical ensemble from path integral formulation
2.2 Dynamical properties from quantum correlation functions
2.3 Empirical quasi-classical methods: RPMD and CMD
2.3.1 Outline of RPMD
2.3.2 Outline of CMD
2.4 Linearized methods
2.4.1 The Linearized Path Integral (LPI) representation of quantum correlation functions
2.4.2 The FK-LPI method
2.4.3 The Linearized Semiclassical Initial Value Representation (LSC-IVR) method with a Local Gaussian Approximation (LGA)
3 Wigner densities with the Phase Integration Method
3.1 Wigner densities via the phase integration method (PIM)
3.1.1 The PIM expression for the thermal Wigner density
3.1.2 Edgeworth expansion for the Wigner density
3.2 Structure of the PIM algorithm
3.2.1 Structure of the algorithm
3.2.2 Definition and evaluation of the numerical estimators
3.3.1 Harmonic oscillator
3.3.2 Morse potential
3.3.3 Proton transfer model
4 Application of PIM to the Infrared spectroscopy
4.1 Symmetrised correlation function for Infrared spectroscopy
4.2 Kubo correlation functions for Infrared spectroscopy
4.2.1 Operators linear in positions
4.2.2 Operators linear in momentum
4.2.3 Application to the Infrared spectroscopy
5 Application of PIM to the calculation of rate constants
5.1 General expression for chemical rate constant
5.2 PIM Flux-side correlation function
5.3 Free energy calculations
5.4.1 Free energy profiles
5.4.2 Rate constants
A Modification of the algorithm for the multidimensional case
B Monte Carlo sampling of the polymer chains
B.1 Staging variables
C Infrared spectroscopy
C.1 Absorption coefficient from the Fermi Golden rule
D Alternative demonstration for the Kubo momentum autocorrelation function
D.1 Notations of Kubo
D.2 Momentum autocorrelation function
D.3 Equipartition of the energy
E Eckart transformation