Introduction to the Chaos library

Calling Sequence:

readlib(chaos)
function ( args )

Description:

To use a chaos function, you must first execute the readlib(chaos) command. Most of the terminology used can be found in Devaney's undergraduate text, Chaotic Dynamical Systems (ISBN: 0-201-55406-2), Frame and Peak's Chaos Under Control (ISBN:0-7167-2429-4), or Barnsley's Fractals Everywhere (ISBN: 0-12-079062-9).

The functions available are:

Basic Dynamics

Orbit( f:function,x:initial value,n:positive integer )

- produces a sequence of the first n exact values for the f -orbit of x. Returns a list.

Orbitf( f:function , x:initial value , n:number of iterations )

- computes the first values of the orbit of under using floating point approximations. Returns a list.

GraphicalAnalysis( f:function,a:lower bound,b:upper bound, :seed,n:iterations )

- displays a graphical analysis of n iterations of the orbit of under f over the interval

AnimatedGraphicalAnalysis(same arguments as GraphicalAnalysis)

- performs a graphical analysis one step at a time and then animates the sequence. (The animation window operates like a VCR)

FixedPointAnalysis ( F:polynomial function )

- tries to find and classifies all fixed points and 2 cycles for real valued polynomial function F

Bifurcation

Bifurcation ( f:function, :lower limit, :upper limit,hpoints:horizontal points, n:number of points, :real number )

- plots the bifurcation diagram of for c in with increment . It plots n points after throwing away the first 50 points. It uses an initial value of for the individual orbits.

Bif3d ( , , , , xgrid,ygrid,toss,keep )

- Computes an extension of the bifurcation diagram in above the Mandelbrot set to the rectangular grid in the complex plane whose corners are and over a grid which is xgrid rows and ygrid columns. It throws away the first toss elements of each orbit and plots the next keep values using the absolute value of the complex values as the third coordinate.

Data Analysis

TimeSeriesPlot ( f , t , n )

- plots the time series plot of the first n iterates of the f -orbit of t , i.e. plots the points [0, t ],[1, f ( t )],...,[ n , ( t )].

ShowMixing ( f , t , n )

-visually tests the mixing property of first n iterates the f -orbit of t by plotting the points [0.5, t ],[0.5, f ( t )],...,[0.5, ( t )]

FirstDelayPlot ( List )

-plots the first delay plot of a list of data values, List, i.e. plots [ , ],[ , ],...,[ , ].

ClosePairsPlot ( DATA , tolerance )

- plots the close-pairs plot of DATA (a list of real values) with tolerance tolerance.

ChaosGameTest ( DATA )

- plays the 4 corner chaos game driven by DATA , which is a list of real values between 0 and 1.

CircleChaosGameTest ( DATA )

- plays a variant of the chaos game driven by DATA , which is a list of real values between 0 and 1. Instead of using a fixed number of goal points and clustering the data into groups, the points on the unit circle are all goal points and are labeled from 0 to 1 CCW direction and each value is plotted by moving half way from the previous point to the point on the unit circle that is labeled by the current data value.

BarnsInt ( DATA,d )

- returns an IFS object (see below) whose attractor is the Barnsley Fractal Interpolation IFS for the data points in DATA, which is a list of n ordered pairs [[ , ],...,[ , ]]. d is a list of n -1 parameters that influence the fractal dimension of the graph.

Particular Maps

Dbl ( x:real )

- The doubling function

LogisticChart ( :real,n:pos. integer )

- makes a nice table of the first n iterates of various orbits for logistic function

LogisticGrowthPlot ( c , , n , L )

-Plots the time series plot for the first n iterates of the f -orbit of logistic map with nice formatting. If L='lines' then it connects the data points, otherwise it doesn't.

BoxLogisticGrowthPlot ( c)

-Plots the time series plots for the first 50 iterates of the f -orbits of 0.1,0.2,...,0.8, and 0.9 under the logistic map all on the same graph to see if they all have the same long term behavior or not.

ChaosGame ( PtList:List of Points,Start:starting point,NumPts:number of points,Ratio:ratio to move )

- plays the chaos game starting with point Start and continuing by moving Ratio percent of the distance toward a randomly chosen point from PtList . It does this NumPts times. Points are input as lists, [ x,y ].

IFS's and AFFINE maps

In the following routines, we define two data types: AFFINE and IFS.

An IFS is a single AFFINE or a list of one or more AFFINE's.

An AFFINE represents an affine map of the plane and comes in four flavors:

affine(a,b,c,d,e,f) - standard form
Affine(R,S,theta,phi,E,F) - geometric form
affineC( , , ) - complex form
affineM(M,B) - matrix form

These are inert forms representing the following maps:

affine(a,b,c,d,e,f) <-> T( x , y )=( ax + by + e , cx + dy + f )

Affine(R,S,theta,phi,E,F) <-> T( x,y )=( , )

affineC( , , ) <-> T(z)= where are complex numbers

affineM( M,b ) <-> where M is a 2x2 matrix and b,v are column vectors

IFS routines

`convert\affine`(A::AFFINE) - converts an AFFINE or IFS object to its affine() form

`convert\Affine`(A::AFFINE) - converts an AFFINE or IFS object to its Affine() form

`convert\affineC`(A::AFFINE) - converts an AFFINE or IFS object to its affineC() form

`convert\affineM`(A::AFFINE) - converts an AFFINE or IFS object to its affineM() form

affineFromPoints ( , , , , , )

- Computes the AFFINE object affine( a,b,c,d,e,f) for the unique affine map T( x , y )=( ax + by + e , cx + dy + f ) which maps points to , to , and to .

IFSFromList ( List )

-creates an IFS data structure from a list of m lists of the form [ R, S, , , E , F ]

Map (IFS)

- converts the AFFINE or IFS object to the map T( x , y )=( ax + by + e , cx + dy + f )

Transform(IFS)

- converts the AFFINE or IFS object to a transform which can be applied to plot objects

DrawDetIFS( figure,IFS,n )

- plots the n th iteration of the deterministic method for IFS starting with figure

DrawIFS( IFS,n )

- plots the attractor of IFS by the random iteration method using n points. Probability weights can be passed as an optional third argument. (Note: don't use a value of n that is larger than 30000)

IFScurve ( t,IFS )

- computes the point in the attractor of an m transformation IFS whose address begins with the digits in the base m expansion of t

DrawIFScurve ( IFS,n )

- plots the attractor defined by IFS using n evenly spaced points in [0..1] for addresses

ContractionFactor ( a,b,c,d,e,f )

-returns the contraction factor for the affine transformation T( x , y )=( ax + by + e , cx + dy + f )

Built-in IFS's

The following IFS data types are built-in:

SierpinskiCarpet SierpinskiTriangle SierpinskiTriangle2

KochCurve CantorSet

In addition IFS's can be quickly produced with the following routine:

GridIFS()

- this routine takes arguments, each of which is from the set . It returns the following IFS: Divide the unit square [0..1]x[0..1] into an n x n grid. Number the grid squared from 1 to starting in the lower left corner and moving from left to right and then bottom to top. The i th argument describes an affine transformation which maps the unit square into the i th grid cell in a manner analogous to one of the symmetry operations in the dihedral group of the square (e.g. Lt means it is rotated 90 degrees CCW, Rt means 90deg CW, -Up is a reflection along a vertical axis, -Lt is reflection along the vertical axis followed by 90deg CCW rotation, etc). If the argument is none, there is no affine map associated with that grid cell. For example, the right Sierpinski Triangle can be produced with GridIFS(Up,Up,Up,none)

Built-in starting figures

The following plot structures are built-in, and can be used as the starting figure for plotting the iterations of the deterministic method for an IFS (see DrawDetIFS()).

MrFace MrGrid MrPoint MrLine

MrSquare MrTriangle MrsTriangle

Mandelbrot and Julia type fractals
JuliaFractal( Fx,Fy,a..b,c..d,options )

- computes a Julia-set type fractal for a function F : C -> C where C is the complex plane. The arguments are as follows:

Fx - the real part of F
Fy - the real part of F
a..b - an optional range indicating the horizontal range to be displayed
(the default is -2..2)
c ..d - an optional range indicating the vertical range to be displayed
(the default is -2..2)

Options:

The remaining options given as equations of the form option = value. The following options are supported:

maxiter - a positive integer indicating the maximum number of iterations to compute for each point before declaring that point to be in the filled in Julia set associated with F (default is 500)

radius - a positive real indicating the escape radius (default is 2)

grid - a list, [n,m], of integers indicating the number of grid points horizontally and vertically to use (default is [100,100])

scheme - an integer from 1..6 indicating which color scheme to use (default is 1)

rate - a number from 1 to maxiter indicating the rate at which colors should change. Low numbers will produce a striped coloration while higher values produce a gradual transition (default is 25)

timer - a boolean indicating whether the routine should report how long it took to compute the fractal (default is false )

Any additional arguments are passed along to plots[display] before rendering.

MandelFractal( Fx,Fy,a..b,c..d,options )

- computes a Mandelbrot set type fractal for a family of functions F : C -> C where C is the complex plane. The arguments are as follows:

Fx - the real part of F (must be a proc(x,y,c0,c1) returning a real number)
(default value is (x,y,c0,c1)->x*x-y*y+c0 )
Fy - the real part of F (must be a proc(x,y,c0,c1) returning a real number)
(default value is (x,y,c0,c1)->2*x*y+c1 )
a..b - an optional range indicating the horizontal range to be displayed
(the default is -2..2)
c ..d - an optional range indicating the vertical range to be displayed
(the default is -2..2)

Options:

The remaining options given as equations of the form option = value. The following options are supported:

maxiter - a positive integer indicating the maximum number of iterations to compute for each point before declaring that point to be in the Mandelbrot set associated with F
(default is 500)

radius - a positive real indicating the escape radius (default is 2)

grid - a list, [n,m], of integers indicating the number of grid points horizontally and vertically to use (default is [100,100])

scheme - an integer from 1..6 indicating which color scheme to use (default is 1)

rate - a number from 1 to maxiter indicating the rate at which colors should change. Low numbers will produce a striped coloration while higher values produce a gradual transition (default is 25)

timer - a boolean indicating whether the routine should report how long it took to compute the fractal (default is false )

Any additional arguments are passed along to plots[display] before rendering.

MakePalette( n,rate,scheme )

- creates a color palette with n colors and period rate. There are six schemes available numbered 1..6.

ViewPalette( palette )

- displays the color palettes that are output by MakePalette()

ComplexToRealFunc( F )

- returns the real and imaginary parts of a complex function F : C -> C

MandComplexToRealFunc( F )

- returns the real and imaginary parts of a complex family of functions F : C -> C parameterized by c0+c1*I

QuadraticFunc( c )

- returns the map z->z^2+c where c is a complex number

JuliaSet( c )

- shorthand for JuliaFractal(ComplexToRealFunc(QuadraticFunc( c ))). Additional args are passed along to JuliaFractal().

MandelbrotSet()

- shorthand for MandelbrotFractal(). Additional args are passed along to MandelbrotFractal().

Other Useful Utilities

Table()

- Takes lists as args and prints them in a table.

BaseN ( t,n )

-Converts a real number t in [0..1] to base n ... returns a list of Digits digits.

FatPoint ( x:real,y:real,n:pos. integer,siz:pos. real,col:color constant )

- plots a big fat point at (x,y) as a regular polygon with n sides and radius siz with color col. Used to make special points more visible in a plot.

Examples

restart;

readlib(chaos); #assuming you have installed it in your library

f:=x->x^2-2;
Q:=(x,c)->x^2+c;
L:=(x,c)->c*x*(1-x);

Orbit(f,1/2,10);

Orb:=Orbitf(f,1/2,6);

Table(Orb,numbered);

1   .50
2 -1.75
3  1.06
4  -.87
5 -1.24
6  -.46
7 -1.79 

Orb:=Orbit(x->Q(x,-0.78+0.04*I),0,10):
Table([Orb,opts(digits=4)],map(abs,Orb));

                            0   0 
                 -.78+.4e-1*I .78 
             -.1732-.224e-1*I .17 
     -.75050352+.4775936e-1*I .75 
-.2190254229-.3168713559e-1*I .22 
-.7330319387+.5388057655e-1*I .74 
-.2455672934-.3899236697e-1*I .25 
-.7212171091+.5915050004e-1*I .72 
-.2633446632-.4532070528e-1*I .27 
-.7127035547+.6386993174e-1*I .72 
-.2761330113-.5104065478e-1*I .28 
GraphicalAnalysis(f,-2,2,.65,250);

AnimatedGraphicalAnalysis(f,-2,2,.65,50);

FixedPointAnalysis(x->x^2-1);

*******************************************************************************
x^2-1
-------------------------------------------------------------------------------
  Fixed Point      F'(p)         F''(p)          F'''(p)      Classification
-------------- -------------- -------------- -------------- -------------------
      1.618034       3.236068       2.000000       0.000000 repelling
      -.618034      -1.236068       2.000000       0.000000 repelling

(x^2-1)^2-1
-------------------------------------------------------------------------------
   2-cycle         F'(p)         F''(p)          F'''(p)      Classification
-------------- -------------- -------------- -------------- -------------------
      0.000000       0.000000      -4.000000       0.000000 attracting
     -1.000000       0.000000       8.000000     -24.000000 attracting
      1.618034      10.472136      27.416408      38.832816 repelling
      -.618034       1.527864        .583592     -14.832816 repelling
*******************************************************************************

Bifurcation(Q,-2,0.25,500,100,0);

Bif3d(-2,-2,2,2,40,40,150,10);

TimeSeriesPlot(x->x^2-0.75,0,100);

ShowMixing(x->3.95*x*(1-x),0.5,200);

n:=1000:
Random:=rand(1..10^11):
List1:=[seq(evalf(Random()/10^11),i=1..n)]: #random data
List2:=Orbit(x->3.99*x*(1-x),.328763415,n-1): #chaotic data

# output omitted: Table(List1,List2,numbered):

FirstDelayPlot(List1);

FirstDelayPlot(List2);

ClosePairsPlot(List1[1..100],0.05);

ClosePairsPlot(List2[1..100],0.05);

ChaosGameTest(List1);

ChaosGameTest(List2);

CircleChaosGameTest(List1);

CircleChaosGameTest(List2);

mb:=BarnsInt([[1,2],[2,3],[3,1],[4,2]],[.4,.4,.4]);

IFSCurve(0,mb);
IFSCurve(1/3,mb);
IFSCurve(2/3,mb);
IFSCurve(1,mb);

IFSCurve(0.5,mb);

DrawIFSCurve(mb,400);

plot(Dbl(x),x=0..1,discont=true);

LogisticChart(3.2,25):

Iter: 0 0.000 .100 .200 .300 .400 .500 .600 .700 .800 .900 1.000 
Iter: 1 0.000 .288 .512 .672 .768 .800 .768 .672 .512 .288 0.000 
Iter: 2 0.000 .656 .800 .705 .570 .512 .570 .705 .800 .656 0.000 
Iter: 3 0.000 .722 .513 .665 .784 .800 .784 .665 .513 .722 0.000 
Iter: 4 0.000 .642 .799 .713 .541 .513 .541 .713 .799 .642 0.000 
Iter: 5 0.000 .735 .513 .655 .795 .799 .795 .655 .513 .735 0.000 
Iter: 6 0.000 .623 .799 .723 .522 .513 .522 .723 .799 .623 0.000 
Iter: 7 0.000 .752 .513 .641 .798 .799 .798 .641 .513 .752 0.000 
Iter: 8 0.000 .598 .799 .737 .515 .513 .515 .737 .799 .598 0.000 
Iter: 9 0.000 .770 .513 .621 .799 .799 .799 .621 .513 .770 0.000 
Iter:10 0.000 .567 .799 .753 .513 .513 .513 .753 .799 .567 0.000 
Iter:11 0.000 .785 .513 .595 .799 .799 .799 .595 .513 .785 0.000 
Iter:12 0.000 .539 .799 .771 .513 .513 .513 .771 .799 .539 0.000 
Iter:13 0.000 .795 .513 .565 .799 .799 .799 .565 .513 .795 0.000 
Iter:14 0.000 .521 .799 .787 .513 .513 .513 .787 .799 .521 0.000 
Iter:15 0.000 .799 .513 .537 .799 .799 .799 .537 .513 .799 0.000 
Iter:16 0.000 .515 .799 .796 .513 .513 .513 .796 .799 .515 0.000 
Iter:17 0.000 .799 .513 .520 .799 .799 .799 .520 .513 .799 0.000 
Iter:18 0.000 .513 .799 .799 .513 .513 .513 .799 .799 .513 0.000 
Iter:19 0.000 .799 .513 .515 .799 .799 .799 .515 .513 .799 0.000 
Iter:20 0.000 .513 .799 .799 .513 .513 .513 .799 .799 .513 0.000 
Iter:21 0.000 .799 .513 .513 .799 .799 .799 .513 .513 .799 0.000 
Iter:22 0.000 .513 .799 .799 .513 .513 .513 .799 .799 .513 0.000 
Iter:23 0.000 .799 .513 .513 .799 .799 .799 .513 .513 .799 0.000 
Iter:24 0.000 .513 .799 .799 .513 .513 .513 .799 .799 .513 0.000 
Iter:25 0.000 .799 .513 .513 .799 .799 .799 .513 .513 .799 0.000 

LogisticGrowthPlot(3.5,.3,100,lines);

BoxLogisticGrowthPlot(3.5);

ChaosGame([[0,0],[0,1],[1,0]],[0,0],10000,0.5);

ChaosGame([[0,0],[1,1],[1,0],[0,1]],[0.3,0.3],10000,.55);

T:=affineFromPoints([0,0],[1,0],[0,1],[0,1],[0,1/2],[1/2,1/2]);

convert(T,Affine);

convert(T,affineC);

convert(T,affineM);

Map(T);

Mapf(T);

P:=plot(sin(x),x=-6..6):P;

display(Transform(T)(P));

ContractionFactor(T);

MyIFS:=GridIFS(
-Up,none,none,-Lt,
none,Up,Lt,none,
none,Rt,Dn,none,
-Rt,none,none,-Dn
);

Map(MyIFS);

Map(MyIFS)[1];

Map(MyIFS)[1](2,3);

convert(MyIFS,affine);

convert(MyIFS,affineC);

convert(MyIFS,affineM);

ContractionFactor(MyIFS);

DrawDetIFS(MrFace,MyIFS,1);

DrawDetIFS(MrFace,MyIFS,2);

DrawIFS(MyIFS,20000);

NiceCurveIFS:=IFSFromList(
[[-1/2,1/2,-120,-120,0,0],
[1/2,1/2,0,0,1/4,sqrt(3)/4],
[-1/2,1/2,120,120,3/4,sqrt(3)/4]]);

IFSCurve(0.5,NiceCurveIFS);

DrawIFSCurve(NiceCurveIFS,3^6);

DrawIFS(SierpinskiCarpet,30000);

DrawIFS(SierpinskiTriangle,20000);

DrawIFS(SierpinskiTriangle2,20000);

DrawIFS(KochCurve,10000);

DrawIFS(CantorSet,5000);

DrawDetIFS(MrFace,GridIFS(Up,none,-Up,none,Lt,none,Rt,none,Lt),1,scaling=constrained);

DrawDetIFS(MrGrid,GridIFS(Up,none,-Up,none,Lt,none,Rt,none,Lt),1,scaling=constrained);

DrawDetIFS(MrPoint,GridIFS(Up,none,-Up,none,Lt,none,Rt,none,Lt),1,axes=none,scaling=constrained);

DrawDetIFS(MrLine,GridIFS(Up,none,-Up,none,Lt,none,Rt,none,Lt),1,axes=none,scaling=constrained);

DrawDetIFS(MrSquare,GridIFS(Up,none,-Up,none,Lt,none,Rt,none,Lt),1,scaling=constrained);

DrawDetIFS(MrTriangle,GridIFS(Up,none,-Up,none,Lt,none,Rt,none,Lt),1,scaling=constrained);

DrawDetIFS(MrsTriangle,GridIFS(Up,none,-Up,none,Lt,none,Rt,none,Lt),1,scaling=constrained);

JuliaSet(-1);

JuliaSet(-1,scheme=6,rate=2);

JuliaSet(-1,1.3..1.75,-0.15..0.15,maxiter=50,grid=[50,50],scheme=2,timer=true);

That took 0 minutes and .6 seconds to compute.

JuliaSet(-1,1.3..1.75,-0.15..0.15,maxiter=400,grid=[150,150],scheme=2,timer=true);

That took 0 minutes and 36.0 seconds to compute.

JuliaSet(-0.7519+0.03523*I);

ComplexToRealFunc(z->sin(z));

JuliaFractal(ComplexToRealFunc(z->sin(z)),-6..6,-6..6,radius=20,maxiter=200,grid=[150,150],scheme=6,rate=15,timer=true);

That took 0 minutes and 44.0 seconds to compute.

MandelbrotSet();

MandelbrotSet(timer=true);

That took 0 minutes and 10.0 seconds to compute.

MandelbrotSet(-0.7519..-0.75158,0.03523..0.035563,maxiter=10000,grid=[300,300],scheme=1,rate=1000,timer=true,axes=normal);

That took 6 minutes and 44.7 seconds to compute.

MandComplexToRealFunc(z->sin(z)+c0+c1*I);

MandelFractal(MandComplexToRealFunc(z->sin(z)+c0+c1*I),-6..6,-6..6,radius=20,maxiter=400,grid=[150,150],scheme=6,timer=true,axes=framed);

That took 2 minutes and 6.8 seconds to compute.

MakePalette(10,10,1);

ViewPalette(MakePalette(10,10,1));

BaseN(1/5,4);

The next three commands are not part of the chaos library but are useful for saving fractals images to postscript or jpeg files instead of displaying them on the screen.

#plotsetup(ps,plotoptions="colour=rgb,width=8in,height=6in",plotoutput="C:/testpic.ps");

#plotsetup(jpeg,plotoptions="width=800,height=600",plotoutput="C:/testpic.jpg");

#plotsetup(default,plotoutput=default);

See Also:

utils

Many mathematics files on this site are in pdf format. If your browser does not display pdf files, click here for assistance.
This page was last  updated on Sunday, May 28, 2000 10:13:52 PM. © Ken Monks