really simple cellular automata to illustrate fire breaks

automata.py

import numpy as np
import matplotlib.pyplot as plt

u = np.array([[4], [2], [1]])

def step(x, rule_b):
    """Compute a single stet of an elementary cellular
    automaton."""
    # The columns contains the L, C, R values
    # of all cells.
    y = np.vstack((np.roll(x, 1), x,
                   np.roll(x, -1))).astype(np.int8)
    # We get the LCR pattern numbers between 0 and 7.
    z = np.sum(y * u, axis=0).astype(np.int8)
    # We get the patterns given by the rule.
    return rule_b[7 - z]

def rand_bin_array(K, N):
    """Compute a random binary array of length N and 
    K ones."""
    arr = np.zeros(N)
    arr[:K]  = 1
    np.random.shuffle(arr)
    return arr

def generate(rule, size=100, steps=100):
    """Simulate an elementary cellular automaton given
    its rule (number between 0 and 255)."""
    # Compute the binary representation of the rule.
    rule_b = np.array(
        [int(_) for _ in np.binary_repr(rule, 8)],
        dtype=np.int8)
    x = np.zeros((steps, size), dtype=np.int8)
    # Random initial state.
    x[0, :] = np.random.rand(size) < .5
    # Apply the step function iteratively.
    for i in range(steps - 1):
        if 40 <= i < 60:
            x[i + 1, :] = step(x[i, :], rule_b)
            x[i + 1, 10:30] = rand_bin_array(12, 20)
            x[i + 1, 70:90] = rand_bin_array(14, 20)
        else:
            x[i + 1, :] = step(x[i, :], rule_b)
    return x
  
x = [generate(154) for _ in range(100)]
y = np.dstack(x)
img = y.mean(axis=2)
plt.figure(figsize = (10,10))
plt.axis('off')
plt.imshow(img, cmap="gray")