Simple Mandelbrot Set In Python

I've looked at other questions regarding this, but I can't seem to figure out where I am going wrong. The goal is to 'Write a program to make an image of the Mandelbrot set by perf

Solution 1:

There are a couple of problems with your code.

Firstly, you are using xvalues and yvalues to index M, but those indices should be pixel index integers in the range 0..(n-1). Numpy will convert the floats in xvalues and yvalues to integers, but the resulting numbers will be in -2..2, so there won't be many pixels plotted, the image will be tiny, and you'll get wrapping due to the way negative indices work in Python.

A simple way to get the required pixel indices is to use the built-in Python function enumerate, but there may be a way to re-organize your code to do this with Numpy functions.

The other problem is that you need to reset z to zero for each pixel. Currently, your code reuses whatever the last z was for the previous pixel, and if that pixel was in the Mandelbrot set then z will be too big.

Here's a repaired version of your code. I don't have pylab, so I wrote a simple bitmap viewer using PIL. You can save the image to a file by calling img.save(filename) inside my show function; PIL will figure out the correct file format from the filename extension.

import numpy as np
from PIL import Image

def show(data):
    img = Image.frombytes('1', data.shape[::-1], np.packbits(data, 1))
    img.show()

n = 100
maxiter = 100

M = np.zeros([n, n], np.uint8)
xvalues = np.linspace(-2, 2, n)
yvalues = np.linspace(-2, 2, n)

for u, x in enumerate(xvalues):
    for v, y in enumerate(yvalues):
        z = 0
        c = complex(x, y)
        for i in range(maxiter):
            z = z*z + c
            if abs(z) > 2.0:
                M[v, u] = 1
                break

show(M)

Here's the output image:

---

Of course, whenever you find yourself iterating over Numpy array indices, that's a sign that You're Doing It Wrong. The main point of using Numpy is that it can perform operations on whole arrays at once by internally iterating over them at C speed; the above code might as well be using plain Python lists instead of Numpy arrays.

Here's a version that gets Numpy to do most of the looping. It uses more RAM, but it's roughly 2.5 times faster than the previous version, and it's somewhat shorter.

This code uses several 2D arrays. c contains all the complex seed numbers, we carry out the core Mandelbrot calculation in z, which is initialsed to zeros. mask is a boolean array that controls where the Mandelbrot calculation needs to be performed. All its items are initially set to True, and on each iteration True items in mask that correspond to items in z that have escaped from the Mandelbrot set are set to False.

To test if a point has escaped we use z.real**2 + z.imag**2 > 4.0 rather than abs(z) > 2.0, this save a little bit of time because it avoids the expensive square root calculation and the abs function call.

We can use the final value of mask to plot the Mandelbrot set, but to make the points in the Mandelbrot set black we need to invert its values, which we can do with 1 - mask.

import numpy as np
from PIL import Image

def show(data):
    img = Image.frombytes('1', data.shape[::-1], np.packbits(data, 1))
    img.show()
    img.save('mset.png')

n = 100
maxiter = 100

a = np.linspace(-2, 2, n)
c = a + 1.j * a[:, None]
z = np.zeros((n, n), np.complex128)
mask = np.ones((n, n), np.bool)

for i in range(maxiter):
    mask[mask] = z[mask].real**2 + z[mask].imag**2 < 4.0
    z[mask] = z[mask]**2 + c[mask]

show(1 - mask)

Solution 2:

from pylab import imshow,show,gray
from numpy import zeros,linspace
N=1000
A=zeros([N+1,N+1],float)
defmandler(x,y):
    c=complex(x,y)
    z=c
    n=0while n<N:
        z=z*z+c
        ifabs(z)>2:
            return0
        n+=1else:
        return1for x in linspace(-2,2,N):
    for y in linspace(2,-2,N):
        a=mandler(x,y)
        s=int(((x+2)*N)/4)
        t=int(((y+2)*N)/4)
        A[s,t]=a
imshow(A)
gray()
show()