Lorenz attractor amusement with Three.js…

First I should admit the two main reasons why I made this project:

  • I wanted to get some experience using Three.js; that is experience in making a 3D web elements. In the past, I have experimented a bit with WebGL, however I wanted to try out a library that could make 3D development a bit easier (that is not needing to use write code for the lower level elements).
  • I wanted to get a better appreciation of the Lorenz attractor. Lately with various books I have been reading, I just accept whatever concepts/ideas and move. Sometimes to only get a complete understanding and appreciation of what’s happening, you need to do the math, make one yourself, or write a program. In doing so, I wasn’t satisfied with plain 2D static pictures; instead I wanted a 3D animation/simulation. ote

A quick comment… Lorenz is not  Lorentz. They are different people and did different things.

Overview…

The Lorenz equations were originally created to model atmospheric convection. They are based on a model introduced that Barry Saltzman introduced in 1962. The original paper is fairly complicated, Lorenz was able to reduce the essence of the model to just 3 equations:

\begin{aligned} \frac{dx}{dt} &= \sigma(y\ -\ x) \\ \frac{dy}{dt} &= x(\rho\ -\ z)\ -\ y \\ \frac{dz}{dt} &= xy\ -\ \beta z \\ \end{aligned}

Where x is proportional to the intensity of convective motion, y is proportional to the temperature difference between ascending and descending currents, and z is proportional to the distortion of the vertical temperature profile. The constants are: \sigma is the Prandtl number which is the ratio of the kinematic viscosity to the thermal conductivity, \rho is the ratio of the Rayleigh number to a critical value, and \beta is a constant related to the given space. Note that in most books, r and b are used for \rho and \beta.

There is a fair bit of math that can be used in the analysis of the above equations, that is enough to fill a chapter or two of a book. To drastically shorten everything up only some of the final results has been listed below. 

The critical points when \rho \gt 1 are defined to be:

c_1 = \left[ \begin{array}{c} 0\\ 0\\ 0 \end{array} \right] \quad c_2 = \left[ \begin{array}{c} \sqrt{\beta(\rho – 1)}\\ \sqrt{\beta(\rho – 1)}\\ \rho – 1 \end{array} \right] \quad c_3 = \left[ \begin{array}{c} -\sqrt{\beta(\rho – 1)}\\ -\sqrt{\beta(\rho – 1)}\\ \rho – 1 \end{array} \right]
  • Case 1: 0 \lt \rho \lt 1
    0 is the asymptotically stable critical point… that is, the system converges to zero.
  • Case 2: \rho = 1
    0 is the neutrally stable solution… again, the system converges to zero.
  • Case 3: 1 \lt \rho \lt 13.926
    Depending on the initial values, the solution will quickly converge to one of the critical points.
  • Case 4: 13.926 \lt \rho \lt 24.74
    Depending on the size of \rho, the system will ‘dance’ around the critical points for a while, and then eventually will converge to one of the critical points
  • Case 5: \rho \gt 24.74
    The famous Lorenz attractor!!!
  • Case 6: \rho \approx 99.65
    There is a stable 6-period orbit in the form a2ba2b.
  • Case 7: \rho \approx 100.75
    There is a stable 3-period orbit a2b, that is it spirals around one critical point once and the other twice.

Note that the book ‘Encounters with Chaos’ was used heavily for the math.

Another thing to note that the Lorenz attractor does look like a butterfly to some people (depending on the viewing angle). However, the main reason why it is associated to the ‘butterfly effect‘ is due to a paper/presentation in 1979 to the American Association for the Advancement of Sciences called ‘Predictability: Does the Flap of a Butterfly’s Wings in Brazil Set Off a Tornado in Texas’. 

Python Implementation…

An implementation was first made in python since it is perhaps the easiest language to implement anything. Two things to note are: 1) there is a package that has various ODE integrators; 2) matplotlib can create 3D plots with some minor encouragement. 

In scipy, there is an ODE integrator similar to MatLabs ODE45. First, the package needs to be imported.

from scipy.integrate import odeint

After a function is defined containing the ODE, then ‘odeint’ function can be used to get some results. Note that initial values contains the starting values for x, y, and z, while timePoints is an array of values representing the times where a solution is to be found. 

def Lorenz (currentValues, t):

    x = currentValues[0]
    y = currentValues[1]
    z = currentValues[2]
    
    dxDt = sigma * (y - x)
    dyDt = x * (rho - z) - y
    dzDt = x * y - beta * z
    
    return [dxDt, dyDt, dzDt]

results = odeint (Lorenz, initialValues, timePoints)

For 3D plotting, an extra package needs to be imported with matplotlib, that is:

from mpl_toolkits import mplot3d
import matplotlib.pyplot as plt

Then to draw the graph an extra few steps need to be made.

ax = plt.axes (projection='3d')
ax.plot3D (x, y, z, 'blue')

Javascript Implementation…

Javascript’s canvas tag is limited to 2D graphics by default.  WebGL allows one to render 3D graphics within the canvas window, however it is a bit of work in order to make a simple application.  Three.js is a set classes and functions that mask the complexities of WebGL and therefore make it much easier to develop something like a Lorenz equation viewer.

The overall steps to get something working when using Three.js are:

  • Initialize… 
    • Create a camera which establishes the viewpoint, clipping planes, drawing area size, clear colour
    • Create a renderer which is how things will be drawn. For the vast majority of the time, THREE.WebGLRenderer is the renderer to use
    • Create a scene which is a container used to hold all the geometry and lights that are to be drawn. 
      • Three.js contains a bunch of functions for simple 3d geometry/objects (such as cubes and spheres) along with the primitives (points, lines, and triangles). 
      • It also contains a variety of lights (including a point light, spotlight, ambient light).
  • Animate/Update…
    • To view the current state of scene, the renderer’s render function is used. 
    • To make an animation/simulation/game, the render function needs to be called repeatedly while the appropriate updates are made to the scene. Instead of making a while loop, javascript has two options: interval timers and requestAnimationFrame. The requestAnimationFrame method is preferred since it is newer and more efficient. 

The following is some example code that could be used to make the camera, scene, and renderer:

// camera
_3DCamera = new THREE.PerspectiveCamera (75, window.innerWidth / window.innerHeight, 0.1, 300);
_3DCamera.position.set ( defViewX, defViewY, defViewZ );

// renderer
_3DRenderer = new THREE.WebGLRenderer ({ antialias: true });
_3DRenderer.setClearColor (BackGroundColour);

// scene
_3DScene = new THREE.Scene ();

To make a simple line, the following can be used:

var pnts = [];
pnts = [ new THREE.Vector3 (0, 0, 0), new THREE.Vector3 ( 10, 20, 30) ];
const mat = new THREE.LineBasicMaterial ({ color: 0x000000 });
const geoLine = new THREE.BufferGeometry().setFromPoints (pnts);
const axisLine = new THREE.Line (geoLine, mat);

_3DScene.add (axisLine);

To make a line that changes during the simulation, the following needs to be done during initialization:

var positions = new Float32Array (3 * _nMaxPoints);

var geometry = new THREE.BufferGeometry();
geometry.setAttribute ('position', new THREE.BufferAttribute (positions, 3));
geometry.setDrawRange (0, 0);

const material = new THREE.LineBasicMaterial( { color: _CurveColour } );

const solLine = new THREE.Line (geometry, material);

Then during each update, 

var positions = solLine.geometry.attributes.position.array;

for (var i = 0; i < xValues.length; i++)
{
    positions [3*i  ] = xValues [i];
    positions [3*i+1] = yValues [i];
    positions [3*i+2] = zValues [i];
}

solLine.geometry.setDrawRange (0, xValues.length);
solLine.geometry.attributes.position.needsUpdate = true;

Results…

1) Variety of Parameters…

The first example shows the results when with the following parameter values:\sigma = 10 \quad \beta = \frac{8}{3} \quad \rho = \{ 0.5, 13, 22.3, 28, 99.65, 100.75\}

To keep things simple, the following features are available:

  • Show critical points: used to toggle is the critical points are shown or not
  • Rotation speed: how big of an angular step is used during each update. The math… the update angle is the product of the slider value with 0.001 radians
  • Computation rate: how often the calculations are performed. This option results in the speed that the solution curve is updated
  • Curve length: the length of the curve. Sometimes shorter curves are better and other times its all about longer ones.

Also note, that the view window has been set to 450 x 450 to reduce any performance load, which is necessary when dealing with 8 year old hardware.

Different Initial Conditions…

Since the Lorenz attractor is a ‘chaotic system’ where the two main elements are unpredictability of future results and sensitivity of initial conditions. Therefore, for this example the initial starting location can be modified and multiple starting points are derived from that one point. By slightly altering the initial positions, we can see that the solutions are the same at the beginning but eventually diverge and become totally different. The different solutions are still in the same ‘Lorenz attractor domain’.

The few extra features included this time around:

  • Point geometry: which direction that the new points will be created; either in the x, y, or x directions
  • Number of points: The number of points used in the simulation. For now, 1 to 4 additional points can be created which corresponds to 2 to 5 total points
  • Point spacing: Is the delta between each new point
  • Initial position sliders: allow the user to alter the starting point of the initial position of the base point

A few things to note… first the point (0, 0, 0) converges to zero and nothing happens. And second, after any of the initial position parameters are altered, any of the resulting solutions are erased to avoid any confusion. 

Other Lorenz Attractors/Systems…

On the web, there are a few other Lorenz attractors/systems out there:

  • This video has a lot of other details about the butterfly effect and the Lorenz attractor. It doesn’t consider other values of \rho
  • This guy made an analog circuit to solve the Lorenz equations and displays the results on an old CRT style oscilloscope 
  • A lot of cool pictures here. Even includes a cool statue.
  • There are some simulation here and here

Code…

All of the code is available here. It was written using vi (which wasn’t as much fun as I hoped; need to move to a real code editor).

 

No Comments

Add your comment