Ising on the Cake

Subtitle: Going through a difficult phase transition

The Ising Model


Wilhelm Lenz was a professor of theoretical physics at the University of Hamburg from 1921 to 1956. Ernst Ising was a student of Lenz’ and solved what is now known as the one-dimensional Ising Model in his thesis Beitrag zur Theorie des Ferromagnetismus published in 1925. The model represents magnetic dipole moments of atomic spin that may be in one of two states, up and down, and modeled as +1+1 or 1-1 over a grid of dipoles. Lenz invented the concept in 1920, but it was Ising who first solved it. Later, Lars Onsager solved the two-dimensional version in 1944 and P. F. Barth solved the problem for arbitrarily large dimensions in 1981.

Imagine you have a grid of tiny magnets, each one pointing either up or down. These magnets represent the “spins” of particles in the material. In the Ising Model, each magnet interacts with its neighbors, meaning it wants to align its direction with them. The key idea is that each magnet has a preference to either align with its neighbors (trying to lower its energy) or go against them (which increases its energy). This preference is determined by the temperature of the system and the strength of the interaction between magnets.

At higher temperatures, the magnets are more likely to be randomly oriented because they have enough energy to overcome the interaction between them. As the temperature decreases, the magnets start to align more and more, leading to what we call “order” in the system. This is similar to how ice crystals form when water cools down.

The Ising Model helps scientists understand phase transitions, where a material changes from one state to another, like from a magnetized state to a non-magnetized state (or vice versa). It’s also used in various fields beyond physics, such as biology and computer science, to model systems with interacting components.

Details of the Model


Schematic representation of a configuration of the 2D Ising model on a square lattice.

In the schematic shown above, dipoles spinning in the +1+1 (up) direction are shown in red, while the 1-1 (down) dipoles are in blue. Interactions occur between the four nearest neighbors connected by the gray gridlines. The spin at node ii is σi\sigma_i and the energy of the entire grid is

E(σ)=J<ij>σiσjhiσiE({\sigma}) = -J \sum_{<ij>} \sigma_i \sigma_{j} - h \sum_i \sigma_i

where <ij><ij> is the interaction between nearest neighbors on the grid. For example, a down spin surrounded by four up spins (either of the two blue spins in the second line from the bottom) would have a sum of 4J4J for the first term. The second term is the sum of all the spins, so in this example, there are 1414 up spins and 1111 down spins for a total of 3h-3h. The parameter JJ measures the strength of the coupling between adjacent magnetic dipoles, and hh is the effect of any external magnetic field.

The Metropolis Algorithm

The Metropolis (or Metropolis-Hastings) Algorithm is a Markov Chain Monte Carlo (MCMC) method that generates a sequence of random numbers from a probability distribution. A Markov Chain is a sequence of states where the current state depends only on the previous state, and the Monte Carlo method is a way to create random samplings.

The Metropolis algorithm is a computational method used to simulate the behavior of systems described by statistical mechanics, such as the Ising Model. It’s particularly useful for systems with many interacting components, where exact solutions are difficult or impossible to obtain.

Here’s how the Metropolis algorithm works within the context of the Ising Model:

  1. Initialize the system: Start with an initial configuration of spins (magnets) on a grid. Each spin can be randomly assigned an orientation, either up or down.
  2. Select a spin to flip: Randomly choose one spin from the grid. This spin will be considered for flipping its orientation.
  3. Calculate energy change: Compute the change in energy that would result from flipping the selected spin. This involves considering the interactions between the selected spin and its neighboring spins according to the rules of the Ising Model.
  4. Accept or reject the flip: If the energy change is negative (i.e., flipping the spin lowers the total energy of the system), accept the flip. If the energy change is positive, accept the flip with a probability determined by the Boltzmann factor: eΔE/kTe^{-\Delta E / kT}, where ΔE\Delta E is the energy change, TT is the temperature of the system, and kk is the Boltzmann constant. This probabilistic acceptance allows the system to explore higher energy configurations temporarily, which is essential for simulating systems at finite temperatures.
  5. Update the system: If the flip is accepted, update the spin configuration accordingly. If it’s rejected, keep the current configuration.
  6. Repeat: Repeat steps 2-5 for a large number of iterations or until the system reaches equilibrium, where its properties stabilize.

By iterating through these steps, the Metropolis algorithm allows us to simulate the behavior of the Ising Model at different temperatures and study phenomena such as phase transitions. It’s a powerful tool for understanding the collective behavior of many interacting spins in a material.

A NetLogo Simulation

The Ising Model can be simulated in NetLogo either by downloading NetLogo or running the model online. In the code section of NetLogo, you might want to change the color of the up spin from blue + 2 to red to match the colors of the schematic. See the yellow highlighted line below.


Ising NetLogo Code

To run the simulation, choose an initial temperature and click setup-random, then go. Watch the “Magnetization” graph evolve as spins flip from +1+1 to 1-1 or vice-versa. Next, try starting with all spins in the same direction, setup -1 or setup +1.


Ising Model Simulation

Initially, the average spin is +1+1 because all dipoles are up, but over time random flips bring the average back to near 00. What happens if the temperature is set very high or very low?

The Ising Model for Phase Transitions

The Ising Model can be used to simulate phase transitions like the one from ice to water or water to steam. Click on setup -1 to orient all of the spins in the same direction, which could be interpreted as all of the molecules locked in the same “ice” condition. Hydrogen bonding between adjacent molecules provides an energy bonus that tends to lock them in a single orientation, the first term of E(σ)E(\sigma).

At higher temperatures, the molecules are more free to reorient as they transition from the ice phase to the water phase. A larger value for the parameter Ediff leads to a sharper transition between phases. Try plotting the sum of all grid couplings (the sum of Ediff for all molecules):

to-report energy
let total-energy 0
ask patches [
let Ediff 2 * molecule-orientation * sum [molecule-orientation] of neighbors
set total-energy total-energy + Ediff
report total-energy

You can modify the model by adding a new plot for total-energy. For a phase transition, rather than ‘Magnetization’ you could rename the plot to ‘Order Parameter’ which measures the average angle between each molecule’s orientation and a reference angle (like a crystal axis). Lower values indicate more alignment. Change ‘average spin’ to ‘average alignment’ or something similar.

Jakub Nowosad is a computational geographer working at the intersection between geocomputation and the environmental sciences, and is an Assistant Professor at Adam Mickiewicz University in Poland. He used the Ising Model to simulate the transition from forest to open areas in two recently published articles, Simulating spatial patterns with the spatial kinetic Ising model and Optimizing the parameters of the spatial kinetic Ising model to simulate spatial patterns. He used R rather than NetLogo for his simulations, but the results are similar.

The Ising Model is an easily understood representation of many real-world interactions where nearest-neighbor states influence surrounding states. It may be useful for improving global warming models such as described in the paper Ising model for melt ponds on Arctic sea ice, and who knows, maybe you can use it to solve this difficult maze problem.


Maze Proof Establishes a ‘Backbone’ for Statistical Mechanics

Image credits

Hero: Ivy Sheet Cake, Homestyle Bakery Nashville.

Wilhelm Lenz: Wilhelm Lenz at the Institute of Physics

Ernst Ising: Biography of Ernst Ising

Schematic representation of a configuration of the 2D Ising model on a square lattice. Sascha Wald

Maze Proof Establishes a ‘Backbone’ for Statistical Mechanics - Quanta Magazine

Code for this article

Ising NetLogo

Further reading

Ising - NetLogo Modeling Commons

Monte Carlo, Metropolis and the Ising Model

Ising model

The Ising model

Optimizing the parameters of the spatial kinetic Ising model to simulate spatial patterns - Jakub Nowosad

Simulating spatial patterns with the spatial kinetic Ising model - Jakub Nowosad

Numerical Solutions to the Ising Model using the Metropolis Algorithm - Danny Bennett

Phase diagram of the square 2D Ising lattice with nearest neighbor and next-nearest neighbor interactions - Harold J. W. Zandvliet

Metropolis–Hastings algorithm - Wikipedia

The Ising Model - Jeffrey Chang, Notes on Statistical Mechanics II

Ising model for melt ponds on Arctic sea ice - Yi-Ping Ma, Ivan Sudakov, Courtenay Strong and Kenneth M Golden, New Journal of Physics, 2019


NetLogo is a multi-agent programmable modeling environment.

Posts using NetLogo


RStudio is an integrated development environment (IDE) for R. It includes a console, syntax-highlighting editor that supports direct code execution, as well as tools for plotting, history, debugging and workspace management.

Posts using RStudio

See all software used on wildpeaches →