The Monte-Carlo method
As we know, it is impossible to calculate statistical averages directly.
One could generate spins configurations at random and approximate the real
thermal averages by Monte-Carlo averages. The next problem here is, most
randomly chosen states will make a negligible contribution to the sum.
Random walks are used to take into account only important spins configurations
by the introduction of a fictitious dynamics (so-called Markov chain of
configurations). This wandering generates states, which are the most probable
from energy-entropy point of view (importance sampling).
The Metropolis algorithm
In the Metropolis algorithm we try to turn over a single spin direction with
where E1, E2 are energies of the old and new
configurations (see details in the Gould and Tobochnik book).
Statistical averages may then be computed as simple arithmetic means.
The Metropolis algorithm can by applied easy e.g. to the XY model simulation.
if E1 < E2
W12 = 1
if E1 > E2
There are many possible choices for W12.
To get equilibrium final spins distribution it is enough that
Wij satisfy the detailed balance conditions
W12 exp(-E1/T) =
W21 exp(-E2/T) or
W12 / W21 =
This Java applet is based on the Gould and Tobochnik book.
A successive loop over all sites (referred as one Metropolis sweep or iteration)
is used. Then "momentary" (in Markov's time) Et, Mt
values are calculated
<M>M = 1/T
∑ t=1,T Mt ,
<E>M = 1/T
∑ t=1,T Et .
Susceptibility χ and specific heat C
of the system can be found by the thermal fluctuations of magnetization
M and energy E
C = d<E>/dT = n/T2
(<E2> - <E>2) ,
χ = limH→ 0 dM/dT
= n/T (<M2> - <M>2) .
Controls Click mouse into lattice to get a new spins configuration.
Drag by mouse "thermometer" to the right (the black bar corresponds to
Tc) or press "Enter" to set a new T value from
the text field.
High non-equilibrium spin flip dynamics
You see below two features in the Metropolis algorithm iterations.
First of all, as since transition probability W1->2 =
exp[(E1-E2)/T] -> 1 for high
temperature T >> (E1 - E2) ~ J , therefore
spins oscillate in big solid clusters. Second, you see growth of dendrite
structures in the left-bottom part of the pictures,
because all spins are flipped by lines from this corner
(to avoid this feature one could choose flipped spins at random).
The thermostat algorithm
In the thermostat algorithm we get a single spin (all the rest are
fixed) into contact with big thermostat at temperature T. Then we
get for transition probabilities (the Glauber formula)
P+ = exp(-E+ /T) /
[exp(-E+ /T) + exp(-E- /T)] =
exp(-(E+-E- )/T) /
[exp(-(E+-E- )/T) + 1] ,
P- = exp(-E- /T) /
[exp(-E+ /T) + exp(-E- /T)] = 1 - P+ .
The Glauber transition probabilities satisfy the detailed balance conditions
too. For high temperature P+ , P- → 1/2
and it leads to fast and smooth "cluster melting" (see the Ising applet with
the thermostat algorithm below).
But for low temperature (e.g. T = 1) and random initial
configurations you can see formation of metastable clusters. Then they are
decreasing slowly and one single cluster will appear. In the Metropolis method
"freezing" is faster, but sometimes we get two clusters with a big metastable
Due to ambiguity of choice of transition probabilities this dynamics
seems me a bit artificial, but it influences on spins relaxation greatly.
G.T.Barkema, M.E.J.Newman New Monte Carlo algorithms for classical
spin systems arXiv:cond-mat/9703179
See also Kenji Harada's
applet demonstrating Metropolis, Swendsen and Wang, Wolff algorithms
for Ising model.
Previous: 2D Ising model
Next: Relaxation time
updated 27 Dec 2001