Random number generation is often done by applying a function recursively to a seed number. The function needs to be bounded to a certain domain. Obtaining a function that is bounded is easily done by applying the modulo operation on the base function. The function must be deterministic. For example the following function
\[ f(x) = (x + 1) % 3 \]
will supply us with the sequence \(0,1,2\) looping indefinitely if we take seed value \(0\). If we want to detect the cycle in a sequence we can use the equation \(X_n = X_{2n}\) for some \(n\). Proofing this relation is left as an exercise. Now we can create some code to detect cycles for a function:
def period_fun(f, x0): """return length and start of recursively applied function 'f' with start value 'x'""" x1 = f(x0) x2 = f(f(x0)) n = 1 while x1 != x2: x1 = f(x1) x2 = f(f(x2)) n += 1 xl = f(x1) l = 1 while xl != x2: xl = f(xl) l += 1 xm = x0 m = 0 while xm != x1: xm = f(xm) x1 = f(x1) m += 1 return l, m, xm
In this code the l
indicates the length of a period, and m
is the
index of the first digit of the cycle. To test if this works, let’s
define a recursive function. A famous random number generating function
by Von Neumann is the Middle-square method. It is here implemented for
2-digit numbers.
def middle_square2(d): return int(str(d * d).zfill(4)[1:3])
Here is a image from the wikipedia, showing the possible cycles:
Let’s see if our results compare with a simple test.
from numpy import zeros n = 100 result = zeros((n, 4)) for i in range(100): l, m, e = period_fun(middle_square2, i) result[i] = [i, l, m, e]
Now visualize our results with matplotlib, the default Python workflow.
import maptlotlib.pyplot as plt fig, ax = plt.subplots(1, 3) ax[0].hist(result[:,1]) ax[0].set_title("cycle length") ax[1].hist(result[:,2]) ax[1].set_title("init sequence length") ax[2].hist(result[:,3]) ax[2].set_title("end number") plt.show()
It looks like it works, you can also rewrite the period_fun
to take
sequences instead of a function. Try to implement it yourself. Also test
some other functions to see what happens.