Attractors

View a running version of this notebook. | Download this project.


Visualizing Attractors

An attractor is a set of values to which a numerical system tends to evolve. An attractor is called a strange attractor if the resulting pattern has a fractal structure. This notebook shows how to calculate and plot two-dimensional attractors of a variety of types, using code and parameters primarily from Lázaro Alonso, François Pacull, Jason Rampe, Paul Bourke, and James A. Bednar.

Clifford Attractors

For example, a Clifford Attractor is a strange attractor defined by two iterative equations that determine the x,y locations of discrete steps in the path of a particle across a 2D space, given a starting point (x0,y0) and the values of four parameters (a,b,c,d):

\begin{equation} x_{n +1} = \sin(a y_{n}) + c \cos(a x_{n})\\ y_{n +1} = \sin(b x_{n}) + d \cos(b y_{n}) \end{equation}

At each time step, the equations define the location for the following time step, and the accumulated locations show the areas of the 2D plane most commonly visited by the imaginary particle.

It's easy to calculate these values in Python using Numba. First, we define the iterative attractor equation:

In [1]:
import numpy as np, pandas as pd, datashader as ds
from datashader import transfer_functions as tf
from datashader.colors import inferno, viridis
from numba import jit
from math import sin, cos, sqrt, fabs

@jit(nopython=True)
def Clifford(x, y, a, b, c, d, *o):
    return sin(a * y) + c * cos(a * x), \
           sin(b * x) + d * cos(b * y)

We then evaluate this equation 10 million times, creating a set of x,y coordinates visited. The @jit here and above is optional, but it makes the code 50x faster.

In [2]:
n=10000000

@jit(nopython=True)
def trajectory_coords(fn, x0, y0, a, b=0, c=0, d=0, e=0, f=0, n=n):
    x, y = np.zeros(n), np.zeros(n)
    x[0], y[0] = x0, y0
    for i in np.arange(n-1):
        x[i+1], y[i+1] = fn(x[i], y[i], a, b, c, d, e, f)
    return x,y

def trajectory(fn, x0, y0, a, b=0, c=0, d=0, e=0, f=0, n=n):
    x, y = trajectory_coords(fn, x0, y0, a, b, c, d, e, f, n)
    return pd.DataFrame(dict(x=x,y=y))
In [3]:
%%time
df = trajectory(Clifford, 0, 0, -1.3, -1.3, -1.8, -1.9)
CPU times: user 1.44 s, sys: 132 ms, total: 1.58 s
Wall time: 1.57 s
In [4]:
df.tail()
Out[4]:
x y
9999995 1.816108 -1.518004
9999996 2.198861 0.040714
9999997 1.675459 -2.176648
9999998 1.334090 0.987108
9999999 -0.665912 -1.525518

We can now aggregate these 10,000,000 continuous coordinates into a discrete 2D rectangular grid with Datashader, counting each time a point fell into that grid cell:

In [5]:
%%time

cvs = ds.Canvas(plot_width = 700, plot_height = 700)
agg = cvs.points(df, 'x', 'y')
print(agg.values[190:195,190:195],"\n")
[[ 34  38  32  43  24]
 [ 25  29  30  34  34]
 [117  30  37  36  29]
 [136 180 117  63  44]
 [ 59  86 132 130  78]] 

CPU times: user 534 ms, sys: 11.5 ms, total: 546 ms
Wall time: 546 ms

A small portion of that grid is shown above, but it's difficult to see the grid's structure from the numerical values. To see the entire array at once, we can turn each grid cell into a pixel, using a greyscale value from white to black:

In [6]:
ds.transfer_functions.Image.border=0

tf.shade(agg, cmap = ["white", "black"])
Out[6]: