Computational Physics

(Rick Simeone) #1
10.3 Importance sampling through Markov chains 309

the trial configuration is accepted or rejected as usual, with an acceptance rate
exp[−βE(X→X′)].
The displacement volume remains to be specified. Any cube should lead to
the correct behaviour as the requirements of the Markov chain (connectedness
and aperiodicity) and the detailed balance condition are fulfilled, but there will
obviously be a particular cube size yielding optimum efficiency, where efficiency
can be defined as the number of statistically independent configurations that we can
generate in a fixed amount of computer time. Obviously, if the cube is chosen very
small, particles are allowed to move over small distances only and it will take many
MC steps to arrive at a statistically independent configuration. On the other hand,
if the cube is chosen too large, the particles will on average make large moves.
This changes the configuration of the system to such an extent that the energy will
increase strongly in the large majority of cases. The probability that the trial move
is accepted will then be vanishingly small and we will spend most of our time
generating configurations which are subsequently rejected, so the configuration
changes very slowly in that case as well.
A widely adopted rule of thumb is to choose the cube such that the acceptance
rate of trial states is on average somewhere between 0.4 and 0.6. Although for hard
sphere systems it should perhaps be lower (around 0.1) [11], this rate is generally
relied upon as being reasonably efficient [2].
When programming the Metropolis method for the case of argon, much of the
MD code (see Section 8.3) can be copied. The particles are again released from face-
centred cubic (fcc) lattice positions. Every trial displacement should be performed
respecting the periodic boundary conditions. Calling half the linear size of the cube
dmax, the trial displacement (without correcting for PBC) in thex-direction is given
in terms of a random number 0≤r<1as


xnew−xold=dmax( 2 r− 1 ) (10.28)

and similarly foryandz. A neighbour list may be kept as in the MD case, but the
list must be constructed such that we have at our disposal all the neighbours (at
a distance smaller than the neighbour list cut-offrcut−off) for every particle. This
means that we need twice as much storage as in the MD method, where only the
neighbours with an index higher than the particle under consideration are stored in
the list. As in the Ising case, the total potential energy can be updated after every
acceptance of a trial configuration. The same can be done with the virial. Physical
properties to be calculated are the same as in MD simulations: the pressure, the
configurational energy, and the pair correlation function. Any quantity dependent
on the positions only can be determined – remember the momenta are not considered
in the MC method.

Free download pdf