mopdl2_generator_modulo2.py

import numpy as np
from PIL import Image

def digital_reverse(n, length, base):
  r = 0
  for _ in range(length):
    r = base*r + n % base
    n /= base
  return r

assert digital_reverse(123, 3, 10) == 321
assert digital_reverse(54321, 5, 10) == 12345
assert digital_reverse(0x123, 3, 16) == 0x321
assert digital_reverse(0x54321, 5, 16) == 0x12345

def save_matrix(m, fn, modulo):
    image = Image.new("L", (m.shape[1], m.shape[0]), 0)
    pix = image.load()
    for y in xrange(m.shape[0]):
        for x in xrange(m.shape[1]):
             pix[x, y] = 255*m[y, x]/(modulo-1)

    image.save(fn, "png")

if __name__ == '__main__':
    # state size
    modulo = 2
    power = 4

    print "modulo&base =", modulo

    half = modulo**power
    size = 2*half
    print "state size", size

    # known state after applying rule 15
    reversible_state = np.zeros((size, size), np.uint16)
    for i in xrange(half):
        r = digital_reverse(i, power, modulo)
        reversible_state[i, r] = 1
        reversible_state[i, r+half] = 1
        reversible_state[i+half, r] = 1
        reversible_state[i+half, r+half] = 1

    # partially known state
    reversed_state = np.empty((size, size), np.uint16)

    # fill initial (known) values
    reversed_state[0,::2] = 0 # first row
    reversed_state[0,1::2] = 1 # first row
    reversed_state[::2,0] = 0 # first column
    reversed_state[1::2,0] = 1 # first column

    # find rest values of the state
    for i0 in xrange(size-1):
        i1 = i0 + 1
        for j0 in xrange(size-1):
            j1 = j0 + 1
            # solve c4[t-1] from (c1[t-1] + c2[t-1] + c4[t-1] + c8[t-1]) == c1[t] (mod modulo)
            reversed_state[i1, j1] = (reversible_state[i0, j0] + reversed_state[i0, j0] + reversed_state[i0, j1] + reversed_state[i1, j0]) % modulo

    # produce image
    save_matrix(reversed_state, "mopdl2-modulo%d-power%d.png"%(modulo, power), modulo)