1549380323-Statistical Mechanics Theory and Molecular Simulation

Problems 313

a. Show that the sequence of distributionsπn(x) satisfies the recursion

πn+1(x) =

∫ 1

x

x

y

πn(y)dy+

∫x

0

πn(y)dy+πn(x)

∫x

0

(

1 −

y

x

)

dy.

b. Show, therefore, thatπn(x) =cxis a fixed point of the recursion, where

cis an arbitrary constant. By normalization,cmust be equal to 2.

c. Now suppose that we start the recursion withπ 0 (x) = 3x^2 and that at

thenth step of the iteration

πn(x) =anx+cnxn+2,

whereanandcnare constants. Show that asn→∞,cngoes asymptot-

ically to 0, leaving a distribution that is purely linear.

∗7.4. In this problem, we will compare some of the Monte Carlo schemesintroduced

in this chapter to thermostatted molecular dynamics for a one-dimensional

harmonic oscillator for which

U(x) =

1

2

mω^2 x^2

in the canonical ensemble. For this problem, you can take the massmand

frequencyωboth equal to 1.

a. Write a Monte Carlo program that uses uniform sampling ofxwithin the

M(RT)^2 algorithm to calculate the canonical ensemble average〈x^4 〉. Try

to optimize the step size ∆ and average acceptance probability to obtain

the lowest possible variance.

b. Write a hybrid Monte Carlo program that uses the velocity Verlet inte-

grator to generate trial moves ofx. Use your program to calculate the

same average〈x^4 〉. Try to optimize the time step and average acceptance

probability to obtain the lowest possible variance.

c. Write a thermostatted molecular dynamics program using the Nos ́e-Hoover

chain equations together with the integrator described by eqns. (4.11.8)-

(4.11.17). Try usingnsy= 7 with the weights in eqn. (4.11.12) andn= 4.

Adjust the time step so that the energy in eqn. (4.10.3) is conserved to

10 −^4 as measured by eqn. (3.14.1).

d. Compare the number of steps of each algorithm needed to converge the

average〈x^4 〉to within the same error as measured by the variance. What

are your conclusions about the efficiency of Monte Carlo versus molecular

dynamics for this problem?

∗7.5. Consider Hamilton’s equations in the form ̇x =η(x). Suppose x = (q,p) is

a two-dimensional phase space vector, and let the equations of motion be

integrated using the following numerical solver:

x(∆t) = x(0) + ∆tη

(

x(0) +

∆t

2

η(x(0))

)

.