Hyperchaos

Section 6.7.3 of my increasingly well-thumbed copy of Sprott [1] introduces coupled chaotic systems. For example we have two Lorenz systems:

Paper patching - this needs all four multipliers and
 six of the eight integrators.

In the above, k is the coupling between the two systems (k = 0 means no coupling - and the systems are independent).

And here's the yx-plot of the attractor (initial conditions as per Fig 6.27, [1]):

k = 0.6, R = 1, y vertical, x horizontal

Compares nicely with Sprott's figure 6.27 [1]. Given it's a coupled system, it's interesting to look at the projection of the attractor with (say) x (from 1st system) horizontal and v (from 2nd system) vertical:

xv projection with k = 0.6

As k gets bigger, there is evidently a tendency for the two systems to become synchronised [1] - for example with k = 0.68 we have the rather attractive (aesthetically speaking) attractor:

xv projection with k = 0.68

With k = 0.7 we have the limit cycle attractor:

xv projection with k = 0.7

Reference

[1] Elegant Chaos: Algebraically Simple Chaotic Flows, Julien Clinton Sprott, World Scientific, 2010, p. 157.

Götthans-Petržela Oscillator

Previously (Beyond Rössler, 1st February 2015), I looked at:

dx/dt = - ax - cos(by),

dy/dt = -ay - cos(bz)

dz/dt = -az - cos(bx).

An interesting variant [1] is to take the sign of the trigonometric function - in effect replacing the (co)sinusoidal right hand side terms with square waves:

dx/dt = - ax - sign[cos(by)],

dy/dt = -ay - sign[cos(bz)]

dz/dt = -az - sign[cos(bx)].

So as to compare results with those in reference [1], I went for the sine function:

dx/dt = - ax + sign[sin(by)],

dy/dt = -ay + sign[sin(bz)]

dz/dt = -az + sign[sin(bx)].

Here's my patching for the above:

(My computer only has two x > y functions - these are used for two of the signums. The remaining one is cooked up using the output x < 0 of the spare special function and an electronic switch. In hindsight, perhaps a minimum of three of each functionality would be nearer the mark - given that we live in a three-dimensional world.)

a = 1, b = 10. y axis is vertical, x horizontal, 2 V /division.

The above compares with (the xy-plane) of Figure 4 of reference [1]. Initial conditions are x = 0.1, y = 0, z = 0.

Setting a = 2 gives the following:

(I've fiddled with the colours in the above to agree with [1]: I am guessing that, in Figure 8 of reference [1], the red plot is the xy-projection, the blue is the xz projection.)

Finally, setting a = 0.1 gives an obvious clipping problem - but in the region for which |x| and |y| < 10 V, I get the same sort of structure as in Figure 9 of [1]:

y vertical, x horizontal, 5 V / div.

Interesting stuff. And probably the most complicated problem I've patched up:

Postscript 

A benefit of the analog computer approach is that it's entirely trivial to try out ideas like this:

where n = 1, 2, ...,6 (at least that's what's available on my computer!).

I (naively) thought that the above might put a cat amongst the pigeons, so to speak - i.e. yield a more 'complicated' attractor - because - at least for odd values of n - we have things like:

...but in an annoying face palm moment it's obvious that:

durhhhh - and so I got the same attractors as before.

(I had a look at these attractors sans the sign functions - but they looked pretty much like the ones for n = 1. I guess the coefficients in front of the 'higher frequency' terms in the above expansions weaken the effect / contribution of those bits.)

Reference

[1] Petrzela, J., Gotthans, T., Hrubos, Z. Analog implementation of Gotthans-Petrzela oscillator with virtual equilibria, Radioelektronika (RADIOELEKTRONIKA), 21st International Conference, 19-20 April 2011, pp. 1-4.

Lorenz 84 System

A model for the long-term atmospheric circulation closest to the equator (the Hadley cell), proposed by Lorenz [1] in 1984, is given by:

I really like this sort of pragmatic physics - i.e. simple low order models which still retain the importance of what's actually going on in Nature. What amazes me also is how easy / quick it is to sketch out how to patch these up on the analog computer:

Adopting a = 0.25, b = 4.0, F = 8.0 and G = 1.0 (as used by Lorenz in reference [1]), yields the following attractors. These show nice agreement with those given in [2]:

Lorenz 84 System Attractor - lower plots from [2]. Upper plots 1 V / div except left plot X axis at 0.5 V / div.

The following plots show the variations of x, y and z as per reference [1], Figure 5. It seems that the initial conditions used by Lorenz are different to mine (I start off with x = y = z = 0); unfortunately I do not know what Lorenz used as initial conditions - the left hand side of his Figure 5 suggests values other than zero. Interestingly, my x, y and z build up into what seems to be Lorenz's values after a period of time - represented by the right hand vertical bar below:

Variation of x, y and z. Left hand side is Figure 5 from reference [1].

The braced regions look similar - about 72 days of Figure 5 corresponding to 14.5 ms of the oscilloscope plots.

Finally, I had a go at producing a plot of where the attractor intersects the y = 0 plane. To do this I connected the output of the integrator which dealt with the 2nd (i.e. dy/dt) equation to one of the two x = y functions on the computer. If x = y then the scope's z-axis is given -10 V else it gets +10 V. On the Tektronix 465 the latter decreases the intensity of the beam, the former brightens it. The result was pretty impressive - and compares nicely with Figure 7 of reference [1]. The distinctive (Cantor-set) gaps between 'fingers' is evident.

Intersection with y = 0 plane; lower plot from reference [1]. Upper plot 0.5 V / div.

(The scope's brightness is not truly on / off - hence some 'ghost' lines showing through.) Nice that I pick up (some of) the more solid (brighter) regions around z = 1, x = 0.5 V.

Here's what the computer looked like for this problem (note sides / back still need fitting):

...and the xy-projection of the attractor drawn out at slow speed on the plotter (x-axis at 0.2 V / cm and y-axis at 0.5 V / cm):

References

[1] E. N. Lorenz, Irregularity: a fundamental property of the atmosphere, Tellus A, 36, 98-110, 1984.

[2] Christophe Letellier, (Lorenz) 1984 A 3D model for global atmospheric circulation, Université & INSA de Rouen, 3rd March 2008.

Scroll Attractors via ⌈x⌉

My exploration of things chaotic continues...

Jinhu Lü et al present multi-scroll chaotic attractors based on a set of equations of the form [1]:

where, for n-scroll attractors, the so-called saturated function series f(x) has a staircase-like appearance - albeit with finite slopes between each horizontal step.  Implementation of Lü's f(x) on my computer is not feasible...

...but, I happen to have a ceiling function available: ceiling(x) = ⌈x⌉ gives the smallest integer not less than x. This is not quite what reference [1] uses, but I imagined (correctly) worth a go. In fact, a bit of fiddling on the computer showed me that constants a, b and d can be set to unity, leaving the rather compact:

where c is a positive constant. (Subtracting the 1/2 from x keeps the attractor centered.) The following results are all y vertical, x horizontal (0.5 V/cm).

c = 0.5

c = 0.2

c = 0.1

The c = 0.2 case on paper (with the computer running 1000 times more slowly than when doing the oscilloscope plots) looks like this:

Here's the patching:

Not surprisingly, the floor function ⌊x⌋ also works. In this case adding 1/2 to x, i.e. ⌊x + ½⌋ keeps the resulting attractor centered.

Reference

[1] Jinhu Lu, Guanrong Chen, Xinghuo Yu, H. Leung, ' Generating multi-scroll chaotic attractors via switching control', 5th Asian Control Conference, 20-23 July 2004 pp. 1753 - 1761.

He ain't Heavi (side)

...he's my neural network...

From Sprott's book (page 162), we have the four dimensional (minimal conservative) artificial (!) neural network:

dx/dt = tanh(y)

dy/dt = tanh(z - x)

dz/dt = tanh(w)

dw/dt = tanh(-z) .

Here the hyperbolic tangent function seems to have little to do with trigonometry, let alone relativity.

Apparently, tanh(x) is popular with neural network modelling because of the shape of its graph - it runs smoothly from between -1 and +1. Other functions are available.

After not too much effort I ended up with the following patching arrangement, for the above equations, on my computer (which just happens to have the four hyperbolic tangent functions needed):

And here's a plot of y (vertical) vs. x (horizontal) (both 1 V / division):- if nothing else it does show that the road is indeed long / with many a winding turn...

(Initial conditions were as per Figure 6.31(b) of the book - note that the figure caption has 'v' for the vertical axis - I only got the similar looking plot when I plotted y vertically - rather than v (which I tried interpreting as dx/dt).)

Smorgasbord #2 Nosé–Hoover

Nosé–Hoover oscillator.

Here she blows (with initial condition y₀ = 0.1, all others zero):

Smorgasbord #1 Rikitake Dynamo

The Rikitake system, which apparently serves as a model for the reversal of polarity of the Earth's electromagnetic field, is written:

Interestingly, I believe this predates the Lorenz attractor.

I wired up the computer to model the above; phase portraits of the Rikitake attractor for μ = 2, α = 5, x₀ = 0, y₀ = 0.1 and z₀ = 0 are shown below (1 V/cm, scope in XY-mode), together with results from reference [1]. I've tried to adjust my results to match the scales used in [1]. Not sure why my results show a preponderance for the right hand side - I think it may be because the analog computer solution decays (over time) due to the integrators not being completely error free: by the time the solution is ready to move across to spiral around on the left hand side, errors have possibly built up, destroying the solution.

[1] Ihsan Pehlivan and Yulmaz Uyaroglu, 2007. Rikitake Attractor and its Synchronization Application for Secure Communication Systems. Journal of Applied Sciences, 7: 232-236.

Bessel Attractors

I was looking at multi-scroll attractors - based on the 3-D linear autonomous system mentioned in [1, equation (11)],

More about the multi-scroll attractors in a later post.

In the meantime I wondered - as you do - what would happen if the function on the right hand side of dz/dt in the above were a Bessel function (first kind - and indeed the only kind on my computer!)...as an aside I am fascinated by an age in which each differential equation effectively had its own named function (what we would now just call a solution)...years ago I managed to get incomplete Anger functions into my thesis...

So, in jerk form we have:

A trivial re-arrangement gives

which has the form of Sprott's 'memory oscillator'...[2, page 82]:

although in my case the first derivative has a non-unity coefficient, b. After integrating the above, term by term, with respect to time, I think it would be reasonable to associate b with some measure of the square of the undamped angular frequency of the oscillator. (It would be interesting to vary b...).

Below are my analogue computer results obtained for a = 0.4, b = c = 0.8, q = 8, and n = 0, 1, and m = 1, 2 and 5:

(Here the top row is at 5 V/div; bottom row is 1 V/div; all others 2 V/div).

The J₁(5x) are worthy of being plotted on actual paper:

XY-plane

XZ-plane

YZ-plane

The last one looks almost composed of cylindrical orbits about a straight axis. Bessel functions are cylinder functions after all...

References

[1] Jinhu Lu, Guanrong Chen, Xinghuo Yu, H. Leung, ' Generating multi-scroll chaotic attractors via switching control', 5th Asian Control Conference, 20-23 July 2004 pp. 1753 - 1761.

[2] Julien C. Sprott , Elegant Chaos: Algebraically Simple Chaotic Flows, World Scientific Publishing Company, 2010.

Beyond Rössler

I came across an interesting paper [1] about something called Labyrinth chaos.

The underlying form of the equations (in three dimensions, at least) which yields the chaos looks like:

dx/dt = f(x,y,z)

dy/dt = f(y,z,x)

dz/dt = f(z,x,y),

where the functions are all identical, except the variables are rotated. This, I believe, is called a circulant system [2]. The above-mentioned Labyrinth aspect occurs when f(x,y,z) is (I assume) periodic e.g. sine, cosine.

Götthans and Petržela give results, amongst others, for the case in which

f(x,y,z) = -ax - cos(by),

where a and b are constants.

Now, as it happens I now have an analog computer with trigonometric (amongst other) functions...

Interestingly, my understanding of the paper by Götthans and Petržela is that their approach is similar to mine: they have used a circuit implementation of the equations, with their integrations being dispatched using op amps with capacitor feedback, and their f(x,y,z) parts dealt with digitally.

(My 'special' functions are similarly dealt with digitally (urghh!!) by looking up MATLAB-calculated values  from an EEPROM, sandwiched between an ADC and DAC - see my post of November 22, 2014...)

My first attempt at implementing Götthans and Petržela's equations directly,

dx/dt = - ax - cos(by),

dy/dt = -ay - cos(bz)

dz/dt = -az - cos(bx),

ended up with the integrators overloading (i.e. outputs greater than 10 V or less than -10 V). To overcome this problem I made the substitutions x = 2u, y = 2v, z = 2w, giving

du/dt = - au - cos(2bv)/2,

dv/dt = -av - cos(2bw)/2

dw/dt = -aw - cos(2bu)/2,

which fixed up the overloading issue. (The above factor of 2 was arrived at by guesswork. This whole area of scaling remains a bit of a black art...)

Here's my thoughts on paper as to how to connect the bits of the computer up:

...and here's how it was patched on the actual thing:

Integrator initial conditions were all set to zero volts.

Et voilà, results for three different values of constant b:

(In all cases the scope is 1 V/div in XY mode. I connected the scope via two times-2 amplifiers on the computer (not shown in above patch diagram), so x = 2u is connected to horizontal axis and y = 2v is connected to vertical axis.)

b = 1, a = 0.4

b = 5, a = 0.4

b = 2, a = 0.4

These look pretty much like the results in [1]: their figure 6 (negative cosine function) - although theirs are projections of the rotated 3D space, and there may well be some issues with signs (in my part - my b = 1 case seems to be in a different quadrant, compared to [1]).

I particularly like the b = 5 case with the computer in FAST mode - the result on the scope looks like wisps of smoke; in normal mode (1 µF / 1 MΩ), the spot darting around on the scope looks very life-like - a fly or something...

And, finally...on paper it looks like (b = 5 case):

x is horizontal, y vertical, b = 5.

References

[1] Experimental Study of the Sampled Labyrinth Chaos. Tomáš Götthans, Jiří Petržela, Radioengineering, Vol. 20, No. 4, December 2011.

[2] Elegant Chaos - Algebraically Simple Chaotic Flows. Sprott, Julien. C. World Scientific, 2010, p. 101.