Tutorials‎ > ‎

Week 6: Fractals




In this tutorial we learn how to create beautiful images of fractals. Along the way, we also will learn a little about complex numbers.








Introduction

In this week's tutorial we will use Python to generate images of fractals. The images below were made  using the pygame package.

 
 
 
 
 
 
 
 
 
(The colour scheme here has been "borrowed" from a nice pygame code developed by Ian Mallett.)

The top left image shows a visualization of the full Mandelbrot set. The subsequent images show the effect of "zooming in on" (magnifying) a certain region. Hopefully you can see that smaller versions of the initial shape in yellow reappear in the later images.  This is called self-similarity, and self-similarity is a typical property of fractals. (Note that the little copies of the body of the Mandelbrot set are all slightly different, mostly because of the thin threads connecting them to the main body of the set.)


Complex Numbers


To understand the algorithm for generating fractals like the Mandelbrot set, we first need to know a little bit about complex numbers. Let's take a quick mathematical detour.

So far in these tutorials, we have mostly been working with the set of natural numbers, 1, 2, 3, etc., (and 0). The natural numbers are part of a larger set: the integers, ..., -2, -1, 0, 1, 2, ..., etc.. This set  includes both positive and negative whole numbers. The integers are, in turn, a subset of the rational numbers: numbers which can be expressed as fractions, e.g. 3/4, 355/113, -2/5, etc. The rational numbers are, in turn, a subset of the real numbers. A real number is a value that represents a quantity along a continuous line. The set of real numbers includes even those numbers which cannot be expressed as fractions. For example, both pi and the square root of 2 are real but not rational. In fact, there are an infinite number of real numbers between every pair of distinct rational numbers.

But can we go further? Yes. An imaginary number is a number that can be written as a real number multiplied by the imaginary unit i, which is defined by the property i2 = -1. We may say that "i is the square root of -1".  Clearly, i cannot be a real number, as the square of any real number is positive.

A complex number is any number that can be expressed in the form

a + i b

where a and b are real numbers, and i is the imaginary unit. (Note that we don't bother writing the "times" symbol between the i and b).  We say that "a" is the real part, and "b" is the imaginary part, of this complex number.

Here is a Venn diagram showing the hierarchy of the number sets (see mathsisfun.com for explanations):


  • N, the set of natural numbers, 1, 2, 3, etc., is contained within
  • Z, the set of integers, ... -2, -1, 0, 1, 2, ..., which is contained within
  • Q, the set of rational numbers, e.g. 1/2, 2/5, 355/113, etc., which is contained within
  • R is the set of real numbers.
  • I is the set of imaginary numbers, which does not intersect with R (except possibly at 0).
  • C is the set of complex numbers.
Complex numbers have many beautiful properties, and turn out to be very useful for solving problems in physics and engineering (e.g. for understanding resonances in atoms, electronic circuits, or bridges).

Addition and Multiplication

For this tutorial, all we need to know is how to add and multiply. Let's suppose we have two complex numbers, a + ib and c + id.  Then their sum is given by:

(a + ib) + (c + id) = (a + c)  +  i (b + d)

In other words, the real part of the sum is a + c and the imaginary part is b + d

Their product is given by:

(a + i b)(c + i d) = a c + i b c + i a d + i2 b d = (a c - b d) + i (b c + a d)

In other words, the real part is a c - b d, and the imaginary part is b c + a d.

Let's try a simple example. Let U = 1 + 2 i and V = 3 - 4 i. Then

 U + V = (1 + 3) + i (2 - 4)   =  4 - 2 i.

and

 U V = (3 + 8) + (6 - 4) i     = 11 + 2 i.


Coordinates

Complex numbers are often represented as points in the complex plane, with the real part giving the x coordinate, and the imaginary part giving the y coordinate.

Modulus

The modulus (also known as "absolute value") of a complex number z, denoted |z|, is found from the square root of the sum of squares of real and imaginary parts. If z = x + i y, then

|z|2 = x2 + y2

In other words, to find the modulus of the complex number we simply apply Pythagoras theorem in the xy plane.

The Mandelbrot Set

There are various of types of fractals. We will only look at the most famous example, described by (and named after) Benoit Mandelbrot.

Let's start by considering a point in the xy plane, which is described by two coordinates, x and y, which are real numbers. We can make a complex number c from the coordinates as follows:

    c = x + i y

Now let's introduce a simple function "f", which acts on complex numbers:

    f(z) = z2 + c.

This function squares its argument (z) and then adds c. Here z and c are both complex numbers. We can square a complex number by multiplying it by itself, according to the rule above (see "Addition and Multiplication").

Now let's think about what happens when we iterate this function -- in other words, what happens when we repeatedly apply this function to some starting value. Mandelbrot considered the sequence generated by repeatedly applying the function to the starting value z=0:

[ f(0), f(f(0)), f(f(f(0))), ....] = [ c, c2 + c, (c2 + c)2 + c, ... ]

He asked a simple question: How does this sequence behave as we iterate many times? There are two possibilities, either:
  1. the sequence remains bounded (the absolute value of the terms in the sequence remains smaller than a certain finite value), or,
  2. the sequence is unbounded (the absolute value of the terms increases without bound).

The Mandelbrot set consists of only those points (x, y) which generate bounded sequences [i.e. condition (1)].  It turns out that the boundary of this set has truly remarkable fractal properties. 

A Compact Set

The Mandelbrot set is only strictly defined for points within a disk of radius 2 around the origin, i.e., for

x2 + y2 < 22.

In fact, if the absolute value of any term in the sequence exceeds 2, then the sequence is guaranteed to increase without bound.

Colouring the fractal

To generate a fractal image we will iterate the function for each pixel on the screen, and use the number of iterations to set the colour of the pixel. We can do this by applying the recurrence relation:

zn+1 = (zn)2 + c

where n is the iteration number, and c = x + i y is generated by the pixel position (x,y). Obviously, we can't carry out an infinite number of iterations in our program. Instead, we will:
  • choose a maximum number of iterations nmax, and a maximum absolute value d = 2.
  • for each pixel, apply the recurrence relation until either:
  •      (1)  |zn| > d,          or   (2)  n > nmax .
  • Now we will use the number of iterations n to choose the colour of the pixel.  To do this, we should apply a colour map: a way of converting an integer (in this case, the number of iterations) into a colour for the pixel.

In my code, I used nmax = 255. The choice of colour map is up to you. I have used a simple but effective colour map which I "borrowed" from this code. It is simply a list of colours, given as Red/Green/Blue values:

ColourPalette = [(241, 233, 191), (248, 201, 95), (255, 170, 0), (204, 108, 0), (153, 87, 0), (106, 52, 3), (66, 30, 15), (25, 7, 26), (25, 7, 26), (9, 1, 47), (4, 4, 73), (0, 7, 100), (12, 44, 138), (24, 82, 177), (57, 125, 209), (134, 181, 229), (211, 236, 248)]

As the number of function iterations may exceed the number of colour values, we can use the remainder operation, i.e.:

numberOfColours = len(ColourPalette)
pixelColour = ColourPalette[ iterations % numberOfColours ]
Note that for the interior regions of the fractal, the number of iterations will be nmax=255, and the number of colours is 17, and since 255 % 17 = 0, the interior will be coloured with the first entry in the palette, which is yellow.

The Program

A step-by-step algorithm for generating an image from the Mandelbrot set is given on the Wikipedia page. If you are feeling confident, you may wish to try implementing the algorithm yourself now.

I have written a pygame code that implements the Wikipedia algorithm. It can be found at the bottom of this page (mandelbrot1.py). I have added many comments to the code, so I hope you can understand how it works. Try running it : Once the image appears on screen, you may "zoom in" by clicking the left mouse button anywhere on the image. Wait a while, and a new image appears. You may zoom in as many times as you wish.

Please try it for yourself! Let me know if it doesn't work for you (pi3challenge@sheffield.ac.uk).

Speed

Calculating fractals is rather slow. On my laptop, it takes around 6 seconds to render each image, and it is sure to take longer on the Raspberry Pi.  To increase the speed, you may wish to :

  • Reduce the resolution, i.e., reduce res = 100 to a smaller value
  • Reduce the maximum number of iterations (decrease max_it = 255 to a smaller value).

(If you have any tips for optimizing the code, please let me know.)

Exercises

  1. What interesting features can you find hidden in the fractal?
  2. Try changing the colour palette -- can you make the images more aesthetically appealing?
  3. Modify the code so that is saves the image after each magnification (look for the commented line). Can you create an animated gif like this one, from the sequence of images?
  4. Can you modify the code so that it draws each pixel on screen as soon as it is computed?
ċ
mandelbrot1.py
(5k)
Sam Dolan,
8 May 2013, 00:33