488 The Feynman path integral
〈α〉f=
CP
2 QP
∫
dx 1 ···dxP
∑P
k=1
xk
∂α
∂xk
e−βα(x^1 ,...,xP)e−βγ(x^1 ,...,xP)
=−
CP
2 βQP
∫
dx 1 ···dxP
∑P
k=1
(
xk
∂
∂xk
e−βα(x^1 ,...,xP)
)
e−βγ(x^1 ,...,xP). (12.6.38)
whereCP= (mP/ 2 πβ ̄h^2 )^1 /^2. Integrating eqn. (12.6.38) by parts yields
〈α〉f=
CP
2 βQP
∫
dx 1 ···dxPe−βα(x^1 ,...,xP)
∑P
k=1
∂
∂xk
[
xke−βγ(x^1 ,...,xP)
]
=
CP
2 βQP
∫
dx 1 ···dxP
[
P−
β
P
∑P
k=1
xk
∂U
∂xk
]
e−βγ(x^1 ,...,xP)e−βα(x^1 ,...,xP)
=
P
2 β
−
〈
1
2 P
∑P
k=1
xk
∂U
∂xk
〉
f
, (12.6.39)
from which eqn. (12.6.33) follows.
Using the virial estimator, we now present a comparison (Fig. 12.13)of path-
integral molecular dynamics with no variable transformations (top row), path-integral
molecular dynamics with a staging transformation (middle row), and staging path-
integral Monte Carlo (bottom row). The system is a one-dimensional harmonic os-
cillator withβ ̄hω= 15.8,mω/ ̄h= 0.03, andP = 400. With these parameters, the
thermodynamic energy is dominated by the ground-state value ̄hω/2. The figure shows
the instantaneous value of the virial estimator (left column), the cumulative average
of the virial estimator (middle column), and the error bar in the valueof the esti-
mator. The error bar is calculated by grouping individual samplings from molecular
dynamics or Monte Carlo into blocks of sizen, computing the average over each block,
and then computing the error bar from these block averages with respect to the global
average (Cao and Berne, 1989). The purpose of this type of “block averaging” is to
remove unwanted correlations between successive samplings. As the right column in-
dicates, the error bar starts off small and then reaches a plateauwhen the blocks are
large enough that correlations are no longer present. This exampledemonstrates that
without variable transformations, path-integral molecular dynamics performs rather
poorly, but when a staging transformation is employed, molecular dynamics and stag-
ing Monte Carlo are equally efficient, as evidenced by the fact that they converge to
the same error bar after the same number of steps.
As written, the virial estimator in eqn. (12.6.34) is only valid for boundsystems
because the termxk(∂U/∂xk) is not translationally invariant. This problem can be
circumvented by defining the virial part of the estimator in terms ofthe path centroid.
The generalization of eqn. (12.6.33) forNparticles inddimensions is
dNP
2 β
−
〈N
∑
i=1
∑P
k=1
1
2
miωP^2
(
r
(k)
i −r
(k+1)
i
) 2 〉
f
=
〈
1
2 P
∑N
i=1
∑P
k=1
r
(k)
i ·
∂U
∂r(ik)
〉
f