Visual Bak-Tang-Wiesenfeld Sandpile Model for Matlab

Matlab library to visualize the avalanche and power law characteristics of the sandpile model Edit

This Matlab library helps to visualize the avalanche and power law characteristics of the popular sandpile model.
The library is meant for visualizing the Bak-Tang-Wiesenfeld sandpile model and plotting the sandpile's statistics on the fly. The visualization has been optimized for performance.

--- How This Sandpile Model Works ---
In this sandpile model, the sandpile is built onto a square grid of a specified "pile_width" The total number of positions in the square grid is (pile_width * pile_width), so for a sandpile with a "pile_width" of 12, 144 grid positions are available.

On every grid position, an arbitrary amount of sand grains can be stacked. However, if the number of sand grains on a grid position is greater than or equal to 4, a toppling event takes place, originating at that grid position.

In a toppling event, four grains get taken away from the originating grid position and are distributed equally among the neighboring grid positions north, east, south, and west of the originating grid position. As a result, the originating grid position possesses four grains less, while the number of grains on all its four neighboring positions has increased by one. This behavior eventually leads to a chain reaction, in which a further toppling event on one of the neighboring positions is caused. Such a chain reaction is called "avalanche". If, through a toppling event, a sand grain is to be added to a position beyond the boundary of the pile, it is ignored and no longer considered in the model. This can be imagined as a sand grain "falling off" the surface on which the sandpile is created.

In this model, new sand grains are added to random grid positions. The sand grains are added one by one, except during toppling or avalanche events. In this case, no grains are added until these events are over. The model is initialized with a random number of up to 3 grains at every grid position.

--- Usage ---
The code is thoroughly documented, the documentation can be accessed through the "doc simulateSandpile" command. To run the model, call the "simulateSandpile" function, the only user-facing function in the script.

Configurable input parameters of "simulateSandpile" are:

- "pile_width": Side length of the square pile

- "no_of_grains": Number of grains that are added one by one to the sandpile as the model runs. The simulation stops after the specified number of grains has been added.

- "draw_speed": Speed of animation of the plots. A speed of 0 skips the sandpile animation entirely, and will only plot a chart of avalanche sizes once the simulation has finished running. A value of 0.5 yields a relatively slow animation, smaller values yield a faster one.

The following code is a simple usage example:

pile_width = 40; % Pile represented through 40x40 grid
no_of_grains = 4000; % Run the simulation until 4000 new sand grains have been added
draw_speed = 0.2; % Animate the pile relatively slow
avalanche_output = simulateSandpile(pile_width, no_of_grains, draw_speed);

The "simulateSandpile" function returns a single output variable, "avalanche_output". This variable is a matrix with the shape (max. avalanche size, 2). The matrix contains information about how many avalanches with a specific size occurred during the run.

Here is an annotated example of a typical "avalanche_output":

[1 50; % 50 avalanches causing 1 toppling event
2 32; % 32 avalanches causing 2 toppling events
3 11; % 11 avalanche causing 3 toppling events
...];