Article 4N2BV Fixed points of the Cliff random number generator

Fixed points of the Cliff random number generator

by
John
from John D. Cook on (#4N2BV)

I ran across the Cliff random number generator yesterday. Given a starting value x0 in the open interval (0, 1), the generator proceeds by

xn+1 = | 100 log(xn) mod 1 |

for n > 0. The article linked to above says that this generator passes a test of randomness based on generating points on a sphere.

Real numbers

While the long term distribution of the generator may be good, it has a problem with its sequential behavior, namely that it has an infinite number of fixed points. If the generator ever reaches one of these points, it gets stuck forever.

Here's a proof. Since x is between 0 and 1, log(x) is always negative. So we can replace the absolute value above with a negative. A fixed point of the generator is a solution to

x = -100 log(x) - k

for some integer k. Define

f(x, k) = -100 log(x) - x - k.

For each non-negative k, the limit of f(x, k) as x goes to 0 is a and the limit as x goes to 1 is negative, so somewhere in between it is zero.

Floating point numbers

So in exact operations over the reals, there is a fixed point for each non-negative integer k. However, when implemented in finite precision arithmetic, it manages to get itself unstuck as the following Python code shows with k arbitrarily chosen to be 20.

 from numpy import log from scipy.optimize import bisect r = bisect(lambda x: -100*log(x) - x - 20, 0.4, 0.999) print(r) for i in range(10): r = abs(-100*log(r)) % 1 print(r)

This produces

 0.8121086949937715 0.8121086949252678 0.8121087033605505 0.8121076646716787 0.8122355649800639 0.7964876237426814 0.7543688054472710 0.1873898667258977 0.4563983607272064 0.4389252746489802 0.3426097596058071

In infinite precision, r above would be a fixed point. In floating point precision, r is not. But it does start out nearly fixed. The first iteration only changes r in the 11th decimal place, and it doesn't move far for the next few iterations.

Update: See the next post for how the random number generator does on the DIEHARDER test suite.

Plotting fixed points

The kth fixed point is the solution to f(x, k) = 0. The following code plots the fixed points as a function of k.

 t = arange(300) y = [bisect( lambda x: -100*log(x)-x-k, 0.001, 0.999) for k in t] plt.plot(t, y) plt.xlabel("k") plt.ylabel("fixed point")

cliff_fps.svg

The fixed points cluster toward zero, or they would in infinite precision arithmetic. As we showed above, the Cliff random number generator performs better in practice than in theory. Maybe the generator does well after wandering close to zero, but I wouldn't be surprised if it has a bias toward the low end of the interval.

Related posts8NWefxR304w
External Content
Source RSS or Atom Feed
Feed Location http://feeds.feedburner.com/TheEndeavour?format=xml
Feed Title John D. Cook
Feed Link https://www.johndcook.com/blog
Reply 0 comments