Posted by: Silvio Capobianco | 2014/05/15

## Random Settings in Cellular Automata Machines

In cellular automata machines such as the CAM8 and SIMP/STEP (which is CAM8’s software incarnation) memory is organized into planes: each plane represents a certain type of information. A single cell spans through planes.

Consider, for example, the HPP (Hardy, de Pazzis and Pomeau) lattice gas model. In HPP, gas particles travel along the lines of a square grid and collide in the nodes of the grid. The collision rule states that, if exactly two particle meet in a node, and they come from opposite direction, then they bounce with a $90^\circ$ angle; in all other cases, particles go straight ahead. In SIMP/STEP, we can construct the HPP lattice gas model by defining interaction as follows:

    def hpp_rule():
if north == south and east == west:
north._ = east
east._  = south
south._ = west
west._  = north
hpp_step = Rule(hpp_rule)



Suppose, however, that the update rule is not deterministic but randomized. For example, in the FHP lattice gas model, which has six direction on a triangular lattice, the principle is the same as for HPP, but the direction (clockwise or counter-clockwise) of the $60^\circ$ rotation is chosen at random. This can be done by initializing a single plane with $50\%$ zeros and ones, and have the cell read the corresponding bit and interpret it as the direction of the rotation. To ensure non-repetition, during the update the plane will be moved by a random offset. In SIMP/STEP, it would look more or less like the following:

    def fhp_rule():
if north == south and east == west and ne == sw:
if r: # rotate clockwise
north._ = west
west._  = sw
sw._    = south
south._ = east
east._  = ne
ne._    = north
else: # rotate counter-clockwise
north._ = ne
west._  = north
sw._    = west
south._ = sw
east._  = south
ne._    = east
fhp_step = Sequence([Stir([r]), Rule(fhp_rule)])



It might be the case, however, that we need different probability distributions in different points. In such a case, it would be anti-economical to set as many random planes as needed distributions: instead, one would want to set as few planes as possible with random distributions, and exploit the flexibility given by the application of boolean functions to said planes. Different cells might then have different types, and interpret differently the state of the random planes.

Let us make a quick example with a single plane, initialized with probability $p_0$. If we pass a bit, whose value is $1$ with probability $p_0$, to a boolean function $f$ of a single argument, then the output bit will have a Bernoulli distribution with parameter $p$ which depends on both $p_0$ and $f$. Specifically, $p$ will be 0 if $f$ is the constant function $0$, $p_0$ if $f$ is the identity, $1-p_0$ if $f$ is the negation, and 1 if $f$ is the constant function 1: which exhausts the possibilities of a single bit.

In general, $N$ random bits with independent Bernoulli distrinutions of parameters $p_0, p_1, \ldots, p_{N-1}$ can yield up to $2^{2^N}$ possible distributions for the output bit of an $N$-to-1 boolean function. It then becomes of interest to understand how to set the parameters $p_0, p_1, \ldots, p_{N-1}$ in order to obtain a convenient range of values for $p$ just by changing the specific $f$. In particular, we are interested in obtaining the uniform range

$0, \dfrac{1}{2^{2^N}-1}, \ldots, \dfrac{2^{2^N}-2}{2^{2^N}-1}, 1$

which minimizes the maximum distance between the value desired and the one which is actually obtained. A general solution to this problem was described by Mark Smith in his 1994 PhD thesis under the supervision of Tommaso Toffoli.

Before examining Smith’s method, let us observe that the denominator in the formula above is the $2^N$-th Mersenne number. The sequence of Mersenne numbers is defined by the formula

$M_n = 2^n - 1 \;\; \forall n \geq 1$

In general, Mersenne number are not primes: indeed, it is a standard exercise in elementary number theory that $M_m$ is a factor of $M_n$ if $m$ is a factor of $n$, so that a Mersenne number can be prime only if its index is prime: such condition is not sufficient, as $M_{11} = 2047 = 23 \cdot 89$.

Let us go back to the case $N=1$: we know that, in this case, given the value of the parameter $p_0$ we can get any of the values $0$, $p_0$, $1-p_0$, $1$. To have these values coincide with $0$, $1/3$, $2/3$, $1$ it is thus sufficient to set $p_0 = 1/3$.

Let us now consider the case $N=2$. We know that, in this case, the denominator should be $2^{2^2}-1 = 15$: so we are strongly tempted to set $p_0 = 1/3$ again, then $p_1 = 1/5$, which would give us $p = 1/15$ when $f(b_0, b_1) = b_0 \wedge b_1$ is the boolean AND. But how do we know that we may get the remeining values as well?

To prove this, let’s do a thought experiment. Suppose we have a rectangle, whose height and width are in $3:5$ ratio: let us make two cuts, one on the short side at $1/3$ of the heigth, the other on the long side at $1/5$ of the width. The sizes of the four pieces obtained are in proportion with each other as the numbers $1 = (1)_2$, $2 = (10)_2$, $4 = (100)_2$, and $8 = (1000)_2$. But every integer number between $0$ and $15$ has a unique writing as a sum of one or more of the fourn numbers above: and indeed, we can also choose the output probability $p$ with any given value of the form $k/15$ with $k$ an integer between $0$ and $15$. To decide which function $f$ we have to use to get probability $p = n/15$ for the output b, we proceed as follows:

1. We write $(n)_{10} = (b_3 b_2 b_1 b_0)_2$.
2. We choose $f$ so that $f(i) = 1-b_i$ for every $i = 0, 1, \ldots, n-1$.

For example, if we want $p = 7/15$, we take the function $f$ that takes value $1$ only on the pair $(1,1)$.

To generalize to arbitrary many bits, we introduce Fermat numbers, defined by the formula

$F_n = 2^{2^n} + 1 \;\; \forall n \geq 0$

We easily get

$F_n \cdot M_{2^n} = M_{2^{n+1}}$

that is, the product of the $n$-th Fermat number and the $2^n$-th Mersenne number is the $2^{n+1}$-th Mersenne number. This leads us to

Theorem. If the probability of bit $b_i$ is $1/F_i$, then all and only the probabilites $k/M_{2^{n+1}}$ with $k$ between $0$ and $M_{2^{n+1}}$ can be obtained. Such arrangement is the only one that allows recovering all said probabilities in output, up to a rearrangement of the bits.

That such arrangement of probabilities for the input bits allows to obtain every possible value in the range $0, 1/M_n, \ldots, (M_n-1)/M_n, 1$, is easily proven via a generalization of the reasoning we made for the case $n=2$. Uniqueness follows inductively from the relation $M_{2^{n+1}} / M_{2^n} = F_n$.