The Sandpile Model - Part 2

Intro

In this workshop, we continue our exploration of the sandpile model. Unlike last time, where we started from a preconfigured initial state and focused on the final configuration, here we capture and visualize the system’s full evolution over time.

Let’s recap the rules of the sandpile model (with one adjustment: grains are dropped one by one):

  1. We have an empty grid (for example, 51x51) where each cell can hold grains of sand.
  2. We drop grains onto the grid, typically at the center cell.
  3. If any cell has 4 or more grains, it topples, distributing grains to its neighbors (up, down, left, and right).
  4. This can trigger a chain reaction of topplings until all cells are stable (below the threshold of 4 grains).

We run this process for a specified number of drops. To visualize the dynamics, we capture the grid state at regular drop intervals, together with the full toppling dynamics of each sampled drop (the avalanche).

This visualization highlights self-organized criticality (SOC): the system naturally evolves toward a critical state in which events of many sizes can occur, from tiny relaxations to system-spanning avalanches.

Simulation

Last time, we illustrated two sandpile simulation algorithms:

  • a 4-grain-at-a-time toppling method, where in each wave all unstable cells topple by 4 grains and distribute grains to neighbors
  • a more efficient multi-4-grain-at-a-time toppling method, where in each wave all unstable cells topple by multiples of 4 grains at once

We chose the efficient approach because we were interested in the final state only.

This time, because we want to capture the full dynamics, we use the 4-grain-at-a-time toppling algorithm. It lets us track each wave and its effect on the grid during every sampled avalanche.

Along the way, we also record avalanche statistics (topples, waves, avalanches, and boundary loss). These statistics will support later analysis of SOC behavior.

We begin with topple_wave(), which applies one wave across the whole grid: all currently unstable cells topple once. The function returns the number of topples in that wave and the number of grains lost at the boundaries.

We then define run_sandpile_experiment(), which runs the full simulation by dropping grains one by one and capturing grid states and avalanche statistics along the way.

We track the following primary statistics:

  • boundary grain loss (grains that fall off the grid)
  • avalanche sizes in topples (one value per avalanche)
  • avalanche sizes in waves (one value per avalanche)

From the two avalanche size arrays, we derive:

  • number of topples: sum of avalanche sizes in topples
  • number of waves: sum of avalanche sizes in waves
  • number of avalanches: length of the avalanche arrays
  • largest avalanche topples count: max of avalanche sizes in topples
  • largest avalanche waves count: max of avalanche sizes in waves

In addition to numpy, we use xarray for labeled multidimensional data handling. xarray.DataArray lets us store grid snapshots and related metadata in a structured form, which makes downstream analysis and visualization easier.

Now we are ready to run the simulation and collect the results. We set the grid size, total number of grains to drop, and snapshot interval, then call run_sandpile_experiment() to generate the data used for analysis and visualization.

Visualization

Evolution of the sandpile

We visualize the sandpile’s evolution by plotting the grid state dynamics at each sampled drop. Each plot shows the distribution of grains across the grid, with color intensity representing grain count. By arranging these plots in a sequence, we can see how the sandpile grows and evolves over time, including the formation of avalanches and the emergence of complex patterns.

Among the sampled drops, we observe many small avalanches and a few large ones. In particular, drops 3000 and 9000 trigger large avalanches that spread across much of the grid, while most other sampled drops trigger small avalanches near the center. This is a hallmark of self-organized criticality (SOC): the system naturally evolves to a critical state in which events of many sizes can occur.

Power law distribution of avalanche sizes

The size distribution of avalanches in the sandpile model follows a power-law form. A power law means the probability of observing an avalanche of size \(s\) scales as \(P(s) \sim s^{-\tau}\) for some exponent \(\tau\). Unlike a normal distribution, where values cluster around a typical scale (i.e., its mean), a power law has a long tail: small events are common, while very large events are rare but important.

A standard way to visualize this is a log-log plot of avalanche-size probability. If the distribution is close to a power law over a range of sizes, the points in that range appear approximately linear, and the slope is related to the exponent \(\tau\).

In the plot below, size is defined as the number of topples per avalanche. The x-axis is avalanche size, and the y-axis is the estimated probability of observing that size; both axes are shown on logarithmic scales.

The scatter points follow an approximately straight, downward trend over part of the range, consistent with power-law-like behavior: small avalanches are frequent, while large avalanches are rare but still present.

If we plot avalanche sizes in waves instead of topples, we should see a similar heavy-tailed pattern, with a potentially different exponent \(\tau\).

Exercise

  • Experiment on other grid sizes and total number of grains. What interesting dynamics do you observe? How does the grid size affect the distribution of avalanche sizes? (Note that setting a larger grid size and a larger number of grains will require more computational resources and time to run the simulation.)

  • Plot the distribution of avalanche sizes in waves and compare it to the distribution of avalanche sizes measured in topples. Does it also follow a power law distribution? If so, what are the estimated exponents \(\tau\) for each distribution?

References

  1. Sandpile Model - Wikipedia
  2. Self-Organized Criticality - Wikipedia
  3. Power Law - Wikipedia
  4. Sandpile Avalanche Simulation - Veritasium and its associated video about power law and self-organized criticality: You’ve (Likely) Been Playing The Game of Life Wrong