Ising Model

Definition

The Ising model is a lattice of binary spins that interact only with their immediate neighbours. Each site i carries a variable that can take just two values:

s_i ∈ { +1, −1 }   (spin up or spin down)

The total energy of a configuration is a sum over nearest-neighbour pairs ⟨i,j⟩, with an optional external magnetic field h:

E = −J · Σ_⟨i,j⟩ s_i·s_j  −  h · Σ_i s_i

The coupling constant J sets the sign of the preference. For the ferromagnetic case J > 0, aligned neighbours (s_i·s_j = +1) lower the energy, so the system wants order. For the antiferromagnetic case J < 0, neighbours prefer to point opposite ways. That is the entire model — and from these two lines of arithmetic emerges one of the deepest phenomena in physics: a phase transition.

How it works — energy versus entropy

At temperature T the probability of any configuration follows the Boltzmann distribution:

P(config) = exp(−E/k_B·T) / Z,   Z = Σ_all configs exp(−E/k_B·T)

The whole drama is a tug-of-war between two terms:

• Energy wants order. Aligning every spin minimizes E. At T = 0 the system sits in a perfectly aligned ground state — all up or all down.
• Entropy wants chaos. There are vastly more disordered configurations than ordered ones, and at high T the Boltzmann factor flattens, so the system samples them democratically. The net magnetization averages to zero.

The free energy F = E − T·S decides the winner. At low T the energy term dominates and the lattice is magnetized; at high T the −T·S term dominates and it is disordered. In two or more dimensions these two regimes are separated by a sharp critical temperature T_c where the system switches abruptly from one to the other. That switch is the phase transition, and the magnetization M = ⟨s⟩ is its order parameter: nonzero below T_c, exactly zero above it.

Onsager's exact solution (2D)

The 1D model was Ernst Ising's 1925 thesis problem. He solved it and found, to his disappointment, that it never orders at finite temperature — and wrongly extrapolated that conclusion to all dimensions. It took until 1944 for Lars Onsager to crack the 2D square lattice with a tour-de-force transfer-matrix calculation, proving that it does undergo a transition. The critical temperature is exact:

k_B·T_c = 2J / ln(1 + √2)   ⟹   T_c ≈ 2.269 J/k_B

Onsager (and later C. N. Yang) also derived the spontaneous magnetization below T_c in closed form:

M(T) = [ 1 − sinh⁻⁴(2J/k_B·T) ]^(1/8)   for T < T_c

The exponent 1/8 is the famous critical exponent β for the 2D Ising universality class. Near the transition the specific heat diverges logarithmically — a much gentler singularity than the older mean-field theories predicted, and exactly what experiments on real systems show. This was the first rigorous proof that a microscopic model could produce a genuine thermodynamic singularity, and it reshaped the entire field of critical phenomena.

Why dimension decides everything

The single most surprising fact about the Ising model is that the same energy rule behaves completely differently depending on how many dimensions the lattice lives in. The deciding question is: how much does a domain wall — the boundary between an up-region and a down-region — cost?

Dimension	Transition?	T_c (units of J/k_B)	Why
1D chain	No (only at T = 0)	0	A wall costs fixed 2J but gains k_B·T·ln(N) entropy — always wins for large N
2D square	Yes (continuous)	2.269 (exact, Onsager 1944)	Wall energy grows with lattice size; order survives finite T
2D triangular	Yes	≈ 3.641	More neighbours (6) → stronger ordering tendency
2D honeycomb	Yes	≈ 1.519	Fewer neighbours (3) → weaker ordering
3D cubic	Yes (continuous)	≈ 4.511 (numerical)	Not exactly solvable; exponents from Monte Carlo + bootstrap
4D+ / mean field	Yes	≈ 6.0 (cubic) → ∞ as d grows	Above 4D, mean-field exponents become exact

In 1D, a single domain wall costs only the fixed energy 2J, but it can sit at any of the N − 1 bond positions, contributing entropy k_B·ln(N − 1). The free-energy change ΔF = 2J − k_B·T·ln(N − 1) is negative for any T > 0 once N is large — so walls proliferate and shred any long-range order. That is why the 1D Ising model has no finite-temperature transition. In 2D the energy of a wall instead scales with its length, which grows with system size, so order can win at finite T. This dimensional cliff is known as the Peierls argument.

Worked example — a single Metropolis flip

Take a 2D square lattice with J = 1 and set k_B = 1 so temperatures are measured directly in units of J. Onsager's critical point is then T_c ≈ 2.269. Suppose we pick a spin that is currently +1 and its four neighbours are +1, +1, −1, +1 (three agree, one disagrees). The local energy contribution before the flip is:

E_before = −J · s · Σ neighbours
         = −1 · (+1) · (+1 + 1 − 1 + 1) = −2

If we flip the spin to −1, the sign reverses: E_after = +2. The energy cost of the flip is:

ΔE = E_after − E_before = 2 − (−2) = 4J = 4

Metropolis accepts an energy-raising flip with probability exp(−ΔE/T). Compare two temperatures:

• Cold (T = 1.0, below T_c): exp(−4/1.0) = exp(−4) ≈ 0.018. The flip is rejected ~98% of the time. The lattice stays ordered.
• Hot (T = 4.0, above T_c): exp(−4/4.0) = exp(−1) ≈ 0.368. The flip succeeds more than a third of the time. Thermal noise scrambles the order.

Now run this single rule over millions of sites and sweeps and the macroscopic magnetization emerges. At T = 1.0 the net magnetization settles near |M| ≈ 0.998 (almost perfectly aligned); at T = 4.0 it fluctuates around 0; and right at T = 2.269 the lattice shimmers with domains of every size — the fingerprint of a critical point.

JavaScript — a 2D Ising Monte Carlo

// Metropolis Monte Carlo for the 2D square-lattice Ising model.
// Units: J = 1, k_B = 1, so temperature T is in units of J/k_B.
// Onsager critical point: T_c = 2 / Math.log(1 + Math.SQRT2) ≈ 2.269.

function makeLattice(N) {
  // Random ±1 spins (a "hot start")
  return Array.from({ length: N }, () =>
    Array.from({ length: N }, () => (Math.random() < 0.5 ? 1 : -1))
  );
}

function magnetization(s) {
  let sum = 0;
  for (const row of s) for (const v of row) sum += v;
  return sum / (s.length * s.length); // ⟨s⟩ ∈ [−1, 1]
}

// One Metropolis sweep = N² attempted single-spin flips.
function sweep(s, T) {
  const N = s.length;
  for (let k = 0; k < N * N; k++) {
    const i = (Math.random() * N) | 0;
    const j = (Math.random() * N) | 0;
    // Sum of the four neighbours with periodic (wrap-around) boundaries
    const nb =
      s[(i + 1) % N][j] +
      s[(i - 1 + N) % N][j] +
      s[i][(j + 1) % N] +
      s[i][(j - 1 + N) % N];
    const dE = 2 * s[i][j] * nb; // ΔE = 2·J·s·Σneighbours, with J = 1
    if (dE <= 0 || Math.random() < Math.exp(-dE / T)) {
      s[i][j] = -s[i][j]; // accept the flip
    }
  }
}

const Tc = 2 / Math.log(1 + Math.SQRT2); // 2.2691853...

function run(N, T, sweeps) {
  const s = makeLattice(N);
  for (let t = 0; t < sweeps; t++) sweep(s, T);
  return Math.abs(magnetization(s));
}

// Below, at, and above the critical temperature
console.log("T = 1.5  |M| ≈", run(64, 1.5, 400).toFixed(3)); // ~0.95+ (ordered)
console.log("T = Tc   |M| ≈", run(64, Tc, 400).toFixed(3));  // intermediate, large fluctuations
console.log("T = 3.5  |M| ≈", run(64, 3.5, 400).toFixed(3)); // ~0.05  (disordered)

The whole physics lives in one line: dE = 2 * s[i][j] * nb and the acceptance rule Math.random() < Math.exp(-dE / T). Cool the system slowly and you will watch |M| jump from near-zero to near-one as T drops through 2.269.

Variants and regimes

• Ferromagnetic (J > 0). Neighbours align; the ground state is uniform. The default case and the one Onsager solved.
• Antiferromagnetic (J < 0). Neighbours anti-align; the ground state is a checkerboard. On frustrated lattices (e.g. the triangular antiferromagnet) no configuration satisfies every bond — a frustration that suppresses ordering entirely in 2D.
• External field (h ≠ 0). Breaks the up/down symmetry and smears the transition. Adding a field makes even the 2D model non-exactly-solvable.
• Potts model. Generalizes spins to q states instead of 2; q = 2 reduces to Ising. For q > 4 in 2D the transition turns first-order.
• XY and Heisenberg models. Replace the ±1 spin with a continuous angle (XY) or a 3D unit vector (Heisenberg). Continuous symmetry forbids ordering in 2D (Mermin–Wagner theorem) — the Ising discreteness is what lets 2D order at all.
• Spin glass (Edwards–Anderson). Random ± couplings J_ij; a notoriously hard, frustrated, disordered cousin central to optimization and neural-network theory.

Performance — critical slowing down

Single-spin Metropolis is cheap per step but breaks down precisely where you most want data: near T_c. As T → T_c the correlation length ξ diverges, the lattice fills with domains of every size, and successive configurations become highly correlated. The number of sweeps needed to produce a statistically independent sample — the autocorrelation time τ — scales as a power of the correlation length:

τ ∝ ξ^z,   z ≈ 2.17  (2D single-spin Metropolis)

This is critical slowing down. Because ξ can grow as large as the lattice size L at the critical point, τ ∝ L^z — so doubling the lattice can more than quadruple the time to equilibrate. The cure is a cluster algorithm that flips an entire correlated region in a single move:

Algorithm	Update	Dynamical exponent z (2D)	Best for
Metropolis	One spin, accept/reject	≈ 2.17	Away from T_c; simplest to code
Glauber / heat-bath	One spin, resampled from local field	≈ 2.17	Dynamics studies
Swendsen–Wang (1987)	All clusters at once	≈ 0.25	Equilibrium near T_c
Wolff (1989)	One random cluster	≈ 0.25	High precision near T_c
Worm (directed loop)	Bond/world-line update	≈ 0.25	Correlation functions
Parallel tempering	Swap replicas at different T	—	Spin glasses, rough landscapes

Wolff and Swendsen–Wang cut z from ~2.17 to ~0.25, which near the critical point can mean the difference between minutes and weeks of compute. They are the workhorses of any serious finite-size scaling study.

Where the Ising model shows up

• Ferromagnetism. Its original and namesake application — uniaxial magnets that spontaneously magnetize below the Curie temperature.
• Liquid–gas critical point. The lattice-gas mapping (occupied/empty ↔ up/down) puts simple fluids in the same universality class as the Ising magnet — they share critical exponents.
• Binary alloys. Order–disorder transitions in mixtures like brass map directly onto Ising spins (atom A / atom B).
• Neuroscience & AI. The Hopfield network and Boltzmann machine are Ising models in disguise; neurons are spins and synaptic weights are couplings J_ij.
• Image processing. Ising-style priors denoise binary images by penalizing mismatched neighbours.
• Social and economic models. Opinion dynamics, voter models, and market herding behaviour borrow the spin-alignment metaphor.
• Quantum annealing & optimization. Many NP-hard problems (max-cut, partitioning) map onto finding the Ising ground state — the basis of Ising-machine hardware and D-Wave quantum annealers.

Common pitfalls and misconceptions

• Believing Ising's own conclusion. Ising solved only the 1D chain and declared the model never orders. He was wrong — in 2D and up it absolutely does. Dimension is everything.
• Confusing "exactly solved" with "solved everywhere." Onsager solved only the zero-field 2D case. Turn on a field, or go to 3D, and there is no closed form — you need numerics or the conformal bootstrap.
• Using ⟨s⟩ instead of ⟨|s|⟩ in finite simulations. A finite lattice below T_c flips between the all-up and all-down states over long times, so the raw average magnetization is ~0. Track ⟨|M|⟩ or ⟨M²⟩ to see the order.
• Forgetting periodic boundary conditions. Free edges introduce surface effects that distort T_c badly on small lattices; wrap the lattice into a torus.
• Sampling near T_c with Metropolis and trusting the error bars. Critical slowing down makes naive error bars far too optimistic — successive samples are correlated. Use cluster updates or correct for the autocorrelation time.
• Reading T_c off a single finite lattice. The apparent transition shifts and rounds with system size. The true T_c emerges only from finite-size scaling (e.g. the crossing of the Binder cumulant across lattice sizes).