Random Sequence Generator based on Chua Circuit

Introduction

First of all, one may ask: "What is a Chua Circuit?"

The Chua Circuit is the simplest electronic circuit exhibiting chaos, and many well-known bifurcation phenomena, as verified from numerous laboratory experiments, computer simulations, and rigorous mathematical analysis. It was introduced in 1983 by Leon O. Chua, who was a visitor at Waseda University in Japan at that time. The ease of construction of the circuit has made it a ubiquitous real-world example of a chaotic system, leading some to declare it "a paradigm for chaos." (from [1] and [2]).

The Chua Circuit can be accurately modeled by means of a set of three nonlinear ordinary differential equations in the variables x(t), y(t) and z(t):

\left\{ \begin{array}{l} \dot{x} = \alpha (y - \varphi(x)) \\ \dot{y} = x - y + z \\ \do{z} = -\beta y\end{matrix} \right

where α and β are two real numbers, and φ(x) is a nonlinear function. By varying α and β one can observe the period-doubling bifurcation route to chaos [2], and create a strange attractor known as double scroll: a beautiful, truly chaotic pattern.

Using a Chua Circuit, in theory, it is possible to build a random sequence generator, as explained in [3]. To do this, one exploits the fact that a chaotic signal has broadband, noise-like, continuous power spectrum, and is very sensitive to initial conditions. A small disturbance to the chaotic circuit will lead to large variation of the output signals, which is named butterfly effect. Therefore, it is actually impossible to foretell the output signal after a time interval long enough, even if we know with precision the initial conditions of the chaotic circuit. That is, if the chaotic signals were sampled at a rate slow enough and the output was set to 1 or 0 according to their amplitudes or phases, a random sequence can be obtained.

The Chua Circuit

The Chua Circuit is indeed easy to build, and someone has even provided a kit [4]. Most implementations make use of an inductor though, like [5], meaning that you have either to buy one or wind it up yourself. There is, however, an implementation based on State Controlled CNN that does not require an inductor. If you're interested on how this is done, or you want to know how to implement a specific system of differential equations with a circuit, I suggest you to read [6] (chapters 2 and 3). I will skip over the details of the design process and give you the result: the schematic of the Chua Circuit implemented using State Controlled CNN (click on the images to enlarge):

chua-schematic

Component values:

R1 20k        R6 1k pot.    R11 1k        R16 1k        R21 8k        TL082 (x2)
R2 13.7k      R7 100k       R12 8.2k      R17 1Meg      R22 1k
R3 5.6k       R8 100k       R13 100k      R18 112k      C1  100n
R4 4k         R9 100k       R14 8.2k      R19 112k      C2  100n
R5 20k        R10 100k      R15 100k      R20 1Meg      C3  100n
			

The circuit isn't particularly complicated, but it uses 4 op-amps, 3 capacitors and quite some resistors. As you might have noticed, not all resistor values belong to the E24 series, so you'll have to use 1% resistors (as I did) or use series/parallel combinations to get close to the required value.
With the resistor values chosen, I have set β to a fixed value, β=R15/R14=12.2, but we are able to change α from R5/R3=3.57 to, theoretically, ∞ by reducing R6 from 1kΩ towards zero. This will permit us to observe the change in behaviour while we vary this resistance. Note that R6 needs to be a multiturn trimmer, because a small variation of its resistance will have quite a big effect on the circuit dynamics. Lastly, you need a dual power supply capable of providing ±9V for this circuit. The current draw shouldn't exceed 20/30mA, so even using batteries could be feasible.

Spice simulations

Before reaching to the breadboard and assemble the circuit, it's always a good practice to run a couple of simulations to validate our design. If something doesn't work as expected it's easy to make a correction at this stage and iterate the process until we're sure everything's in order. As you can see on the schematic, I have labeled the nodes in the circuit with progressive numbers, and these same numbers will be used in the netlist that you can download here. The simulations have been performed using ngspice. As always, to launch the simulation run:

$ ngspice chua.net

I chose to simulate the transient response of the circuit for 100ms. The parameter α was set at 9.4 by setting R6=380Ω, and this should give us the famous double scroll. The first thing to do is looking at the state variables, plotted for the last 20ms of simulation:

X variable:

x

Y variable:

y

Z variable:

z

We can observe that these chaotic signals seem aperiodic, as they should be. This can be confirmed by plotting the FFT of X for example:

fft_x

Now let's look at state variables in pairs. What we can notice is that X and Z seem to have about the same amplitude and opposite phase for this value of α, while the amplitude of Y is much smaller:

X and Y:

xy

X and Z:

xz

Y and Z:

yz

And on to the interesting part: the double scroll. To observe this patter we need to plot one state variable versus another. Here are the pictures of all three combinations:

Y vs X (double scroll):

double-scroll-yx

Z vs X (double scroll):

double-scroll-zx

Y vs Z (double scroll):

double-scroll-yz

This is what we were looking for: the behaviour of this circuit is indeed chaotic and looking at the three attractors their trajectories never cross in phase plane, as expected.

Circuit test

I built the circuit on a breadboard, powered it with a ±9V dual supply and connected the oscilloscope probes to the state variables. But all I could see was...well, nothing. The state variables were constant. The circuit was in one of its equilibrium point. So I started to turn R6 slowly, in the sense that decreased its resistance, and sure enough variable signals started to appear on the oscilloscope. All in all, the circuit is quite robust, and you can generate a whole set of patters, including single and double scrolls of course. So I set R6 at about 400Ω and looked at the signals: what I found was in perfect accordance with the simulations. These are the pictures I took (sorry if they are quite dark, but I had to a lower the shutter speed in order to freeze the waveforms):

X state variable (V=1V/div, H=2ms/div):

osc_x

Y state variable (V=0.2V/div, H=2ms/div):

osc_y

Z state variable (V=1V/div, H=2ms/div):

osc_z

Y vs X (V=0.2V/div, H=0.5V/div):

osc_ds_yx

Z vs X (V=1V/div, H=0.5V/div):

osc_ds_zx

Y vs Z (V=0.1V/div, H=0.5V/div):

osc_ds_yz

It was very interesting to witness the famous period-doubling bifurcation route to chaos. Success, we can affirm that our Chua Circuit works!

Random sequence generator

So, we have managed to create chaos inside an electronic circuit, next questions is: can we use it to generate a true random bit sequence? According to [3], we can. In fact, randomness coming from initial conditions is one of the three mechanisms responsible for random behavior in systems. This aspect is studied by chaos theory and is observed in systems whose behavior is very sensitive to small variations in initial conditions, such as the Chua Circuit. So, the definition of chaos itself gives us the key to understand why the Chua Circuit can be used for this purpose:

Chaos theory studies the behavior of dynamical systems that are highly sensitive to initial conditions, an effect which is popularly referred to as the butterfly effect. Small differences in initial conditions (such as those due to rounding errors in numerical computation) yield widely diverging outcomes for chaotic systems, rendering long-term prediction impossible in general. This happens even though these systems are deterministic, meaning that their future behavior is fully determined by their initial conditions, with no random elements involved. In other words, the deterministic nature of these systems does not make them predictable. This behavior is known as deterministic chaos, or simply chaos. (from [7])

If we look again at the Z state variable waveform, we notice that it has variant signs in different scrolls. Thanks to the butterfly effect of chaos, we can not foretell in which scroll will the trajectory be located after a time interval long enough. So, if we sample the Z signal with a low rate and set output as 1 or 0 according to its sign, we will obtain a random sequence. This is easily done using another op-amp working as a comparator (open loop), as in the following diagram:

chua-comparator

The comparator output looks like this (the values of the voltage divider Rc1/Rc2 have been chosen so that we never exceed 3.3V in the worst case when using a ±9V dual supply):

comp

The only things that remains to be done is to acquire these bits coming out of the comparator for example every 200ms (assuming this is a time interval long enough). For this I used a Raspberry Pi and a few lines of C code, but it can be easily done in bash or Python as well. Just configure GPIO4 as input at the beginning, then read the bit, print its value, wait the delay and repeat. A voltage level of 2.2V at the output of the comparator should be enough to be recognized as high by the RPi.
In particular, I used two similar programs: the first one sits in a endless loop where it reads the bit and prints it to the standard output. If you use this solution and don't want to wear out your SD card with all the writes, you can use netcat and send the output of your program to another computer on the network. Run this on your computer:

$ netcat -v -l PORT | tee chua-sequence

and then start sending the bits from the RPi:

$ sudo ./chua | netcat -v HOST PORT

The second one reads bits for an hour into a memory buffer and then appends this buffer to the output file, thus minimizing writings to the SD card. The process is repeated until you send SIGINT, in which case the partially full buffer is appended to the output file before the program terminates. You can start this one via ssh and leave it running for long acquisitions using the following command:

$ nohup sudo ./chua_hour &

Both programs save the bitstream in ascii format, where each line is either "0\n" or "1\n". I chose this format because it's easy to import in Matlab/Octave for further processing.
Lastly, I added the LED and the transistor because I wanted to see the Raspberry Pi "absorbing" the chaos. The software should configure GPIO9 as output, and then setting or clearing this output based on the bit just read.

Acquiring the sequences

While testing this bit generator, I noticed that it presented a systematic bias (that is the bitstream was prone to be biased, with either 1s or 0s predominating). So before acquiring a bit sequence we need to minimize bias inherent in the operation of the generator. This can easily be done with the following calibration technique: display the ZX double scroll on your oscilloscope and make sure that the two scrolls are equally bright and about the same dimension (both the external dimension and the "hole" at the center of the scroll). To do this you have to slowly adjust either the +9V supply or the -9V, with steps of, say, 10mV. The reason why the shape of the scroll depends on the supply voltage is because in this implementation the nonlinearity depends on the Vsat of the operational amplifiers, which of course change with the supply voltage. Equal brightness means that our trace will spend approximately 50% of the time in each scroll, giving us an equal probability of ones and zeros.
Note that even if the above hardware bias reduction steps have been taken, the bit-stream should still be assumed to contain bias and correlation, so a de-skewing technique for reducing bias and correlation should be used (more on this later).
My initial goal was trying to reproduce the results presented in the paper listed as [3] on my page. In order to do this I chose a time inverval of 200ms and ran every acquisition for 14 hours. This means I was able to collect 252000 samples at 5 Samples/s. So far I have acquired the following sequences:

Basic tests

Yes, I hear you asking: "Are these sequences really random?" Well, I am by no means an expert in this field, but I have been trying to analyze the data collected in several ways to make sure I didn't make any obvious mistake.
The simplest and easiest way to visually check if a sequence looks random is to display it as a bitmap image with each pixel representing one bit. This is not, of course, a serious test, but it's enough to show if something is obviously wrong. The image can be generated with the command:

$ cat seq1_252000_200ms | ./make_bitmap 400 400

The following table shows the images for all the sequences:

Sequence 1

seq1

Sequence 2

seq2

Sequence 3

seq3

Sequence 4

seq4

These bitmaps do not show any patter, and they look indistinguishable from white noise to the human eye, as expected.

Another simple way to test for randomness is evaluating the Monte Carlo for π as explained in [8]. This method is based on this idea: we take a square and inscribe within it a circle that touches each edge of the square. We know that if the radius of the circle is r, then the area of the circle is πr², and the area of the square is 4r². If we calculate the ratio, q, of the circle area to the square area, we get q = π/4; we can therefore calculate π with π = 4q. The same result still holds if we work in the first quadrant.
The ratio q can be found using pairs of random points (x, y) that we extract from our sequence. In particular (x,y) pairs are generated using blocks of consecutive 48-bit, with each coordinate being a 24-bit number. If we count the number of points that fall inside the circle and divide this number by the total number of points we get an estimation of q. For independent random points that are uniformly distributed within the square, this method should give us a sequence that slowly converges to π.
These are the results for all the sequences:

$ cat seq1_252000_200ms | ./montecarlo_pi 10000 2>points | tail -n1
n= 5250  pi=3.163428571  e= 0.695%
$ cat seq2_252000_200ms | ./montecarlo_pi 10000 2>points | tail -n1
n= 5250  pi=3.119238095  e=-0.712%
$ cat seq3_252000_200ms | ./montecarlo_pi 10000 2>points | tail -n1
n= 5250  pi=3.179428571  e= 1.204%
$ cat seq4_252000_200ms | ./montecarlo_pi 10000 2>points | tail -n1
n= 5250  pi=3.116190476  e=-0.809%
			

Full log for this test here.

Sequence 1

seq1

Sequence 2

seq2

Sequence 3

seq3

Sequence 4

seq4

Lastly, I wrote a simple script in Octave to try to reproduce one of the tables in [3] and these are the results:

Sequence number Length Number of 0 Number of 1 NRL-1 NRL-2 NRL-3 NRL-4
1 252000 125378 126622 63054 31587 15664 7851
2 252000 125801 126199 62784 31728 15874 7858
3 252000 127018 124982 62844 31730 15658 7866
4 252000 125858 126142 62653 31317 15850 7907

NRL-n means number of runs with length=n. For a random sequence we expect that: 1) the number of ones is nearly equal to that of zeros, and 2) the number of runs with length=n is about 1/(n+1) of the number of all runs. So, using these properties the numbers look good.

Statistical tests

I did some reading about this problem, and it turns out there are lots of tests to assess the statistical properties of RNGs. A finite set of statistical tests may only detect defects in statistical properties. Furthermore, these tests cannot verify the unpredictability of the random sequence. For us, unpredictability of the sequence is a consequence of having a chaotic signal as a source.

A natural source of random bits may not give unbiased bits as direct output. Many applications, especially in cryptography, rely on sequences of unbiased bits. There are various techniques to extract unbiased bits from a defective generator with unknown bias. These are called de-skewing techniques or whitening algorithms. These techniques also eliminate the correlation in the output of the natural sources of random bits. John von Neumann invented a simple algorithm to fix simple bias, and reduce correlation. It considers successive pairs of consecutive, non-overlapping bits from the input stream, taking one of three actions: when two successive bits are equal, they are discarded; a sequence of 1,0 becomes a 1; and a sequence of 0,1 becomes a 0. It thus represents a falling edge with a 1, and a rising edge with a 0. This eliminates simple bias, and is easy to implement. This technique works no matter how the bits have been generated. It cannot assure randomness in its output, however. What it can do is transform a biased random bit stream into an unbiased one. However, significant numbers of bits are discarded. Thanks to the hardware bias reduction, we can affirm that the probability of 1 is p ~ 1/2, and in this case the discard of input pairs is minimum, giving us an output stream that is 1/4 the length of the input on average.

I have considered the following tools to analyze the sequences:

rngtest: rngtest works on blocks of 20000 bits at a time, using the FIPS 140-2 (errata of 2001-10-10) tests to verify the randomness of the block of data. It will use the first 32 bits of data to bootstrap the FIPS tests (these bits are not tested for randomness). Each block is subjected to five tests: monobit test, poker test, runs test, long run test, continuos run test. If any of these fails the block fails the test. For more informations, visit [9].
I have ran this test on all the sequences before applying de-skewing. This is because some whitener designs can pass statistical tests with no random input. While detecting a large deviation from perfection would be a sign that a true random noise source has become degraded, small deviations are normal and can be an indication of proper operation.
This is what the output of rngtest looks like on the first sequence:

$ cat seq1_252000_200ms | ./ascii_to_bin | rngtest 
rngtest 2-unofficial-mt.12
Copyright (c) 2004 by Henrique de Moraes Holschuh
This is free software; see the source for copying conditions.  There is NO warranty; not even for MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.

rngtest: starting FIPS tests...
rngtest: entropy source exhausted!
rngtest: bits received from input: 252000
rngtest: FIPS 140-2 successes: 12
rngtest: FIPS 140-2 failures: 0
rngtest: FIPS 140-2(2001-10-10) Monobit: 0
rngtest: FIPS 140-2(2001-10-10) Poker: 0
rngtest: FIPS 140-2(2001-10-10) Runs: 0
rngtest: FIPS 140-2(2001-10-10) Long run: 0
rngtest: FIPS 140-2(2001-10-10) Continuous run: 0
rngtest: input channel speed: (min=1.041; avg=4.935; max=6357.829)Mibits/s
rngtest: FIPS tests speed: (min=45.413; avg=48.750; max=49.542)Mibits/s
rngtest: Program run time: 58739 microseconds
			

As stated previously, this test works on block of 20000 bits. To make sure our generator passes the test in all circumstances, I repeated the test after moving the block boundaries from 50 to 19950 bits in steps of 50 bits. The results are summarized in the following table, where PASS means the sequence passed the test for all the possible starting offsets:

Sequence 1 Sequence 2 Sequence 3 Sequence 4
FIPS 140-2 PASS PASS PASS PASS

Full log for this test here.

ent: ent applies various tests to sequences of bytes stored in files and reports the results of those tests. The program is useful for evaluating pseudorandom number generators for encryption and statistical sampling applications, compression algorithms, and other applications where the information density of a file is of interest. See [10] for an explation of its output.
Before running this test, each sequence has been passed through the de-skewing step. ent output for the first sequence:

$ cat seq1_252000_200ms | ./deskewing | ./ascii_to_bin | ent -b
Entropy = 1.000000 bits per bit.

Optimum compression would reduce the size
of this 62992 bit file by 0 percent.

Chi square distribution for 62992 samples is 0.02, and randomly
would exceed this value 75.00 percent of the times.

Arithmetic mean value of data bits is 0.4997 (0.5 = random).
Monte Carlo value for Pi is 3.149390244 (error 0.25 percent).
Serial correlation coefficient is -0.002159 (totally uncorrelated = 0.0).
			

Summary of the results:

Sequence number Length Entropy Chi Square Mean Monte Carlo for π Serial correlation
1 62992 1.000000 0.02 (75%) 0.4997 3.149390244 (e=0.25%) -0.002159
2 62944 0.999996 0.32 (50%) 0.4989 3.148741419 (e=0.23%) -0.002483
3 63072 0.999994 0.53 (50%) 0.4986 3.126331811 (e=0.49%) -0.000769
4 62856 0.999985 1.27 (50%) 0.4978 3.126050420 (e=0.49%) -0.003584

Let's repeat this test for all the data we've acquired:

$ cat seq* | ./deskewing | ./ascii_to_bin | ent -b
Entropy = 0.999995 bits per bit.

Optimum compression would reduce the size
of this 251880 bit file by 0 percent.

Chi square distribution for 251880 samples is 1.64, and randomly
would exceed this value 25.00 percent of the times.

Arithmetic mean value of data bits is 0.4987 (0.5 = random).
Monte Carlo value for Pi is 3.164474938 (error 0.73 percent).
Serial correlation coefficient is -0.002230 (totally uncorrelated = 0.0).
			

After 56 hours of operation and 250000 bits out of the generator (after de-skewing), the statistical properties are aligned with what we'd expect from a truly random sequence. Thus, at least, from the data collected so far we cannot say it is not random.
Full log for this test here.

Download

All the scripts and utilities I have used while testing this generator can be downloaded here.

Conclusions

Our random sequence generator seems to pass all the statistical tests I have tried so far on this small number of sequences, and this is a good sign. Of course it would be interesting to analyze many more, longer sequences, but all in all I am satisfied with the results I got with this simple circuit, considering this is something I started out of curiousity and I've learned quite a lot about random number generators while at it. The only downside seems to be its inherent slowless: at an average of 1.25 bit/s after de-skewing it is not exactly fast, but it's okay for my purposes. So, if you're interested in this kind of stuff, I'd say go ahead and build one of these circuits yourself. It is certainly interesting to experiment with, and maybe let me know your results!

Note: I have tried to be as rigorous as possible while building the circuit, acquiring the sequences and performing the tests, but it is quite possible that I have made some mistake in the process. If so, I'd be glad to be pointed in the right direction. My email is on the homepage of this site.

References

[1] Chua's circuit - Wikipedia, the free encyclopedia
[2] Chua circuit - Scholarpedia
[3] A random sequence generator based on chaotic circuits, J. Chen, L. Ran, K. Chen - link1 link2
[4] Chua's Circuit Kit
[5] Robust Op-Amp Realization of Chua's Circuit, M. P. Kennedy - pdf
[6] Chua’s Circuit Implementations: Yesterday, Today and Tomorrow, L. Fortuna, M. Frasca, M.G. Xibilia - google, amazon
[7] Chaos theory - Wikipedia, the free encyclopedia
[8] Efficient Random Number Generation and Application Using CUDA
[9] FIPS PUB 140-2 Security Requirements for Cryptographic Modules
[10] ENT Pseudorandom Number Sequence Test Program
[11] chuacircuits.com

Valid XHTML 1.0 Strict

Copyright (C) 2012 Giorgio Vazzana