Exercises 257
(a) Using the leap-frog/Verlet algorithm, show that
pi(t+h/ 2 )=pi(t)+hFi/2.
The refreshed momentapi(t)are drawn from a Maxwell–Boltzmann distribution,
and the momenta at timet+h/2, which are needed in the Verlet/leap-frog
algorithm are then calculated using this last formula.
(b) Implement the Andersen update algorithm for argon and compare the results with
the microcanonical program.
(c) Now suppose that the momenta are refreshed ateverystep. Show that in that case
we have
ri(t+h)=ri(t)+h^2 Fi/ 2 +hζi(t),
whereζi(t)is theith random momentum component generated according to the
Maxwell–Boltzmann distribution. This is a kind of Langevin equation. Discuss
the difference with the Langevin equation described in Section 8.8.
8.9 [C] In this problem we consider the implementation of the Nosé–Hoover thermostat
in the microcanonical MD simulation for Lennard–Jones argon described in
Section 8.3. The extension is straightforward – the equations are given in
Section 8.5.1. You can verify now that the behaviour of the Nosé–Hoover thermostat
is often nonergodic. ForT=1.5 andρ=0.8 the behaviour is as it should be for
coupling constantQ=1. You can check that the standard deviation in the
temperature is in accordance withEq. (8.81). For lower temperatures, likeT=0.85,
ρ=1.067, the temperature exhibits large oscillations. The period of these
oscillations depends onQ[26].
8.10 (a) Verify that when we takeg= 3 Ninstead ofg= 3 N+1 in the derivation of the
Nosé–Hoover thermostat, the probability density for configurations(P,R)turns
out to be:
ρ(P,R)=
1
3 N(
2 πQ
kBT) 1 / 2
exp[
−H 0 (P,R)( 3 N+ 1 )
3 NkBT]
.(b) For this choice, verify that quantities sampled in a simulation yield averages as
given inEq. (8.94).
8.11 [C] In this problem, a code for evaluating the potential felt by the particles in a
two-dimensional Coulomb (or gravitational) system is developed, using the
tree-code method ofSection 8.7.2.
Although experienced programmers would be tempted to start building tree
structures using pointers and recursive programming for this problem, it can be dealt
with using more pedestrian methods. The point is that the squares can be coded by
two integersNX, NYwhich are considered as bit-strings. The first of these contains
information about thex-coordinate of the square and the second about the
y-coordinate. They are ordered linearly: the leftmost column of squares hasNX=0,
the rightmost columnNX= 2 n−1 etc., and a similar coding is adopted for the
rows. If squares are neighbours, their respectiveNXandNY-codes should differ at
most by 1 (and they should not be equal). The codes of the parents can be found