Mitch Richling: Lorenz Strange Attractor

Author:	Mitch Richling
Updated:	2023-01-14

1. Strange Attractors
2. The Lorenz Strange Attractor
3. Simulation & Visualization: Digital Computer
4. Simulation & Visualization: Analog Computer
5. References

1. Strange Attractors

First, we need to set the stage with the idea of a dynamical system. Instead of a formal definition, I think a few examples will convey the idea best. Think about the movement of a pendulum, a planet and moon orbiting a star in our solar system, or a bouncing ball on your kitchen floor. An attractor for such a system is the behavior, or set of behaviors, the system will approach after a long enough time. In the case of the bouncing balls, the attractor is quite simple – it is the flat surface of the floor. Attractors don't need to be so simple. Consider for example the orbiting planet. In this case, the attractor might be a complex orbit. Now, a strange attractor is an attractor that is "complicated".

The thing that makes strange attractors so fascinating is that they show the overall patterns of a dynamical system in spite of the fact that the systems described are unbelievably complex. The way a body orbits three or more gravitating bodies is a very good example. Another good example is the strange attractor attributed to Lorenz discussed on this page.

2. The Lorenz Strange Attractor

Lorenz was studying weather prediction, and he developed a rather simple model involving a differential equation in three dimensional Euclidean space. He had great difficulty solving this equation numerically because it is very sensitive to its initial conditions. Each time he would solve the system, even with tiny rounding changes in the 6th digit, he would obtain radically different solutions. This sort of sensitivity has come to be part of the definition of chaotic systems. The differential equation he was working with was:

\[\frac{df}{dt}=[a(y-x),x(b-z)-y,xy-cz]\] \[\frac{df}{dt}=A \circ [y-x,x,-z]-[0,y+xz,-xy]\]

Note \(\circ\) is the Hadamard product (element-wise). Traditional initial conditions are the following:

\[[x, y, z] = [1/10,0,0]\]

While not part of the system proper, this is the most common set of parameters used:

\[A = [a, b, c] = [10,28,8/3]\]

3. Simulation & Visualization: Digital Computer

3.1. Solving & Plotting the Equations with Euler's Method

The picture below is a plot of this equation projected onto the XY-plane. The equation was solved using a short program implementing a very simple version of Euler's method. You can find several example implementations in various languages linked below. Once you have the data, you can generate an image with POV-Ray using any of the example POV-Ray scene files below.

For such a simple graph POV-Ray is overkill, but it demonstrates the basic idea of using a ray tracer to obtain nice plots – like the ones below.

3.2. Solving & Plotting the Equations with Adaptive Euler's Method

Euler's method is the most common numerical ODE algorithm used to solve this equation; however, a more sophisticated solution method is required to obtain high quality plots of the curve. When graphing a curve, the best plots are usually obtained when the sample points used to draw the curve are evenly spaced along the curve. Unfortunately Euler's method generates solution points that are not uniformly spaced along the curve in spite of the fact that the input variable delta is held constant. One option is to modify Euler's method to adaptively change the input variable delta so that the solution points are never too far apart. The resulting programs are a bit more complex – complex enough that I didn't implement a pure POV-Ray version. You can find (C & Perl implementations below. We will still connect the dots with cylinders, but the cylinders will be very short – so the curve will appear to be smooth. The images and movies below were all generated with POV-Ray using one of the example POV-Ray scene files below.

Click on the thumbnail for a larger version!

The adaptive Euler's method above sets tight upper bounds on maximum step size so that we may obtain visually appealing graphs of the Lorenz system. This methodology is, in a way, at direct odds with most modern, advanced numerical ODE algorithms that generally use adaptive step control to make the largest steps possible while still controlling error. By observing the sizes of the steps taken by such algorithms, we can learn a bit about the systems – where they are more and less stable for example. In the image below the solution points produced by a fourth order Runge–Kutta algorithm are visible as blue stripes on the curve. The viewpoint has been changed to more clearly illustrate how tightly packed together the solution points are near the central parts of the spiral.

3.3. Code

3.3.1. Makefile

Makefile used for most of the images on this page.

# -*- Mode:Makefile; Coding:us-ascii-unix; fill-column:132 -*-
####################################################################################################################################
# @file      makefile
# @author    Mitch Richling http://www.mitchr.me
# @Copyright Copyright 1994,1995,1997,2001,2014 by Mitch Richling.  All rights reserved.
# @Revision  $Revision: 1.7 $
# @SCMdate   $Date: 2014/10/11 14:55:27 $
# @brief     Make images and movies.@EOL
# @Keywords  
# @Std       GNUmake
#
#            
#            

#-----------------------------------------------------------------------------------------------------------------------------------

SQUAL = -Q11
#SQUAL = -Q0

SSIZE = -W3840 -H2160
#SSIZE = -W1920 -H1080
#SSIZE = -W960  -H540
#SSIZE = -W480  -H270


MQUAL = -Q11
#MQUAL = -Q5
MSIZE = -W960  -H540
#MSIZE = -W480  -H270
#MSIZE = -W240  -H135

TARGET_L1S = lorenz1_c_s1.png lorenz1_c_s2.png lorenz1_pl_s1.png lorenz1_pl_s2.png lorenz1_pov_s1.png lorenz1_pov_s2.png

TARGET_L2S = lorenz2_c_s1.png lorenz2_c_s2.png lorenz2_pl_s1.png lorenz2_pl_s2.png lorenz2_pl_s3.png lorenz2_pl_s4.png 

TARGET_L2M = lorenz2_pl_s1_movie.mp4 lorenz2_pl_s2_movie.mp4 lorenz2_pl_s3_movie.mp4

TARGET_L2MS = lorenz2_pl_s1_movie_240x135.gif lorenz2_pl_s1_movie_480x270.gif lorenz2_pl_s2_movie_240x135.gif lorenz2_pl_s2_movie_480x270.gif lorenz2_pl_s3_movie_240x135.gif lorenz2_pl_s3_movie_480x270.gif

TARGET_ALL = $(TARGET_L1S) $(TARGET_L2S) $(TARGET_L2M)

################################################################################################################################################################################################################################################################

l1s : $(TARGET_L1S)

l2s : $(TARGET_L2S)

l2m : $(TARGET_L2M)

l2ms : $(TARGET_L2MS)

all : $(TARGET_ALL)

################################################################################################################################################################################################################################################################

lorenz1_c : lorenz1.c
    gcc lorenz1.c -o lorenz1_c

lorenz1_f : lorenz1.f08
    gfortran lorenz1.f08 -o lorenz1_f

lorenz1_c.inc : lorenz1_c
    lorenz1_c > lorenz1_c.inc

lorenz1_c_s1.png : lorenz1_c.inc lorenz_s1.pov
    povray $(SSIZE) $(SQUAL) +A +AM2 -K0.0 +R4 +J3 -P +D  -HIlorenz_s1.pov -Ilorenz1_c.inc -Olorenz1_c_s1.png

lorenz1_c_s2.png : lorenz1_c.inc lorenz_s2.pov
    povray $(SSIZE) $(SQUAL) +A +AM2 -K0.0 +R4 +J3 -P +D  -HIlorenz_s2.pov -Ilorenz1_c.inc -Olorenz1_c_s2.png

lorenz1_pl.inc : lorenz1.pl
    perl lorenz1.pl > lorenz1_pl.inc

lorenz1_pl_s1.png : lorenz1_pl.inc lorenz_s1.pov
    povray $(SSIZE) $(SQUAL) +A +AM2 -K0.0 +R4 +J3 -P +D  -HIlorenz_s1.pov -Ilorenz1_pl.inc -Olorenz1_pl_s1.png

lorenz1_pl_s2.png : lorenz1_pl.inc lorenz_s2.pov
    povray $(SSIZE) $(SQUAL) +A +AM2 -K0.0 +R4 +J3 -P +D  -HIlorenz_s2.pov -Ilorenz1_pl.inc -Olorenz1_pl_s2.png

lorenz1_pov_s1.png : lorenz1.inc lorenz_s1.pov
    povray $(SSIZE) $(SQUAL) +A +AM2 -K0.0 +R4 +J3 -P +D  -HIlorenz_s1.pov -Ilorenz1.inc -Olorenz1_pov_s1.png

lorenz1_pov_s2.png : lorenz1.inc lorenz_s2.pov
    povray $(SSIZE) $(SQUAL) +A +AM2 -K0.0 +R4 +J3 -P +D  -HIlorenz_s2.pov -Ilorenz1.inc -Olorenz1_pov_s2.png

################################################################################################################################################################################################################################################################

lorenz2_c : lorenz2.c
    gcc -lm lorenz2.c -o lorenz2_c

lorenz2_c.inc : lorenz2_c
    lorenz2_c > lorenz2_c.inc

lorenz2_c_s1.png : lorenz2_c.inc lorenz_s1.pov
    povray $(SSIZE) $(SQUAL) +A +AM2 -K0.0 +R4 +J3 -P +D  -HIlorenz_s1.pov -Ilorenz2_c.inc -Olorenz2_c_s1.png

lorenz2_c_s2.png : lorenz2_c.inc lorenz_s1.pov
    povray $(SSIZE) $(SQUAL) +A +AM2 -K0.0 +R4 +J3 -P +D  -HIlorenz_s2.pov -Ilorenz2_c.inc -Olorenz2_c_s2.png

lorenz2_pl.inc : lorenz2.pl
    perl lorenz2.pl > lorenz2_pl.inc

lorenz2_pl_s1.png : lorenz2_pl.inc lorenz_s1.pov
    povray $(SSIZE) $(SQUAL) +A +AM2 -K0.0 +R4 +J3 -P +D  -HIlorenz_s1.pov -Ilorenz2_pl.inc -Olorenz2_pl_s1.png

lorenz2_pl_s2.png : lorenz2_pl.inc lorenz_s2.pov
    povray $(SSIZE) $(SQUAL) +A +AM2 -K0.0 +R4 +J3 -P +D  -HIlorenz_s2.pov -Ilorenz2_pl.inc -Olorenz2_pl_s2.png

lorenz2_pl_s3.png : lorenz2_pl.inc lorenz_s3.pov
    povray $(SSIZE) $(SQUAL) +A +AM2 -K0.0 +R4 +J3 -P +D  -HIlorenz_s3.pov -Ilorenz2_pl.inc -Olorenz2_pl_s3.png

lorenz2_pl_s4.png : lorenz2_pl.inc lorenz_s4.pov
    povray $(SSIZE) $(SQUAL) +A +AM2 -K0.0 +R4 +J3 -P +D  -HIlorenz_s4.pov -Ilorenz2_pl.inc -Olorenz2_pl_s4.png

################################################################################################################################################################################################################################################################

lorenz2_pl_s2_frame120.png : lorenz_s2.pov lorenz2_pl.inc
    povray $(MSIZE) $(MQUAL) +A +R4 +J3 -P +D -KFI1 -KFF120 -KI0 -KF1 -HIlorenz_s2.pov -Ilorenz2_pl.inc -Olorenz2_pl_s2_frame.png

lorenz2_pl_s2_movie.gif : lorenz2_pl_s2_frame120.png
    convert lorenz2_pl_s2_frame*.png lorenz2_pl_s2_movie.gif

lorenz2_pl_s2_movie.mp4 : lorenz2_pl_s2_movie.gif
    avconv -i lorenz2_pl_s2_frame%3d.png lorenz2_pl_s2_movie.mp4

lorenz2_pl_s1_frame120.png : lorenz_s1.pov lorenz2_pl.inc
    povray $(MSIZE) $(MQUAL) +A +R4 +J3 -P +D -KFI1 -KFF120 -KI0 -KF1 -HIlorenz_s1.pov -Ilorenz2_pl.inc -Olorenz2_pl_s1_frame.png

lorenz2_pl_s1_movie.gif : lorenz2_pl_s1_frame120.png
    convert lorenz2_pl_s1_frame*.png lorenz2_pl_s1_movie.gif

lorenz2_pl_s1_movie.mp4 : lorenz2_pl_s1_movie.gif
    avconv -i lorenz2_pl_s1_frame%3d.png lorenz2_pl_s1_movie.mp4

lorenz2_pl_s3_frame120.png : lorenz_s3.pov lorenz2_pl.inc
    povray $(MSIZE) $(MQUAL) +A +R4 +J3 -P +D -KFI1 -KFF120 -KI0 -KF1 -HIlorenz_s3.pov -Ilorenz2_pl.inc -Olorenz2_pl_s3_frame.png

lorenz2_pl_s3_movie.gif : lorenz2_pl_s3_frame120.png
    convert lorenz2_pl_s3_frame*.png lorenz2_pl_s3_movie.gif

lorenz2_pl_s3_movie.mp4 : lorenz2_pl_s3_movie.gif
    avconv -i lorenz2_pl_s3_frame%3d.png lorenz2_pl_s3_movie.mp4

lorenz2_pl_s1_movie_480x270.gif : lorenz2_pl_s1_movie.gif
    convert -sample 480x270 lorenz2_pl_s1_movie.gif lorenz2_pl_s1_movie_480x270.gif 

lorenz2_pl_s2_movie_480x270.gif : lorenz2_pl_s2_movie.gif
    convert -sample 480x270 lorenz2_pl_s2_movie.gif lorenz2_pl_s2_movie_480x270.gif 

lorenz2_pl_s3_movie_480x270.gif : lorenz2_pl_s3_movie.gif
    convert -sample 480x270 lorenz2_pl_s3_movie.gif lorenz2_pl_s3_movie_480x270.gif 

lorenz2_pl_s1_movie_240x135.gif : lorenz2_pl_s1_movie.gif
    convert -sample 240x135 lorenz2_pl_s1_movie.gif lorenz2_pl_s1_movie_240x135.gif 

lorenz2_pl_s2_movie_240x135.gif : lorenz2_pl_s2_movie.gif
    convert -sample 240x135 lorenz2_pl_s2_movie.gif lorenz2_pl_s2_movie_240x135.gif 

lorenz2_pl_s3_movie_240x135.gif : lorenz2_pl_s3_movie.gif
    convert -sample 240x135 lorenz2_pl_s3_movie.gif lorenz2_pl_s3_movie_240x135.gif 

################################################################################################################################################################################################################################################################

clean_intr : 
    rm -rf lorenz1_c lorenz1_c.inc lorenz1_pl.inc lorenz2_c lorenz2_c.inc lorenz2_pl.inc lorenz2_pl_s2_frame*.png lorenz2_pl_s1_frame*.png lorenz2_pl_s3_frame*.png *~ *.pov-state

clean_targets : 
    rm -rf $(TARGET_ALL)

cleanall : clean_targets clean_intr

3.3.2. Euler's Method Solvers

3.3.2.1. Perl

#!/usr/local/bin/perl
# -*- Mode:Perl; Coding:us-ascii-unix; fill-column:132 -*-

####################################################################################################################################
##
# @file      lorenz1.pl
# @author    Mitch Richling http://www.mitchr.me
# @Copyright Copyright 1994,1997,2014 by Mitch Richling.  All rights reserved.
# @Revision  $Revision: 1.1 $
# @SCMdate   $Date: 2014/10/10 19:36:56 $
# @brief     Generate a PovRay include file for the lorenz strange attractor. @EOL
# @Keywords  lorenz strange attractor
# @Std       Perl5
#
# This is a very simple program to generate a curve based on Lorenz's attractor.  It uses Euler's method to solve the equations.
# The "curve" consists of a series of spheres placed at each point found on the curve.  The spheres are connected with cylinders.
# Alternate implimentations: lorenz1.c, lorenz1.inc
#            

#-----------------------------------------------------------------------------------------------------------------------------------

$maxBalls = 10000;
$tDelta = 0.01;
$x = 0.1;
$y = 0.0;
$z = 0.0;
$a = 10.0;
$b = 28.0;
$c = 8.0 / 3.0;
printf(STDERR "Computation starting...\n");
printf("sphere {<%f,%f,%f>, %f texture { crvTex } }\n", $x, $y, $z, 0.2);

for($numBalls=0;$numBalls<$maxBalls;$numBalls++) {
  $xNew = $x + $a*($y-$x)*$tDelta;
  $yNew = $y + ($x*($b-$z)-$y)*$tDelta;
  $zNew = $z + ($x*$y-$c*$z)*$tDelta;

  if(($numBalls % int($maxBalls/20)) == 0) {
    printf(STDERR "Step: %6d  tDelta: %15.3f x: %15.2f y: %15.2f z: %15.2f\n", $numBalls, $tDelta, $x, $y, $z);
  }

  printf("sphere {<%f,%f,%f>, %f texture { crvTex } }\n", $xNew, $yNew, $zNew, 0.2);
  printf("cylinder {<%f,%f,%f>, <%f,%f,%f>, %f texture { crvTex } }\n", $x, $y, $z, $xNew, $yNew, $zNew, 0.2);

  $x=$xNew;
  $y=$yNew;
  $z=$zNew;
}

3.3.2.2. POV-Ray

This implementation is in POV-Ray itself – i.e. we are not using POV-Ray to simply plot some data, but to solve the equations!

// -*- Coding:us-ascii-unix; fill-column:132 -*- 
////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// @file      lorenz1.inc
// @author    Mitch Richling http://www.mitchr.me
// @Copyright Copyright 1994,1997,2014 by Mitch Richling.  All rights reserved.
// @Revision  $Revision: 1.1 $ 
// @SCMdate   $Date: 2014/10/10 19:36:53 $
// @brief     Generate a PovRay include file for the lorenz strange attractor. @EOL
// @Keywords  lorenz strange attractor
// @Std       Povray_3.5
//
// This is a very simple program to generate a curve based on Lorenz's attractor.  It uses Euler's method to solve the equations.
// The "curve" consists of a series of spheres placed at each point found on the curve.  The spheres are connected with cylinders.
// Alternate implimentations: lorenz1.pl, lorenz1.c           

//----------------------------------------------------------------------------------------------------------------------------------

#declare tDelta = 0.01;
#declare xCur = 0.1;
#declare yCur = 0.0;
#declare zCur = 0;
#declare a = 10;
#declare b = 28;
#declare c = 8.0 / 3.0;

#declare numBalls=1;
#while (numBalls < 10000)
    #declare xNew = xCur + a*(yCur-xCur)*tDelta;
    #declare yNew = yCur + (xCur*(b-zCur)-yCur)*tDelta;
    #declare zNew = zCur + (xCur*yCur-c*zCur)*tDelta;

    sphere {<xNew,yNew,zNew>, 0.2 texture { crvTex } }
    cylinder {<xCur,yCur,zCur>, <xNew,yNew,zNew>, 0.2 texture { crvTex } }

    #declare xCur=xNew;
    #declare yCur=yNew;
    #declare zCur=zNew;

    #declare numBalls=numBalls+1;
#end

3.3.2.3. C

/* -*- Mode:C; Coding:us-ascii-unix; fill-column:132 -*- */
/**********************************************************************************************************************************/
/**
   @file      lorenz1.c
   @author    Mitch Richling http://www.mitchr.me
   @Copyright Copyright 1994,1997,2014 by Mitch Richling.  All rights reserved.
   @brief     Generate a PovRay include file for the lorenz strange attractor. @EOL
   @Std       C89

   This is a very simple program to generate a curve based on Lorenz's attractor.  It uses Euler's method to solve the equations.
   The "curve" consists of a series of spheres placed at each point found on the curve.  The spheres are connected with cylinders.
   Alternate implimentations: lorenz1.pl, lorenz1.inc

***********************************************************************************************************************************/

#include <stdio.h>

int main() {
  int maxBalls  = 10000;
  double tDelta = 0.01;
  double x      = 0.1;
  double y      = 0.0;
  double z      = 0.0;
  double a      = 10.0;
  double b      = 28.0;
  double c      = 8.0 / 3.0;
  int numBalls;
  double xNew, yNew, zNew;
  fprintf(stderr, "Computation starting...\n");
  printf("sphere {<%f,%f,%f>, %f texture { crvTex } }\n", x, y, z, 0.2);

  for(numBalls=0;numBalls<maxBalls;numBalls++) {
    xNew = x + a*(y-x)*tDelta;
    yNew = y + (x*(b-z)-y)*tDelta;
    zNew = z + (x*y-c*z)*tDelta;

    /* Display 20 "status" messages. */
    if((numBalls % (int)(maxBalls/20)) == 0)
      fprintf(stderr, "Step: %6d  tDelta: %15.3f x: %15.5f y: %15.5f z: %15.5f\n", numBalls, tDelta, x, y, z);

    printf("sphere {<%f,%f,%f>, %f texture { crvTex } }\n", xNew, yNew, zNew, 0.2);
    printf("cylinder {<%f,%f,%f>, <%f,%f,%f>, %f texture { crvTex } }\n", x, y, z, xNew, yNew, zNew, 0.2);

    x=xNew;
    y=yNew;
    z=zNew;
  }

  return 0;
}

3.3.2.4. Fortan 2008

! -*- Mode:F90; Coding:us-ascii-unix; fill-column:132 -*-
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!
! @file      lorenz1.f08
! @author    Mitch Richling <https://www.mitchr.me>
! @date      2021-05-18
! @brief     Generate a PovRay include file for the lorenz strange attractor. @EOL
! @std       F2003
!
! This is a very simple program to generate a curve based on Lorenz's attractor.  It uses Euler's method to solve the equations.
! The "curve" consists of a series of spheres placed at each point found on the curve.  The spheres are connected with cylinders.
! Alternate implimentations: lorenz1.pl, lorenz1.c, lorenz1.inc
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

PROGRAM lorenz1

  USE, INTRINSIC :: ISO_FORTRAN_ENV

  IMPLICIT NONE

  INTEGER, PARAMETER                    :: maxBalls = 10000
  REAL(REAL64), PARAMETER               :: tDelta   = 0.01_REAL64;
  REAL(REAL64), PARAMETER, DIMENSION(3) :: A        = [10.0_REAL64, 28.0_REAL64, 8.0_REAL64 / 3.0_REAL64 ]
  REAL(REAL64), DIMENSION(3)            :: x        = [ 0.1_REAL64,  0.0_REAL64,              0.0_REAL64 ], xNew
  INTEGER                               :: numBalls

  WRITE(ERROR_UNIT, fmt=*) "Computation starting..."
  WRITE(OUTPUT_UNIT, fmt='("sphere {<", f0.3, ",", f0.3, ",", f0.3, ">, 0.2 texture { crvTex } }")') x(1), x(2), x(3)
  DO numBalls=0,(maxBalls-1)

     ! This is a more natural translation of the equations into vector form.  
     ! This rearrangement of the order of the arithmetic steps leads to very different results by iteration 3000
     ! That's the essence of sensitive dependence (chaos)!!
     ! xNew = x + tDelta * (A * [x(2)-x(1), x(1), -x(3)] - [0.0_REAL64, x(2)+x(1)*x(3), -x(1)*x(2)])

     ! This is precisely the same sequence of arithmetic as used in the C, Perl, & POV-Ray code:
     xNew = x + [ A(1)*(x(2)-x(1)), x(1)*(A(2)-x(3))-x(2), x(1)*x(2)-A(3)*x(3) ] * tDelta

    ! Display 20 status messages.
    IF (MODULO(numBalls, maxBalls/20) == 0) THEN
       WRITE(ERROR_UNIT, fmt='(i6, f15.5, 3(f15.10))') numBalls, tDelta, x
    END IF

    WRITE(OUTPUT_UNIT, fmt='("sphere {<", f0.3, ",", f0.3, ",", f0.3, ">, 0.2 texture { crvTex } }")') xNew(1), xNew(2), xNew(3)
    WRITE(OUTPUT_UNIT, fmt='("cylinder {<", 2(f0.3, ","), f0.3, ">, <", 2(f0.3, ","), f0.3, ">, 0.2 texture { crvTex } }")') &
         x(1), x(2), x(3), xNew(1), xNew(2), xNew(3)

    x=xNew        
  END DO

END PROGRAM lorenz1

3.3.3. Adaptive Euler's Method Solvers

3.3.3.1. Perl

#!/usr/local/bin/perl
# -*- Mode:Perl; Coding:us-ascii-unix; fill-column:132 -*-

####################################################################################################################################
##
# @file      lorenz2.pl
# @author    Mitch Richling http://www.mitchr.me
# @Copyright Copyright 1994,1997,2001,2014 by Mitch Richling.  All rights reserved.
# @Revision  $Revision: 1.1 $
# @SCMdate   $Date: 2014/10/10 19:40:31 $
# @brief     Generate a PovRay include file for the lorenz strange attractor. @EOL
# @Keywords  lorenz strange attractor povray adaptive Euler
# @Std       Perl5
#
# This is an adaptive Euler's method to solve Lorenz's differential equations.  The value of delta f is bounded above and below so
# that the curve looks smooth when we are done.  The bounding is adjusted by adjusting t delta in each leap-frog step.  The value of
# delta t is bounded absolutely above and below to keep things from going nuts.  The method of finding an acceptable t delta is to
# use bisection between the bounds on t delta.  This is a strange adaptation of Euler's method, but it works well for this
# application. Alternate implimentation: lorenz2.c
#            

#-----------------------------------------------------------------------------------------------------------------------------------

$distToGo = 2500;                  # How long should the curve be?

$sphereRad = 0.2;                  # Size of the spheres to use

$maxMovDelta   = $sphereRad / 4;   # The maximum delta for the function.
$minMovDelta   = $sphereRad / 5;   # The minimum delta for the function.

$minTdelta        = 0.000001;      # Minimum value for tDelta
$maxTdelta        = 0.01;          # Maximum value for tDelta 
$tDeltaZeroThresh = 0.00000001;    # Episolon ball for zero.
$maxNumBisect     = 10;            # Max bisections on tDelta.

# Initial conditions, and constants.
$x = 0.1;
$y = 0.0;
$z = 0;
$a = 10;
$b = 28;
$c = 8.0 / 3;

# Solve the equations.....
$tDelta       = $maxTdelta;
while($dist < $distToGo) {
  $numBalls++;
  # Take a step up so that we minimize the number of balls.
  $tDelta = ($maxTdelta + $tDelta)/2;
  if($tDelta > $maxTdelta) { $tDelta = $maxTdelta; }
  $curMaxTdelta = $maxTdelta;
  $curMinTdelta = $minTdelta;
  # Bisect t delta until we have done it too much or
  # until we have a good value for t delta.
  $numBisect = 0;
  $doneBisecting = 0;
  while( !($doneBisecting)) {
    $numBisect++;
    $Xdelta = $a*($y-$x)*$tDelta;
    $Ydelta = ($x*($b-$z)-$y)*$tDelta;
    $Zdelta = ($x*$y-$c*$z)*$tDelta;
    $movDelta = sqrt(abs($Xdelta**2+$Ydelta**2+$Zdelta**2));
    if($numBisect > $maxNumBisect) {
      $doneBisecting = 1;
    } else {
      if($movDelta > $maxMovDelta) {
        if(abs($tDelta-$curMinTdelta)<$tDeltaZeroThresh) {
          $doneBisecting = 1;
        } else {
          $curMaxTdelta = $tDelta;
          $tDelta = ($tDelta + $curMinTdelta)/2;
          if($tDelta < $minTdelta) { $tDelta = $minTdelta; }
        }
      } elsif($movDelta < $minMovDelta) {
        if(abs($tDelta-$curMaxTdelta)<$tDeltaZeroThresh) {
          $doneBisecting = 1;
        } else {
          $curMinTdelta = $tDelta;
          $tDelta = ($curMaxTdelta + $tDelta)/2;
          if($tDelta > $maxTdelta) { $tDelta = $maxTdelta; }
        }
      } else {
        $doneBisecting = 1;
      }
    }
  }
  $dist += $movDelta;
  $xOld = $x;
  $yOld = $y;
  $zOld = $z;
  $x = $x + $Xdelta;
  $y = $y + $Ydelta;
  $z = $z + $Zdelta;
  if($movDelta>$maxMovDelta) {
    $moveCh = '*';
    $numOver++;
  } else {
    $moveCh = ' ';
  }

  #printf(STDERR "Balls: %7d tDel: %9.6f x: %9.2f y: %9.2f z: %9.2f Dist: %7.1f fDel: %9.3f #Over: %7d %1s-%1s\n", $numBalls, $tDelta, $x, $y, $z, $dist, $movDelta, $numOver, $moveCh, ($numBisect>$maxNumBisect?'*':' '));
  print "sphere {<$x,$y,$z>, $sphereRad texture { crvTex } }\n"; 
  print "cylinder {<$x,$y,$z>, <$xOld,$yOld,$zOld>, $sphereRad texture { crvTex } }\n";

}
printf(STDERR "Balls: %7d Dist: %7.1f #Over: %7d\n", $numBalls, $dist, $numOver);

3.3.3.2. C

/* -*- Mode:C; Coding:us-ascii-unix; fill-column:132 -*- */
/**********************************************************************************************************************************/
/**
   @file      lorenz2.c
   @author    Mitch Richling http://www.mitchr.me
   @Copyright Copyright 1994,1997,2001,2014 by Mitch Richling.  All rights reserved.
   @Revision  $Revision: 1.2 $ 
   @SCMdate   $Date: 2014/10/10 19:40:05 $
   @brief     Generate a PovRay include file for the lorenz strange attractor. @EOL
   @Keywords  lorenz strange attractor povray adaptive Euler
   @Std       C89

   This is an adaptive Euler's method to solve Lorenz's differential equations.  The value of delta f is bounded above and below so
   that the curve looks smooth when we are done.  The bounding is adjusted by adjusting t delta in each leap-frog step.  The value
   of delta t is bounded absolutely above and below to keep things from going nuts.  The method of finding an acceptable t delta is
   to use bisection between the bounds on t delta.  This is a strange adaptation of Euler's method, but it works well for this
   application.  Alternate implimentation lorenz2.pl.

***********************************************************************************************************************************/

#include <stdio.h>
#include <math.h>

#define DEBUG 0
#define CYL 1

int main () {  
  double distToGo = 2500.0;             /* How long should the curve be? */

  double sphereRad = 0.2;               /* Size of the spheres to use */

  double maxMovDelta = sphereRad / 4.0; /* The maximum delta for the function. */
  double minMovDelta = sphereRad / 5.0; /* The minimum delta for the function. */


  double minTdelta        = 0.000001;   /* Minimum value for tDelta */
  double maxTdelta        = 0.01;       /* Maximum value for tDelta */
  double tDeltaZeroThresh = 0.00000001; /* Episolon ball for zero. */
  double maxNumBisect     = 10;         /* Max bisections on tDelta. */

  /* Initial conditions, and constants. */
  double x = 0.1;
  double y = 0.0;
  double z = 0;
  double a = 10;
  double b = 28;
  double c = 8.0 / 3;

  double curMaxTdelta, curMinTdelta, tDelta, dist, Xdelta, Ydelta, Zdelta, movDelta, xOld, yOld, zOld;
  int    numBisect, doneBisecting, numBalls, numOver;
  char   moveCh;

 /*  Solve the equations..... */
  numBalls = 0;
  tDelta = maxTdelta;
  dist = 0.0;
  numOver = 0;
  while (dist < distToGo) {
    numBalls++;
    /*  Take a big step up to minimize the number of balls. */
    tDelta = (maxTdelta + tDelta) / 2;
    if (tDelta > maxTdelta) {
      tDelta = maxTdelta;
    }
    curMaxTdelta = maxTdelta;
    curMinTdelta = minTdelta;
    /*  Bisect t delta until we have done it too much or until we have a good value for t delta. */
    numBisect = 0;
    doneBisecting = 0;
    while (!(doneBisecting)) {
      numBisect++;
      Xdelta = a * (y - x) * tDelta;
      Ydelta = (x * (b - z) - y) * tDelta;
      Zdelta = (x * y - c * z) * tDelta;
      movDelta = sqrt (fabs (Xdelta * Xdelta + Ydelta * Ydelta + Zdelta * Zdelta));
      if (numBisect > maxNumBisect) {
        doneBisecting = 1;
      } else {
        if (movDelta > maxMovDelta) {
          if (fabs (tDelta - curMinTdelta) < tDeltaZeroThresh) {
            doneBisecting = 1;
          } else {
            curMaxTdelta = tDelta;
            tDelta = (tDelta + curMinTdelta) / 2;
            if (tDelta < minTdelta) {
              tDelta = minTdelta;
            }
          }
        } else if (movDelta < minMovDelta) {
          if (fabs (tDelta - curMaxTdelta) < tDeltaZeroThresh) {
            doneBisecting = 1;
          } else {
            curMinTdelta = tDelta;
            tDelta = (curMaxTdelta + tDelta) / 2;
            if (tDelta > maxTdelta) {
              tDelta = maxTdelta;
            }
          }
        } else {
          doneBisecting = 1;
        }
      }
    }
    dist += movDelta;
    xOld = x;
    yOld = y;
    zOld = z;
    x = x + Xdelta;
    y = y + Ydelta;
    z = z + Zdelta;
    if (movDelta > maxMovDelta) {
      moveCh = '*';
      numOver++;
    } else {
      moveCh = ' ';
    }
#if DEBUG
    fprintf (stderr, "%d Balls: %7d tDel: %9.6f x: %9.2f y: %9.2f z: %9.2f Dist: %7.1f fDel: %9.3f #Over: %7d %1c-%1c\n",
             numBisect, numBalls, tDelta, x, y, z, dist, movDelta, numOver, moveCh, (numBisect > maxNumBisect ? '*' : ' '));
#endif
    printf ("sphere {<%f,%f,%f>, %f texture { crvTex } }\n", x, y, z, sphereRad);
#if CYL
    printf("cylinder {<%f,%f,%f>, <%f,%f,%f>, %f texture { crvTex } }\n", x, y, z, xOld, yOld, zOld, sphereRad);
#endif

  }
  fprintf (stderr, "Balls: %7d Dist: %7.1f #Over: %7d\n", numBalls, dist, numOver);

  return 0;
}

3.3.4. POV-Ray scene files

////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// @file      lorenz_s1.pov
// @author    Mitch Richling http://www.mitchr.me
// @Copyright Copyright 1994,1997,2014 by Mitch Richling.  All rights reserved.
// @Revision  $Revision: 1.4 $ 
// @SCMdate   $Date: 2014/10/11 03:38:26 $
// @brief     Basic setup for rendering the Lorenz strange attractor.@EOL
// @Keywords  lorenz strange attractor povray movie 
// @Std       Povray_3.7
//
//            
//            

//----------------------------------------------------------------------------------------------------------------------------------

#version 3.7;

#include "colors.inc"
#include "textures.inc"
#include "woods.inc"
#include "metals.inc"
#include "skies.inc"

global_settings { assumed_gamma 1 }

camera {
    location <cos(2*pi*clock)*-25-5, -5, sin(2*pi*clock)*45+30>
    sky      <0,1,0>
    look_at  <0,0,25>
}

light_source { < 40,-4,  40>                                     color White*0.3 }
light_source { <-40,-50,  40>                                    color White*0.6 }  // back shadow
light_source { < 40,-4, -40>                                     color White*0.3 }
light_source { <-40,-4, -40>                                     color White*0.8 } // Shadow on left

background { color Black }

#declare crvTex=texture { pigment { color Red } }

////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// @file      lorenz_s2.pov
// @author    Mitch Richling http://www.mitchr.me
// @Copyright Copyright 1994,1997,2014 by Mitch Richling.  All rights reserved.
// @Revision  $Revision: 1.4 $ 
// @SCMdate   $Date: 2014/10/11 03:38:38 $
// @brief     Basic setup for rendering the Lorenz strange attractor.@EOL
// @Keywords  lorenz strange attractor povray movie 
// @Std       Povray_3.7
//
//            
//            

//----------------------------------------------------------------------------------------------------------------------------------

#version 3.7;

#include "colors.inc"
#include "textures.inc"
#include "woods.inc"
#include "metals.inc"
#include "skies.inc"

global_settings { assumed_gamma 1 }

camera {
    location <cos(2*pi*clock)*-30, -5, sin(2*pi*clock)*45+30>
    sky      <0,1,0>
    look_at  <0,0,25>
}

light_source { < 40,-4,  40>                                     color White*0.3 }
light_source { <-40,-50,  40>                                    color White*0.6 }  // back shadow
light_source { < 40,-4, -40>                                     color White*0.3 }
light_source { <-40,-4, -40>                                     color White*0.8 } // Shadow on left

background { color Black }

#declare crvTex=texture {
   pigment { color Red }
   finish  {
      ambient .15
      diffuse 0.01
      reflection 0.4
      specular 0.9
      roughness 0.04
      phong .1 
      phong_size 200
    }
}

////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// @file      lorenz_s3.pov
// @author    Mitch Richling http://www.mitchr.me
// @Copyright Copyright 1994,1997,2014 by Mitch Richling.  All rights reserved.
// @Revision  $Revision: 1.4 $ 
// @SCMdate   $Date: 2014/10/11 03:38:43 $
// @brief     Basic setup for rendering the Lorenz strange attractor.@EOL
// @Keywords  lorenz strange attractor povray movie 
// @Std       Povray_3.7
//
//            
//            

//----------------------------------------------------------------------------------------------------------------------------------

#version 3.7;

#include "colors.inc"
#include "textures.inc"
#include "glass.inc"
#include "metals.inc"
#include "skies.inc"

global_settings { assumed_gamma 1 }

camera {
    location <cos(2*pi*clock)*-30, -5, sin(2*pi*clock)*45+30>
    sky      <0,1,0>
    look_at  <0,0,25>
}

light_source { < 40,-4,  40>                                     color White*0.3 }
light_source { <-40,-50,  40>                                    color White*0.6 }  // back shadow
light_source { < 40,-4, -40>                                     color White*0.3 }
light_source { <-40,-4, -40>                                     color White*0.8 } // Shadow on left

background { color Black }

#declare crvTex=texture { Candy_Cane }

////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// @file      lorenz_s4.pov
// @author    Mitch Richling http://www.mitchr.me
// @Copyright Copyright 1994,1997,2014 by Mitch Richling.  All rights reserved.
// @Revision  $Revision: 1.4 $ 
// @SCMdate   $Date: 2014/10/11 03:38:50 $
// @brief     Basic setup for rendering the Lorenz strange attractor.@EOL
// @Keywords  lorenz strange attractor povray movie 
// @Std       Povray_3.7
//
//            
//            

//----------------------------------------------------------------------------------------------------------------------------------

#version 3.7;

#include "colors.inc"
#include "textures.inc"
#include "glass.inc"
#include "metals.inc"
#include "skies.inc"

global_settings { assumed_gamma 1 }

camera {
    location <cos(2*pi*clock)*-30, -5, sin(2*pi*clock)*45+30>
    sky      <0,1,0>
    look_at  <0,0,25>
}

light_source { < 40,-4,  40>                                     color White*0.3 }
light_source { <-40,-50,  40>                                    color White*0.6 }  // back shadow
light_source { < 40,-4, -40>                                     color White*0.3 }
light_source { <-40,-4, -40>                                     color White*0.8 } // Shadow on left

background { color Black }

#declare crvTex=texture { Candy_Cane }
#declare crvTex=texture { T_Copper_4E }

////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// @file      lorenz_s4.pov
// @author    Mitch Richling http://www.mitchr.me
// @Copyright Copyright 1994,1997,2014 by Mitch Richling.  All rights reserved.
// @Revision  $Revision: 1.4 $ 
// @SCMdate   $Date: 2014/10/11 03:38:50 $
// @brief     Basic setup for rendering the Lorenz strange attractor.@EOL
// @Keywords  lorenz strange attractor povray movie 
// @Std       Povray_3.7
//
//            
//            

//----------------------------------------------------------------------------------------------------------------------------------

#version 3.7;

global_settings { assumed_gamma 1 }

#include "colors.inc"

#declare lookAt=<0,0,0>;

camera {
    location <0,0,-30>
    sky      <0,1,0>
    look_at  lookAt
}

cylinder { <-1000,0,0>, lookAt, 0.20 pigment { color Pink } }  
cylinder { <0,-1000,0>, lookAt, 0.20 pigment { color Pink } }  
cylinder { lookAt, <1000,0,0>, 0.20 pigment  { color Blue } }   
cylinder { lookAt, <0,1000,0>, 0.20 pigment  { color Yellow } } 

light_source { <0,0,-50> color White }

background { color Black }

#declare crvTex=texture { pigment { color Red } }

////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// @file      lorenz_s4.pov
// @author    Mitch Richling http://www.mitchr.me
// @Copyright Copyright 1994,1997,2014 by Mitch Richling.  All rights reserved.
// @Revision  $Revision: 1.4 $ 
// @SCMdate   $Date: 2014/10/11 03:38:50 $
// @brief     Basic setup for rendering the Lorenz strange attractor.@EOL
// @Keywords  lorenz strange attractor povray movie 
// @Std       Povray_3.7
//
//            
//            

//----------------------------------------------------------------------------------------------------------------------------------

#version 3.7;

#include "colors.inc"
#include "textures.inc"
#include "glass.inc"
#include "metals.inc"
#include "skies.inc"

global_settings { assumed_gamma 1 }

camera {
    location <cos(2*pi*clock)*-30, -5, sin(2*pi*clock)*45+30>
    sky      <0,1,0>
    look_at  <0,0,25>
}

light_source { < 40,-4,  40>                                     color White*0.3 }
light_source { <-40,-50,  40>                                    color White*0.6 }  // back shadow
light_source { < 40,-4, -40>                                     color White*0.3 }
light_source { <-40,-4, -40>                                     color White*0.8 } // Shadow on left

background { color Black }

#declare crvTex=texture { Jade }

3.3.5. Resizing Images

#!/bin/bash
# -*- Mode:Shell-script; Coding:us-ascii-unix; fill-column:132 -*-
####################################################################################################################################
##
# @file      conv.sh
# @author    Mitch Richling <https://www.mitchr.me>
# @date      2021-05-21
# @brief     Convert png & gif files into various sizes for web use.@EOL
# @std       bash
#
####################################################################################################################################

for f in *.png;      do echo $f; convert -quality 100 -resize 240x135   $f `echo $f | sed 's/\.png$/_240x135.jpg/'`;   done
for f in *.png;      do echo $f; convert -quality 100 -resize 480x270   $f `echo $f | sed 's/\.png$/_480x270.jpg/'`;   done
for f in *.png;      do echo $f; convert -quality 100 -resize 960x540   $f `echo $f | sed 's/\.png$/_960x540.jpg/'`;   done
for f in *.png;      do echo $f; convert -quality 100 -resize 1920x1080 $f `echo $f | sed 's/\.png$/_1920x1080.jpg/'`; done
for f in *.png;      do echo $f; convert -quality 100                   $f `echo $f | sed 's/\.png$/_3840x2160.jpg/'`; done
for f in *movie.gif; do echo $f; cp                                     $f `echo $f | sed 's/\.gif$/_960x540.gif/'`;   done
for f in *movie.gif; do echo $f; convert              -resize 480x270   $f `echo $f | sed 's/\.gif$/_480x270.gif/'`;   done
for f in *movie.gif; do echo $f; convert              -resize 120x67    $f `echo $f | sed 's/\.gif$/_120x67.gif/'`;    done

You can find all the images here.

4. Simulation & Visualization: Analog Computer

Given any system of differential equations it is possible to find an electronic circuit governed by the same equations. These circuits are called "analogs" of the system of equations, and that's why we use the word "analog" for most non-digital electronics today. Finding and implementing these analogs in a systematic way is the fundamental idea behind Analog computers.

4.1. Circuit

This circuit is from "The Art of Electronics" by Paul Horowitz.

4.2. SPICE Simulation Results

4.3. Physical Implementation

4.4. Behavior of Physical Circuit

Note the output of the circuit is passed through another analog processing unit that translates, scales, and rotates the results – that's why the image is spinning. Someday I'll do a little write-up of that circuit…

Scope: Tektronix AN/USM488 (military version of the 2235). Equipped with the standard clear plexiglass CRT filter
Capture Device: Olympus OM-D E-M1 MARK II with Olympus M.Zuiko ED 60mm f/2.8 Macro Lens
Post: ffmpeg (Encode to VP9, crop, & gamma correction)

Scope: Tektronix TDS3052B Digital Phosphor Oscilloscope
Capture Device: GRACETOP USB VGA Capture Card connected to a Raspberry PI 4 and converted with VLC
Post: ffmpeg (Encode to VP9)

Scope: Siglent SDS2504X Plus Super Phosphor Oscilloscope
Capture Device: Scope's VNC server output connected to a Raspberry PI 4 and converted with VLC
Post: ffmpeg (Encode to VP9)

5. References

Check out the chaos in electronics section of my reading list, as well as the "The Art of Electronics" by Paul Horowitz. For more about chaos in general, take a look at the fractals section of my reading list.