automata_rule65535.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 rule65535(x, modulo):
  X = np.empty(x.shape, dtype=np.int16)
  I = x.shape[0]
  J = x.shape[1]
  for i in xrange(I):
    ip1 = (i + 1) % I
    ip2 = (i + 2) % I
    ip3 = (i + 3) % I
    for j in xrange(J):
      jp1 = (j + 1) % J
      jp2 = (j + 2) % J
      jp3 = (j + 3) % J
      t  = x[i, j]   + x[ip1, j]   + x[ip2, j]   + x[ip3, j]  
      t += x[i, jp1] + x[ip1, jp1] + x[ip2, jp1] + x[ip3, jp1] 
      t += x[i, jp2] + x[ip1, jp2] + x[ip2, jp2] + x[ip3, jp2] 
      t += x[i, jp3] + x[ip1, jp3] + x[ip2, jp3] + x[ip3, jp3]
      X[i, j] = t % modulo
  return X

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__':
    # Image size
    power = 5
    modulo = 4

    size = modulo**power
    print "size", size

    permuted_diagonal = np.zeros((size, size), np.int16)
    for i in xrange(size):
        ri = digital_reverse(i, power, modulo)
        permuted_diagonal[ri, i] = 1

    state = permuted_diagonal
 
    for step in xrange(size*power):
        print "step =", step
        save_matrix(state, "img-p%d-g%05d.png"%(power, step), modulo) 
        state = rule65535(state, modulo)