mi i

2

24 kB T 2 i

Assuming V can be written as a sum of pair potentials U (r):

1

V= U (rij ) = U (rij ), (3.122)

2

i<j i,j=i

we can evaluate the Laplacian and arrive at

2 2U (rij )

1

fcor = ’ U (rij ) + . (3.123)

2

24 kB T 2 mi rij

i j=i

Next we rewrite the total potential energy on the basis of Feynman“Hibbs

pair potentials:

1

V FH = FH

Vij (rij )

2

i,j=i

2 2U (rij )

1 1 1

cl

=V + + U (rij ) +

2 24 kB T mi mj rij

i,j=i

2 2U (rij )

1

= V cl + U (rij ) + . (3.124)

24 kB T mi rij

i j=i

74 From quantum to classical mechanics: when and how

Expanding exp(’βV FH ) to ¬rst order in the correction term we ¬nd:

⎡ ¤

2 2U (rij ) ¦

1

e’βV = e’βV ⎣1 ’

FH cl

U (rij ) + , (3.125)

24 kB T mi rij

i j=i

which, after integration, gives exactly the fcor of (3.123).

3.5.4 Corrections for high-frequency oscillators

Bond oscillations are often of such a high frequency that ω/kB T > 1 and

order- 2 corrections are not su¬cient to describe the thermodynamics of

the vibrations correctly. A good model for non-classical high-frequency vi-

brations is the harmonic oscillator, which is treated in Chapter 17. Fig-

ure 17.5 on page 478 shows the free energy for the harmonic oscillator for

the classical case, the 2 -corrected classical case and the exact quantum

case. When ω/kB T >≈ 5, the bond vibrational mode is essentially in its

ground state and may be considered as ¬‚exible constraint in simulations.9

For ω/kB T <≈ 2, the 2 -corrected values, which can be obtained by a

proper Feynman“Hibbs potential, are quite accurate. For the di¬cult range

2 < ω/kB T < 5 it is recommended to use the exact quantum correc-

tions. At T = 300 K, this “di¬cult” range corresponds to wave numbers

between 400 and 1000 cm’1 in which many vibrations and librations occur

in molecules. One may also choose not to include the 2 corrections at all

and apply the exact quantum corrections for the full range ω/kB T >≈ 0.5,

i.e., all frequencies above 100 cm’1 .

In a condensed-phase simulation, one does not know all the normal modes

of the system, from which the quantum corrections could be computed. The

best way to proceed is to perform a classical simulation and compute the

power spectrum of the velocities of the particles.10 The power spectrum

will contain prominent peaks at frequencies corresponding to normal modes.

We follow the description by Berens et al. (1983), who applied quantum

corrections to water.

First compute the total mass-weighted velocity-correlation function for an

N -particle ¬‚uid:

3N

C(„ ) = mi vi (t)vi (t + „ ), (3.126)

i=1

9 See the treatment of ¬‚exible constraints in molecular dynamics on page 158.

10 See Section 12.8 for the description of power spectra and their relation to correlation functions.

3.5 Quantum corrections to classical behavior 75

and from this the spectral density of states S(ν):

∞

4

S(ν) = C(„ ) cos 2πν„ d„. (3.127)

kB T 0

Note that the correlation function is the inverse transform of S(ν) (see Sec-

tion 12.8):

∞

C(„ ) = kB T S(ν) cos 2πν„ dν, (3.128)

0

with special case

∞

C(0)

S(ν) dν = = 3N. (3.129)

kB T

0

Now we have the classical density of states, we can compute the quantum

corrections to thermodynamics quantities. The results are (Berens et al.,

1983):

∞

1 ’ e’ξ

’A ln ’ξ/2 ’ ln ξ ,

qu cl

A = kB T dνS(ν) (3.130)

e

0

∞

ξ ξ

’U ’1 ,

qu cl

U = kB T dνS(ν) +ξ (3.131)

2 e ’1

0

∞

ξ

’ ln 1 ’ e’ξ

’S

qu cl

S = kB T dνS(ν) , (3.132)

ξ ’1

e

0

∞

ξ 2 eξ

qu

’ ’1 ,

cl

CV CV = kB T dνS(ν) (3.133)

(1 ’ eξ )2

0

where

hν

ξ= . (3.134)

kB T

One should check if the computed spectral density of states integrate to 3N .

3.5.5 The fermion“boson exchange correction