import numba
import numpy as np
from pylab import show
from time import time
import blz
from shutil import rmtree

@numba.njit
def mandel(x, y, max_iters):
    """
        Given the real and imaginary parts of a complex number,
        determine if it is a candidate for membership in the Mandelbrot
        set given a fixed number of iterations.
    """
    c = complex(x, y)
    z = 0.0j
    for i in xrange(max_iters):
        z = z*z + c
        if (z.real*z.real + z.imag*z.imag) >= 4:
            return i

    return max_iters


def create_fractal(height, width, min_x, max_x, min_y, max_y, image, iters):

    pixel_size_x = (max_x - min_x) / width
    pixel_size_y = (max_y - min_y) / height

    for x in xrange(height):

        imag = min_y + x * pixel_size_y

        for y in xrange(width):

            real = min_x + y * pixel_size_x
            color = mandel(real, imag, iters)
            image[x, y] = color


%matplotlib inline
height = 1024
width = 1536

image = np.zeros((height, width), dtype=np.uint8)

t1 = time()
create_fractal(height, width, -2.0, 1.0, -1.0, 1.0, image, 20)
t2 = time()

elapsed1 = t2-t1
print elapsed1
imshow(image)

@numba.njit
def mandel(x, y, max_iters):
    """
        Given the real and imaginary parts of a complex number,
        determine if it is a candidate for membership in the Mandelbrot
        set given a fixed number of iterations.
    """
    c = complex(x, y)
    z = 0.0j
    for i in xrange(max_iters):
        z = z*z + c
        if (z.real*z.real + z.imag*z.imag) >= 4:
            return i

    return max_iters


def create_fractal(height, width, min_x, max_x, min_y, max_y, image, row, iters):

    pixel_size_x = (max_x - min_x) / width
    pixel_size_y = (max_y - min_y) / height

    for x in xrange(height):

        imag = min_y + x * pixel_size_y

        for y in xrange(width):

            real = min_x + y * pixel_size_x
            color = mandel(real, imag, iters)
            row[y] = color

        image.append(row)


height = 1024
width = 1536

#If the blz already exist, remove it
rmtree('images/Mandelbrot.blz', ignore_errors=True)

image = blz.zeros((0, width), rootdir='images/Mandelbrot.blz', dtype=np.uint8,
                  expectedlen=height*width,
                  bparams=blz.bparams(clevel=9, shuffle=True, cname='zlib'))
row = np.zeros((width), dtype=np.uint8)

t1 = time()
create_fractal(height, width, -2.0, 1.0, -1.0, 1.0, image, row, 20)
t2 = time()

image.flush()
elapsed2 = t2-t1
print elapsed2
imshow(image)

print str(elapsed1/elapsed2) + ' times faster'

print 'Size in memory: %s' % image.nbytes
print 'Size in disk: %s' % image.cbytes
print 'Compress ratio: %f' % (float(image.nbytes)/float(image.cbytes))

height = 20000
width = 30000

#If the blz already exist, remove it
rmtree('images/Mandelbrot.blz', ignore_errors=True)

image = blz.zeros((0, width), rootdir='images/Mandelbrot.blz', dtype=np.uint8,
                  expectedlen=height*width,
                  bparams=blz.bparams(clevel=9, shuffle=True, cname='zlib'))
row = np.zeros((width), dtype=np.uint8)

t1 = time()
create_fractal(height, width, -2.0, 1.0, -1.0, 1.0, image, row, 20)
t2 = time()

image.flush()
elapsed3 = t2-t1
print elapsed3

print 'Size in memory: %s' % image.nbytes
print 'Size in disk: %s' % image.cbytes
print 'Compress ratio: %f' % (float(image.nbytes)/float(image.cbytes))
We can compute images as fast as with regular methods but we get compression and we don't run out of memory.

blz is awesome, numba is awesome you should use them.