Day 7 – Automatic on a Cellular Level

Today’s advent calendar post deals with Cellular Automata.

What are cellular automata? I’m glad you asked! They are systems made up of only a few parts: A kind of field or “world” made up of cells, a set of states that every cell can be in at any point, a “neighborhood” describing what cells are visible to each cell, and a set of rules governing what a cell will change its state into given its own state and the state of all cells in its neighborhood.

That is, of course, a very abstract description, so let me give a few examples for the individual parts to hopefully give you an idea of what you might see in a cellular automaton:

In typical worlds you might find cells arranged like beads on a string, or like fields on a chess or chinese checkers board. There are also more exotic configurations you can make up: Any 2-dimensional field can be mapped onto any surface, like for example the stanford bunny

Sets of states you can find out in the wild would be “the numbers from 0 through n”, “there are bacteria here or not”, “the colors black and white” (or more). Since you can represent basically any piece of information as “a number”, and any number of states is allowed, there can also be states representing “how many particles in this cell are going up, down, left, or right at this moment?” as integer or even floating point numbers.

A neighborhood can be thought of as a pattern in which cells are “wired together”. Typical neighborhoods would be “the one in front and the one behind” for the “string of beads” field, and “north, east, south, west” or “north, northeast, east, southeast, south, southwest, west, northwest” for the chess board field – these two are the Von Neumann neighborhood and the Moore neighborhood respectively.

The set of rules that govern what states in each cell’s neighborhood will cause what state to go to what other state could be considered the heart of a particular cellular automaton.

In today’s advent calendar post we’ll explore what you might call the simplest automatons. We’ll string the fields up like beads on a string, and we’ll try a neighborhood or two and a few different state sets.

To make things interesting for those at home, I’ll be giving you links that let you run example code right here in your browser, or at home with your local perl6 compiler!

Doing for Learning

Let’s get started, then:

First of all, what do we need? There’ll have to be storage for the world, a bit of code to get a cell’s eligible neighbors, and a bit of code to calculate a cell’s next state given its own state and the states of the neighbors. On top of it all, we’ll want to see what’s going on, so we’d have some code for that, too.

After deciding what states each cell can have, we’ll know what storage is appropriate for our world. Going with an array of 8 bit integers will allow us to choose from any set of states that doesn’t have more than 255 individual ones. Let’s go with 3 states in total for now, though. We can initialize the world any way we please, but setting every field to a random valid state is a good starting point. Another is to put a single cell of one of the states in the middle and have every other cell with a different state.

constant number-of-states = 3;
constant field-width = 60;

my int8 @field;

sub init-field {
    @field = (^number-of-states).roll(field-width);
}

Displaying the field is pretty straight-forward, depending on what output we’re using. Here’s a piece of code that works in any console, and below that a link to a 6pad that outputs little colored squares in the HTML part of the pad.

sub output-a-row {

    for @field {
        # Using the unicode characters "Light shade", "Medium shade", and "Dark shade"
        # and printing each field twice so they look square rather than slim and tall.
        print ["\x2591", "\x2592", "\x2593"][$_] x 2
    }

    say "";
}

init-field;
output-a-row;

Run this code in your browser. You will have to wait a few seconds for the currently rather big perl6 compiler in javascript.

Stepping ahead

Getting from just one single row to a simulated run of the cellular automaton needs a whole bunch of parts at once.

Cellular automata by their definition will advance all of their cells at the same time. We’ll of course not go out to the cloud to get a machine with as many cpu cores as there are cells in our field. We’ll settle for “calculate the next step for each cell based on what their neighbor cells were in the previous step” by going through all fields with a simple loop.

The straight-forward approach for this is to have an extra field to put calculation results and copy the results into the “real” field array after each step. Let’s give that a try.

sub simulate-step {
    my int8 @output;

    for ^field-width -> $x {
        # do some calculations here
    }
    
    @field = @output;
}

Let’s see what we need in order for the calculations to happen: The new state will depend on the neighborhood and the cell itself. We’ll go with possibly the most obvious neighborhood: A cell, it’s predecessor and successor in the field. But wait, what happens with the very first and the very last cell? Let’s just pretend they have an extra neighbor that’s just always at state 0. That way

sub simulate-step {
    my int8 @output;
  
    for ^field-width -> $x {
        my $left   = $x - 1 < 0 ?? 0 !! @field[$x - 1];
        my $middle = @field[$x];
        my $right  = $x + 1 >= field-width ?? 0 !! @field[$x + 1];
        
        # do some calculation with $left, $middle, and $right
        # then push the result into @output
    }
    
    @field = @output;
}

So, we’re finally at the spot where we need to decide what our cellular automaton is actually supposed to do in the first place. But how are we supposed to figure out what it should do when we don’t even have a meaning assigned to states “0, 1, and 2”?

The answer is simple! Literally make it all up.

Making stuff up

What I mean by that is, we currently don’t really care what the cellular automaton does, as long as it looks good. So why not decide up front for every possibility what’s supposed to happen by rolling some imaginary dice?

For this purpose it helps to know how many possible “configurations” there even are for the cell and its neighbors to be in. Fortunately, that’s very simple. You can imagine the three cells as a number made up of three digits, and each digit is allowed to be 0, 1, or 2.

000  001  002
010  011  012
020  021  022
100  101  102
110  111  112
120  121  122
200  201  202
210  211  212
220  221  222

If I didn’t mess up, that’s all the possibilities for left, middle, and right cell to make up. Just like a four digit binary number can be one of 2⁴ numbers, this three digit ternary number can be one of 3³. That means that we just have to pick 3³ numbers between 0 and 2, one for each number in the table up above.

It’s actually pleasantly easy to do this!

my int8 @lookup-table = (^number-of-states).roll(3³);

And given our $left, $middle, and $right variables, we can multiply the first one with 9, the second one with 3, and add the three together to get an index into our lookup table:

sub simulate-step {
    my int8 @output;
  
    for ^field-width -> $x {
        my $left   = $x - 1 < 0 ?? 0 !! @field[$x - 1];
        my $middle = @field[$x];
        my $right  = $x + 1 >= field-width ?? 0 !! @field[$x + 1];
        
        my $index = $left * 9 + $middle * 3 + $right;
        
        @output.push(@lookup-table[$index]);
    }
    
    @field = @output;
}

And running this will already give us a shiny to look at. All we need to do is hook up the subs:

constant number-of-states = 3;
constant field-width = 60;

my int8 @field;

sub init-field {
    @field = (^number-of-states).roll(field-width);
}
init-field;

sub output-a-row {

    for @field {
        # Using the unicode characters "Light shade", "Medium shade", and "Dark shade"
        # and printing each field twice so they look square rather than slim and tall.
        print ["\x2591", "\x2592", "\x2593"][$_] x 2
    }

    say "";
}

my int8 @lookup-table = (^number-of-states).roll(3³);

sub simulate-step {
    my int8 @output;
  
    for ^field-width -> $x {
        my $left   = $x - 1 < 0 ?? 0 !! @field[$x - 1];
        my $middle = @field[$x];
        my $right  = $x + 1 >= field-width ?? 0 !! @field[$x + 1];
        
        my $index = $left * 9 + $middle * 3 + $right;
        
        @output.push(@lookup-table[$index]);
    }
    
    @field = @output;
}

for ^100 {
    simulate-step;
    output-a-row;
}

The results already look pretty interesting some of the time! Of course we need to get lucky with the random lookup table. Here’s one I like to get you started if you get lots of uninteresting ones:

my int8 @lookup-table = <0 0 2 0 0 0 1 2 0 0 1 1 2 1 1 2 1 1 1 0 1 2 2 0 2 1 1>;

And here’s the link to the 6pad where you can try it all in your browser.

Thirdly, here’s a screenshot from my machine in case you’re reading on a mobile device or something else that can’t run perl6 yet.

fish -tmp_112

Alterations

Now that our simulator does what it’s supposed to, let’s have some fun with some tweaks.

First, let’s see what it takes to increase the number of different states:

constant number-of-states = 4;

# the size of the lookup table should be based on the number of states
my int8 @lookup-table = (^number-of-states).roll(number-of-states³);

sub output-a-row {

    for @field {
        # add unicode character "Full block" for the fourth state
        print ["\x2591", "\x2592", "\x2593", "\x2588"][$_] x 2
    }

    say "";
}

And the calculation needs to do its computation based on number-of-states as well:

my $index = $left * number-of-states * number-of-states
            + $middle * number-of-states
            + $right;

And that’s already it! Not even terribly hard so far.

Changing up the neighborhood

Now this one’s much more interesting. Changing the neighborhood will require our calculation loop to get more variables for the index calculation, and the lookup table will change its size once more as well.

Let’s go back to 3 states instead of 4 and replace the neighborhood with one that has just one more cell in it: We’ll take a cell’s predecessor and its successor, but leave out the cell itself. We then add the predecessor’s predecessor and the successor’s successor:

# three states, but four neighbors
constant number-of-states = 3;
constant number-of-neighbors = 4;

# ...

# exponentiate number-of-states with number-of-neighbors, like
# you would to get a number-of-neighbors number in base number-of-states.
my int8 @lookup-table = (^number-of-states).roll(number-of-states ** number-of-neighbors);

sub simulate-step {
   my int8 @output;

   for ^field-width -> $x {
       my $leftleft   = $x <= 1 ?? 0 !! @field[$x - 2];
       my $left       = $x == 0 ?? 0 !! @field[$x - 1];

       my $right      = $x == field-width - 1 ?? 0 !! @field[$x + 1];
       my $rightright = $x >= field-width - 2 ?? 0 !! @field[$x + 2];

       # many multiplications later ...
       my $index = $leftleft * number-of-states * number-of-states * number-of-states
                   + $left   * number-of-states * number-of-states
                   + $right  * number-of-states
                   + $rightright;

       @output.push(@lookup-table[$index]);
   }

   @field = @output;
}

fish  -tmp_113.png

Here’s the 6pad to try it out

Sadly, all it seems to do is make the output a little more chaotic.

Optimization opportunities?

At the moment, the code is a compromise between performant and readable. It could also have looked like this:

for (0, |@field, 0).rotor(3 => -2) -> ($left, $middle, $right) {
    my $index = :3[$right, $middle, $left];
}

Though my intuition tells me that would have been noticeably slower.

But we can make the code a bit faster still, without sacrificing too much readability, even!

There is one thing that our calculation loop is doing way too much of: Array accesses! Every cell gets accessed three times in a row: Once when it becomes the $right, again when it becomes the $middle and another time when it becomes $left.

So how can we do this better? The first thing I thought of was to have the variables $left, $middle, and $right stick around between iterations and shifting the cell values through:

my $left   = 0;
my $middle = @field[0];
my $right  = @field[1];

for ^field-width -> $x {
    my $index = $left * number-of-states * number-of-states
            + $middle * number-of-states
            + $right;

    @output.push: @lookup-table[$index];
    $left = $middle;
    $middle = $right;
    $right = $x + 1 >= field-width ?? 0 !! @field[$x + 1];
}

Cool, we’ve even gotten rid of one of the checks of $x vs field-width! But there’s another thing that happens over and over that we could make just a tiny bit simpler. We can have the $left, $middle, and $right variables already hold the exact value needed for the addition:

my $left   = 0;
my $middle = @field[0] * 3;
my $right  = @field[1];

for ^field-width -> $x {
    my $index = $left + $middle + $right;

    @output.push: @lookup-table[$index];
    $left = $middle * 3;
    $middle = $right * 3;
    $right = $x + 1 >= field-width ?? 0 !! @field[$x + 1];
}

I think that looks pretty neat!

Other changes?

One kind of cellular automaton I’ve encountered is one where every cell has a chance to do its calculation on each step, and otherwise just keeps its state for an extra step. Let’s see how that’s implemented:

constant probability = 0.75e0;

my $left   = 0;
my $middle = @field[0] * 3;
my $right  = @field[1];

for ^field-width -> $x {
    if rand < probability {
        my $index = $left + $middle + $right;

        @output.push: @lookup-table[$index];
    }
    else {
        @output.push: $middle;
    }
    $left = $middle * 3;
    $middle = $right * 3;
    $right = $x + 1 >= field-width ?? 0 !! @field[$x + 1];
}

Lower probabilities are easy to spot, because they will cause the resulting image to look vertically stretched. Higher probabilities can cause completely regular patterns to stay mostly intact, but at some point be broken up in a spot or two.

Here’s a screenshot for you!

fish -tmp_114

Is this useful?

Cellular automata in general are extremely versatile, and even very simple ones can already handle universal computation, like “Rule 110”. There are more complex automatons like Von Neumann’s machines capable of self-replication and WireWorld, which has been used to build a little machine that can calculate prime numbers and display them on a little seven-segment display.

Pretty astonishingly, there’s a turing machine with a literal tape built from the pretty popular Game of Life and, potentially more astonishingly, a configuration of Game of Life that can calculate and display a Game of Life.

All in all, I find cellular automata a fascinating topic. Two-dimensional cellular automata were barely mentioned in this post, but there are many interesting ones besides the ones already mentioned in this section.

Implementation wise, you would most probably not use CPU code to simulate cellular automata. At the very least, you wouldn’t use a loop that goes over every individual cell – see the fantastic HashLife algorithm that chops the world up into bigger and bigger blocks that occur often and does many whole-world steps at once. Otherwise, you would most likely simulate a CA on the GPU, which offers extreme degrees of parallelism when the code run for each cell is the same.

Thank you for staying with me through this really rather long post!

I hope I could awaken even an idle interest in the fantastic and wide world of cellular automata!

Have a lovely december, everyone!