5Phase transitions
II Statistical Physics
5.3 The Ising model
The Ising model was invented as a model of a ferromagnet. Consider a
d
-
dimensional lattice with
N
sites. On each of these sites, we have a degree of
freedom called the “spin”
s
i
=
(
+1 spin up
−1 spin down
.
The Hamiltonian is given by
H = −J
X
hiji
s
i
s
j
− B
X
i
s
i
for some
J, B
. We should think of
B
as the magnetic field, and for simplicity,
we just assume that the magnetic moment is 1. The first term describes interac-
tions between the spins themselves, where
P
hiji
denotes the sum over nearest
neighbours.
For example in d = 1, a lattice looks like
and the nearest neighbours are the literal neighbours.
In d = 2, we can have many different lattices. We could have a square grid
The nearest neighbours of the central red dot are the orange ones.
We let
q
be the number of nearest neighbours. So for example, when
d
= 1,
we have q = 2; for a d = 2 square lattice, then q = 4.
What is
J
? It is clearly a measure of the strength of the interaction. But the
sign matters. The case
J >
0 corresponds to the case where neighbouring spins
prefer to align. This is called a ferromagnet; and the
J <
0 case corresponds to
the case where they prefer to anti-align, and this is called an anti-ferromagnet.
It doesn’t really matter a lot, but for the sake of definiteness, we assume
J > 0.
Now if turn on the
B
field, then on the basis of energy, it will make the spins
to try to align with the
B
field, and the interactions will further try to make it
align.
But energy is not the only thing that matters. There is also entropy. If all
spins align with the field, then entropy will be zero, which is not very good. So
there is a competition between energy, which wants to align, and entropy, which
prefers to anti-align, and the end result depends on how this competition turns
out.
To understand this properly, we use the canonical ensemble
Z =
X
{s
i
}
e
−βE[{s
i
}]
.
We define the average spin, or the magnetization, to be
m =
1
N
X
i
hs
i
i =
1
Nβ
∂ log Z
∂B
T
. (†)
It turns out there is a different possible interpretation of the Ising model. It can
also be viewed as a description of a gas! Consider “hard core” particles living on
a lattice. Then there is either 0 or 1 particle living on each site
i
. Let’s call it
n
i
.
We suppose there is an attractive force between the particles, and suppose the
kinetic energy is not important. The model is then defined by the Hamiltonian
H = −4J
X
hiji
n
i
n
j
.
This is a rather crude model for studying gas on a lattice. To understand this
system of gas, we use the grand canonical ensemble
Z
gas
=
X
{n
i
}
e
−β(E[{n
i
}]−µ
P
i
n
i
)
.
We look at the exponent
E[{n
i
}] − µ
X
i
n
i
= −4J
X
hiji
n
i
n
j
− µ
X
i
n
i
.
We can compare this with the Ising model, and they look very similar. Unsur-
prisingly, this is equal to E
Ising
[{s
i
}] if we set
n
i
=
s
i
+ 1
2
,
and let
B =
µ
2
+ qJ.
Then we have
Z
gas
= Z
Ising
.
So if we understand one, then we understand the other.
We now concentrate on
Z
Ising
. Can we actually compute it? It turns out
it is possible to do it for
d
= 1, which we will do some time later. In the case
d
= 2 with
B
= 0, we can still manage, but in general, we can’t. So we want to
develop some approximate method to study them.
We’ll do so using mean field theory. In general, we would expect
hs
i
i
=
m
.
So, being a good physicist, we expand around this equilibrium, and write
s
i
s
j
= [(s
i
− m) + m][(s
j
− m) + m)]
= (s
i
− m)(s
j
− m) + m(s
j
− m) + m(s
i
− m) + m
2
.
We will assume that
X
hiji
(s
i
− m)(s
j
− m)
is negligible compared to other terms in
H
. This is a subtle assumption. It is
not true that
h
(
s
i
−m
)
2
i
is negligible. In fact, since
s
i
= 1 always, we must have
hs
2
i
i = 1. Thus,
h(s
i
− m)
2
i = hs
2
i
i − 2mhs
i
i + m
2
= 1 − m
2
,
and for small magnetization, this is actually pretty huge. However,
hs
i
s
j
i
can
behave rather differently than
hs
i
i
2
. If we make suitable assumptions about the
fluctuations, this can potentially be small.
Assuming we can indeed neglect that term, we have
H ≈ −J
X
hiji
(m(s
i
+ s
j
) − m
2
) − B
X
i
s
i
.
Now the system decouples! We can simply write this as
H =
1
2
JNqm
2
− (Jqm + B)
X
s
i
.
Now this is just the 2-state model we saw at the very beginning of the course,
with an effective magnetic field
B
eff
= B + Jqm.
Thus, we have “averaged out” the interactions, and turned it into a shifted
magnetic field.
Using this, we get
Z ≈ e
−
1
2
βJNqm
2
+B
eff
P
i
s
i
= e
−
1
2
βJNqm
2
e
−βB
eff
+ e
βB
eff
N
= e
−
1
2
βJNqm
2
2
N
cosh
N
(βB + βJqm)
Note that the partition function contains a variable
m
, which we do not know
about. But we can determine m(T, B) using (†), and this gives
m = tanh(βB + βJqm). (∗∗)
It is not easy to solve this analytically, but we can graph this and understand it
qualitatively:
We first consider the case B = 0. Note that for small x, we have
tanh x ≈ x −
1
3
x
3
.
We can then plot both sides of the equation as a function of
m
. If
βJq <
1, then
we see that the only solution is m = 0:
m
On the other hand, if βJq > 1, then we have three possible states:
m
m
0
(T )
Something interesting happens when βJq = 1. We define T
c
by
kT
c
= Jq.
Then we have
βJq =
T
C
T
.
We can now rephrase our previous observation as follows — if
T > T
C
, then
there is only one solution to the system, namely
m
= 0. We interpret this as
saying that at high temperature, thermal fluctuations make them not align. In
this case, entropy wins.
If
T < T
c
, then there are 3 solutions. One is
m
= 0, and then there are
non-trivial solutions
±m
0
(
T
). Let’s look at the
m
= 0 solution first. This we
should think of as the unstable solution we had for the liquid-gas equation.
Indeed, taking the derivative of (∗∗) with respect to B, we find
∂m
∂B
T
B=0
=
β
1 − βJq
< 0.
So for the same reason as before, this is unstable. So only
±m
0
(
T
) is physical.
In this case the spins align. We have
m
0
→ ±
1 as
T →
0, and the two signs
correspond to pointing up and down.
The interesting part is, of course, what happens near the critical temperature.
When we are close to the critical temperature, we can Taylor expand the equation
(∗∗). Then we have
m
0
≈ βJqm
0
−
1
3
(βJqm
0
)
3
.
We can rearrange this, and we obtain
m
0
∼ (T
c
− T )
1/2
.
We can plot the solutions as
T
m
−1
1
T
c
There is a phase transition at
T
=
T
C
for
B
= 0. For
T > T
C
, the stable states
have
m
= 0, while for
T < T
C
the stable states have non-zero
m
. Thus, the
magnetization disappears when we are above the critical temperature T
C
.
We now want to understand the order of the phase transition. For general
B, we can write the free energy as
F (T, B) = −kT log Z =
1
2
JNqm
2
− NkT log (2 cosh(βB + βJqm)) ,
where m = m(T, B), which we found by solving (∗∗).
Restricting back to
B
= 0, for
T ≈ T
C
, we have a small
m
. Then we can
expand the expression in powers of m, and get
F (T, 0) ≈
1
2
NkT
C
1 −
T
C
T
m
2
− NkT log 2.
The second term does not depend on
m
, and behaves in a non-singular way. So
we are mostly interested in the behaviour of the other part. We have
F + NkT log 2 ∼
(
−(T
C
− T )
2
T < T
C
0 T > T
C
If we look at the first derivative, then this is continuous at
T
=
T
C
! However,
the second derivative is discontinuous. Hence we have a second-order phase
transition at B = 0, T = T
C
.
We’ve understood what happens when
B
= 0. Now let’s turn on
B
, and put
B 6= 0. Then this merely shifts the tanh curve horizontally:
m
m(T, B)
We see that we always have a unique solution of
m
of the same sign as
B
. We
call this
m
(
T, B
). However, when we are at sufficiently low temperature, and
B
is not too big, then we also have other solutions
U
and
M
, which are of opposite
sign to B:
m
M
U
m(T, B)
As before,
U
is unstable. But
M
is something new. It is natural to ask ourselves
which of the two states have lower free energy. It turns out
M
has higher
F
, and
so it is metastable. This state can indeed be constructed by doing very delicate
experiments.
This time, at least qualitatively from the pictures, we see that
m
(
T, B
)
depends smoothly on T and B. So there is no phase transition.
We can also see what happens at high temperature. If
βJq
1 and
β|B|
1,
then we can expand (∗∗) and get
m ≈ βB =
B
kT
.
This is Curie’s law . But now let’s do something different. We fix
T
, and vary
B
.
If T > T
c
, then m(T, B) depends smoothly on B, and in particular
lim
B→0
m(T, B) = 0.
However, if we start lower than the critical temperature, and start at
B >
0,
then our magnetization is positive. As we decrease
B
, we have
B →
0
+
, and
so
m → m
0
(
T
). However, if we start with negative
B
, we have a negative
magnetization. We have
lim
B→0
+
m(T, B) = m
0
(T ) 6= −m
0
(T ) = lim
B→0
−
−m
0
(T ).
So we indeed have a phase transition at B = 0 for T < T
C
.
We now want to determine the order of this transition. We note that we have
m
1
Nβ
∂ log Z
∂B
T
= −
1
N
∂f
∂B
T
.
Since
m
is discontinuous, we know this first derivative of the free energy is
discontinuous. So this is a first-order phase transition.
We can plot our rather unexciting phase diagram:
B
T
We can compare this with the liquid-gas phase diagram we had:
T
p
critical point
liquid
gas
It looks roughly the same.
Let’s now compute some critical components. We saw that for
B
= 0 and
T → T
−
C
, we had
m
0
∼ (T
C
− T )
β
,
with β =
1
2
.
We can also set T = T
C
, and let B → 0. Then our equation (∗∗) gives
m = tanh
B
T q
+ m
.
We can invert this equation by taking tanh
−1
, and get
B
Jq
+ m = tanh
−1
m = m +
1
3
m
3
+ ··· .
So we find
m ≈
3B
Jq
1/3
.
So we find B ∼ m
δ
, where δ = 3.
Finally, we can consider the susceptibility
χ = N
∂m
∂B
T
.
We set B = 0, and let T → T
+
C
. We then find that
χ =
Nβ
1 − βJq
∼ (T − T
C
)
−γ
,
where γ = 1.
These were exactly the (incorrect) critical exponents we computed for the
liquid-gas phase transition using the van der Waals equation!
In the case of van der Waals, we saw that our result was wrong because we
made the approximation of ignoring fluctuations. Here we made the mean field
approximation. Is this a problem? There are two questions we should ask —
is the qualitative picture of the phase transition correct? And are the critical
exponents correct?
The answer depends on the dimension.
–
If
d
= 1, then this is completely wrong! We can solve the Ising model
exactly, and the exact solution has no phase transition.
–
If
d ≥
2, then the phase diagram is qualitatively correct. Moreover, the
critical exponents correct for
d ≥
4. To learn more, one should take III
Statistical Field Theory.
–
In
d
= 2, we can indeed solve the Ising model exactly, and the exact
solution gives
β =
1
8
, γ =
7
4
, δ = 15.
So the mean field approximation is pretty useless at predicting the critical
exponents.
–
In
d
= 3, then there is no exact solution known, but there has been a lot
of progress on this recently, within the last 2 or 3 years! This led to a
very accurate calculation of the critical exponents, using a combination of
analytic and numerical methods, which gives
β = 0.326, γ = 1.237, δ = 4.790.
These are actually known to higher accuracy.
These are exactly the same as the measured critical exponents for the
liquid-gas transition!
So this funny model of spins on a lattice is giving us exactly the same number as
real-world liquid-gas! This is evidence for universality. Near the critical point,
the system is losing knowledge of what the individual microscopic description
of the system is, and all systems exhibit the same kind of physics. We say the
critical points of the three-dimensional Ising model and liquid-gas system belong
to the same universality class. This is something described by a conformal field
theory.
One-dimensional Ising model
Let’s now solve the one-dimensional Ising model. Here we just have
H = −J
N
X
i=1
s
i
s
i+1
−
B
2
n
X
i=1
(s
i
+ s
i+1
).
To make our study easier, we impose the periodic boundary condition
s
N+1
≡ s
1
.
We then have
Z =
1
X
s
1
=−1
1
X
s
2
=−1
···
1
X
s
N
=−1
N
Y
i=1
exp
BJs
i
s
i+1
+
βB
2
(s
i
+ s
i+1
)
.
We define the symmetric 2 × 2 matrix T by
T
st
= exp
BJst +
βB
2
(s + t)
,
where s, t = ±1. We can then rewrite this Z using matrix multiplication —
Z =
X
s
1
···
X
s
N
T
s
1
s
2
T
s
2
s
3
···T
s
N
s
1
= Tr(T
N
).
The trace is defined to be the sum of eigenvalues, and if we know the eigenvalues
of T , we know that of T
N
. We can do this directly. We have
λ
±
= e
βJ
cosh βB ±
q
e
2βJ
cosh
2
βB − 2 sinh 2βJ.
This is not terribly attractive, but we can indeed write it down. As expected,
these eigenvalues are real, and we picked them so that
λ
+
> λ
−
. Then we have
Z = λ
N
+
+ λ
N
−
.
We have thus solved the 1
d
Ising model completely. We can now see if there is a
phase transition. We have
m =
1
Nβ
∂ log Z
∂B
T
=
1
βZ
λ
N−1
+
∂λ
+
∂B
T
+ λ
N−1
−
∂λ
−
∂B
T
.
To evaluate this derivative, we see that the
B
-dependence lives in the
cosh
terms.
So we have
∂λ
±
∂B
T
∝ sinh βB.
But this is all we need, so we only need to evaluate it at
B
= 0, and
sinh
vanishes.
So we know that when
B
= 0, we have
m
= 0 for all
T
! So there is no phase
transition at B = 0.
More generally, we can look at the free energy
F = −
1
β
log Z
= −
1
β
log
λ
N
+
1 +
λ
−
λ
+
N
!!
= −
N
β
log λ
+
−
1
β
log
1 +
λ
−
λ
+
N
!
.
Recall that we said phase transitions are possible only in the thermodynamic
limit. So we want to take N → ∞, and see what we get. In this case, we have
λ
−
λ
+
N
→ 0.
So we have
F
N
→ −
1
β
λ
+
.
We have
λ
+
>
0, and this depends smoothly on
T
and
B
(as
√
·
never vanishes).
So there is no phase transition.
In fact, this holds for any 1d system without long-range interactions (Peierls).