# Lorenz Attractor (projected onto XY-plane)
# http://en.wikipedia.org/wiki/Lorenz_attractor
# FB - 201107317
import random
from PIL import Image
imgx = 800
imgy = 600
image = Image.new("RGB", (imgx, imgy))
maxIt = 100000 # number of pixels to draw
size = 30
xa = -size
xb = size
ya = -size
yb = size
# initial state
x = random.random() * size * 2 - 1
y = random.random() * size * 2 - 1
z = random.random() * size * 2 - 1
# dx/dt = delta * (y - x)
# dy/dt = r * x - y - x * z
# dz/dt = x * y - b * z
delta = float(10) # Prandtl number
r = float(28)
b = float(8) / 3
h = 1e-3 # time step
def Lorenz(x, y, z):
dx_dt = delta * (y - x)
dy_dt = r * x - y - x * z
dz_dt = x * y - b * z
x += dx_dt * h
y += dy_dt * h
z += dz_dt * h
return (x, y, z)
for i in range(maxIt):
(x, y, z) = Lorenz(x, y, z)#; print x, y, z
xi = int((imgx - 1) * (x - xa) / (xb - xa))
yi = int((imgy - 1) * (y - ya) / (yb - ya))
if xi >=0 and xi < imgx and yi >= 0 and yi < imgy:
image.putpixel((xi, yi), (255, 255, 255))
image.save("Lorenz_Attractor.png", "PNG")
