The Game of Life, Inheritance, and Cellular Automata¶
Preamble: Run the cells below to import the necessary Python packages
This notebook created by William Gilpin. Consult the course website for all content and GitHub repository for raw files and runnable online code.
import numpy as np
# Import local plotting functions and in-notebook display functions
import matplotlib.pyplot as plt
%matplotlib inline
The Game of Life and Cellular Automata¶
- Cellular automata are systems where all physical laws and fields are discretized: space is a lattice, and time is discrete
- Physical laws take the form of IF-THEN statements
- Basic cellular automata are Markovian: given current field values at some set of lattice points, perform an instantaneous update based on a set of rules
The Game of Life¶
- Conway, 1970
- Complex dynamics emerge from simple rules
- Markovian: next state of universe depends only on current state of universe
- Local: next state at a lattice point depends only on current site and nearest neighbors
- Synchronous: we update all sites at once, rather than raster scanning
- Finite fields: A site can have values of 1 (alive) or 0 (dead)
Game life rules:¶
- Underpopulation: Any live cell with fewer than two live neighbours dies
- Survival: Any live cell with two or three live neighbours lives stays alive
- Overpopulation: Any live cell with more than three live neighbours dies
- Reproduction: Any dead cell with exactly three live neighbours becomes a live cell
These simple rules give rise to surprisingly complex behaviors; in fact, the Game of Life has been shown to support universal computation, given a suitable encoding scheme based on initial conditions. We've included some examples of fixed-point, limit-cycle, and soliton like dynamics below
.
Images from Scholarpedia
Implementing Cellular Automata using inheritance in Python¶
Below, we are going to define a "base" class for cellular automata, and then define subclasses for various specific automata. The base class handles state counting, initialization, and will contain a simulate method for the simulation loop (where we repeatedly call the update rule).
In the base class, we will define a next_state method because the simulate method depends on it. But because the state rule is different for each automaton, we will leave the implementation of next_state to be overwritten by the subclasses (the base class raises an exception)
This convention is an example of the Template Method Pattern in software development.
# A base allows us to impose structure on specific cases. We will use this as a parent
# for subsequent classes
class CellularAutomaton:
"""
A base class for cellular automata. Subclasses must implement the step method.
Parameters
n (int): The number of cells in the system
n_states (int): The number of states in the system
random_state (None or int): The seed for the random number generator. If None,
the random number generator is not seeded.
initial_state (None or array): The initial state of the system. If None, a
random initial state is used.
"""
def __init__(self, n, n_states, random_state=None, initial_state=None):
self.n_states = n_states
self.n = n
self.random_state = random_state
np.random.seed(random_state)
## The universe is a 2D array of integers
if initial_state is None:
self.initial_state = np.random.choice(self.n_states, size=(self.n, self.n))
else:
self.initial_state = initial_state
self.state = self.initial_state
self.history = [self.state]
def next_state(self):
"""
Output the next state of the entire board
"""
raise NotImplementedError
def simulate(self, n_steps):
"""
Iterate the dynamics for n_steps, and return the results as an array
"""
for i in range(n_steps):
self.state = self.next_state()
self.history.append(self.state)
return self.state
ca = CellularAutomaton(10, 2, random_state=1)
ca.simulate(10)
# plt.imshow(ca.initial_state)
--------------------------------------------------------------------------- NotImplementedError Traceback (most recent call last) Cell In[28], line 49 45 return self.state 48 ca = CellularAutomaton(10, 2, random_state=1) ---> 49 ca.simulate(10) 50 # plt.imshow(ca.initial_state) Cell In[28], line 43, in CellularAutomaton.simulate(self, n_steps) 39 """ 40 Iterate the dynamics for n_steps, and return the results as an array 41 """ 42 for i in range(n_steps): ---> 43 self.state = self.next_state() 44 self.history.append(self.state) 45 return self.state Cell In[28], line 36, in CellularAutomaton.next_state(self) 32 def next_state(self): 33 """ 34 Output the next state of the entire board 35 """ ---> 36 raise NotImplementedError NotImplementedError:
Just to show how inheritance works, we will start by defining two simpler automata
An extinction automaton that sends every single input condition to the same output state
A random automaton that assigns each input state to either 0 or 1 randomly
class ExtinctionAutomaton(CellularAutomaton):
"""
A cellular automaton that simulates the extinction of a species.
"""
def __init__(self, n, **kwargs):
super().__init__(n, 2, **kwargs)
def next_state(self):
"""
Output the next state of the entire board
"""
next_state = np.zeros_like(self.state)
return next_state
class RandomAutomaton(CellularAutomaton):
"""
A cellular automaton with random updates
"""
def __init__(self, n, random_state=None, **kwargs):
np.random.seed(random_state)
super().__init__(n, 2, **kwargs)
self.random_state = random_state
def next_state(self):
"""
Output the next state of the entire board
"""
next_state = np.random.choice(self.n_states, size=(self.n, self.n))
return next_state
Let's try running these automata, to see how they work
model = ExtinctionAutomaton(40, random_state=0)
model.simulate(200)
plt.figure()
plt.imshow(model.initial_state, cmap="gray")
plt.title("Initial state")
plt.figure()
plt.imshow(model.state, cmap="gray")
plt.title("Final state")
Text(0.5, 1.0, 'Final state')
model = RandomAutomaton(40, random_state=0)
model.simulate(200)
plt.figure()
plt.imshow(model.initial_state, cmap="gray")
plt.title("Initial state")
plt.figure()
plt.imshow(model.state, cmap="gray")
plt.title("Final state")
Text(0.5, 1.0, 'Final state')
We are now ready to define the Game of Life as a subclass of the base class.
Game life rules:¶
- Underpopulation: Any live cell with fewer than two live neighbours dies
- Survival: Any live cell with two or three live neighbours lives stays alive
- Overpopulation: Any live cell with more than three live neighbours dies
- Reproduction: Any dead cell with exactly three live neighbours becomes a live cell
# This child class inherits methods from the parent class
class GameOfLife(CellularAutomaton):
"""
An implementation of Conway's Game of Life in Python
Args:
n (int): The number of cells in the system
**kwargs: Additional keyword arguments passed to the base CellularAutomaton class
"""
def __init__(self, n, **kwargs):
# the super method calls the parent class's __init__ method and passes the
# arguments to it.
super().__init__(n, 2, **kwargs)
def next_state(self):
"""
Compute the next state of the lattice
"""
# Compute the next state
next_state = np.zeros_like(self.state) # preallocate the next state
for i in range(self.n):
for j in range(self.n):
# Count the number of neighbors
n_neighbors = np.sum(
self.state[i - 1:i + 2, j - 1:j + 2]
) - self.state[i, j]
# Update the next state
if self.state[i, j] == 1:
if n_neighbors in [2, 3]:
next_state[i, j] = 1
else:
next_state[i, j] = 0
else:
if n_neighbors == 3:
next_state[i, j] = 1
else:
next_state[i, j] = 0
return next_state
We can now run this automaton and see how it behaves. We will also try using the timing utility to see how long it takes to run the simulation for one step.
model = GameOfLife(40, random_state=0)
model.simulate(500)
# print timing results
print("Average time per iteration: ", end="")
%timeit model.next_state()
plt.figure()
plt.imshow(model.initial_state, cmap="gray")
plt.title("Initial state")
plt.figure()
plt.imshow(model.state, cmap="gray")
plt.title("Final state")
Average time per iteration: 3.22 ms ± 20.4 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
Text(0.5, 1.0, 'Final state')
We can also access the time series of automaton states, which our base class stores in a list. We can use this to make animations, and to analyze the behavior of the automaton over time.
from matplotlib.animation import FuncAnimation
from IPython.display import HTML
# Assuming frames is a numpy array with shape (num_frames, height, width)
frames = np.array(model.history).copy()
fig = plt.figure(figsize=(6, 6))
img = plt.imshow(frames[0], vmin=0, vmax=1, cmap="gray");
plt.xticks([]); plt.yticks([])
# tight margins
plt.margins(0,0)
plt.gca().xaxis.set_major_locator(plt.NullLocator())
def update(frame):
img.set_array(frame)
ani = FuncAnimation(fig, update, frames=frames, interval=50)
HTML(ani.to_jshtml())