In the preceding discussion, it has been assumed that the N particles
in the system under consideration could be
treated as classical point particles. In many cases, this treatment
is justifiable, however, there is a large class of systems for which such
an approximation is not valid. In general, systems where hydrogen/proton motion
is important, for example, proton transfer processes, often have significant
nuclear quantum effects. The problem of treating nuclear quantum effects
in a system at finite temperature requires the solution of a quantum statistical
mechanical problem. One approach that has been applied with considerable
success is based on the Feynman path integral formalism of
statistical mechanics [94, 95].

Consider the quantum canonical partition
function for a single particle in one spatial dimension. The partition
function is given by the trace:

where the trace is carried out in the coordinate basis.
Assuming H=T+U, where T is the kinetic energy operator, and U is the
potential, the Trotter theorem, Eq. (6.4), allows
to be expressed as
in the limit . The Trotter theorem expression for
is then substituted into Eq. (8.1), and
an identity operator in the form of
is inserted in between each factor of ,
yielding

Then, using the fact that

one obtains the final expression for Q as a function of P

where is an effective potential given by

with .

Equation (8.4) is in the form of a configurational partition function
for a P-particle system in one dimension subject to a potential .
The configurational partition function can also be expressed in
a quasi phase space form by recognizing that the prefactor can be written as a product
of P uncoupled Gaussian integrals:

where

In Eq. (8.6), the constant is an overall
constant that ensures equality of
Eqs. (8.6) and (8.4).
In addition, the mass m', being a fictitious mass, is arbitrary, a fact that can be
exploited in devising an MD scheme for Eq. (8.6) as will be shown below.
As was pointed out by Chandler and Wolynes [96],
Eqs. (8.6) and (8.7) together show that,
for finite P, the path integral of a single quantum particle is
isomorphic to a classical system of P particles subject with a Hamiltonian
given by Eq. (8.7).
Inspection of Eq. (8.5) shows that the P particles form
a closed polymer chain with nearest neighbor
harmonic coupling and are subject to a potential U.
The classical isomorphism allows molecular dynamics to be used to
simulate a finite-temperature quantum system. The extension of the path integral
scheme to N particles in three dimensions is straightforward if it is assumed
that the particles obey Boltzmann statistics, i.e., all spin statistics
are neglected. In this case, the partition function is

where the classical Hamiltonian is given by

In principle, the equations of motion resulting from Eq. (8.9)
could be implemented as a MD procedure, from which quantum equilibrium properties
of a system could be computed [97].
A number of well known difficulties arise in a straightforward
implementation of MD to the path integral. Primarily, since ,
the force constant of the harmonic coupling increases as P increases, giving
rise to a stiff harmonic interaction and a time scale separation. As was shown
by Hall and Berne [98], this time scale separation gives rise to
non-ergodic trajectories that do not sample the available canonical phase space.
A solution to this problem was first presented in Ref. [30]. There,
it was shown that several elements are needed to devise an efficient MD scheme
for path integrals. First, a change of variables that diagonalizes
the harmonic coupling is introduced. This has the effect of isolating the various
time scales present in the Hamiltonian of Eq. (8.9). The
change of variables is linear, having the general form

where the matrix is a constant matrix of unit determinant.
Two different choices of the matrix , discussed in Ref. [99]
and [100], lead to the staging and normal mode
transformations. The transformed coordinates are known
as staging or normal mode variables. If the change of variables is made in
Eq. (8.8), then the corresponding classical Hamiltonian takes the form:

where the s-dependent masses result from the variable transformation.
For a staging transformation, the
masses are , for ,
while for the normal mode transformation, the
masses are proportional to the normal mode eigenvalues.
Thus, it is clear that the fictitious masses should be chosen according to
and . In this way, all modes
will move on the same time scale, leading to maximally efficient exploration
of the configuration space.

In addition to variable transformations, it is necessary to ensure that
a canonical phase space is generated. This can be achieved via one of
the non-Hamiltonian MD schemes for generating the NVT ensemble.
It has been found that maximum efficiency is obtained if each Cartesian
direction of each mode variable is coupled to its own thermostat, as
was clearly demonstrated in Ref. [99], and multiple time scale
integration techniques are employed [30].

It is worth mentioning that the path integral MD scheme outlined
here has been combined with ab initio MD to yield an ab initio
path integral Car-Parrinello method [101, 99].
This allows quantum effects on chemical
processes to be studied. More recently, the ab initio path integral
scheme has been extended to incorporate approximate quantum dynamical
properties [102]
via the so called centroid dynamics method [103, 104]. Finally, the
path integral MD scheme has been modified to allow path integral simulations under
conditions of constant temperature and pressure to be carried out [100].