This is one of the 100 recipes of the IPython Cookbook, the definitive guide to high-performance scientific computing and data science in Python.
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
Our discrete dynamical system is defined by the recursive application of the logistic function:
$$x_{n+1}^{(r)} = f_r(x_n^{(r)}) = rx_n^{(r)}(1-x_n^{(r)})$$def logistic(r, x):
return r*x*(1-x)
n = 10000
r = np.linspace(2.5, 4.0, n)
iterations = 1000
last = 100
x = 1e-5 * np.ones(n)
lyapunov = np.zeros(n)
plt.figure(figsize=(6,7));
plt.subplot(211);
for i in range(iterations):
x = logistic(r, x)
# We compute the partial sum of the Lyapunov exponent.
lyapunov += np.log(abs(r-2*r*x))
# We display the bifurcation diagram.
if i >= (iterations - last):
plt.plot(r, x, ',k', alpha=.04)
plt.xlim(2.5, 4);
plt.title("Bifurcation diagram");
# We display the Lyapunov exponent.
plt.subplot(212);
plt.plot(r[lyapunov<0], lyapunov[lyapunov<0] / iterations,
',k', alpha=.2);
plt.plot(r[lyapunov>=0], lyapunov[lyapunov>=0] / iterations,
',r', alpha=.5);
plt.xlim(2.5, 4);
plt.ylim(-2, 1);
plt.title("Lyapunov exponent");
plt.tight_layout();
The bifurcation diagram brings out the existence of a fixed point for $r<3$, then two and four equilibria... until a chaotic behavior when $r$ belongs to certain areas of the parameter space. We observe an important property of the Lyapunov exponent: it is positive when the system is chaotic (in red here).
You'll find all the explanations, figures, references, and much more in the book (to be released later this summer).
IPython Cookbook, by Cyrille Rossant, Packt Publishing, 2014 (500 pages).