Generating huge Mandelbrot fractals with Python, the quick and memory-safe way¶
Let's see how custom large generated images behave when working with blz and numba. For this purpose we are going to generate some Mandelbrot fractals directly to a blz container.
import numba
import numpy as np
from pylab import show
from time import time
import blz
from shutil import rmtree
We will be using the next code to generate the fractals.
@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
Let's generate a small fractal so we can se how it works.
%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)
2.01292991638
<matplotlib.image.AxesImage at 0x971f350>
That was good but, what happens if I want to generate a really large fractal? Will I run out of memory? Sadly, the answer is yes.
But this is where blz shines, can we modify the code so we generate the image directly to a blz container? Happily, the answer is yes.
@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)
Let's take a look to the previous code, create_fractal now receives a new parameter 'row'. This parameter will allow us to generate a complete row and the append it to the disk blz, so we never run out of memory. This is a really fast method, way better than generating the whole blz and assigning values to the matrix.
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)
2.01285195351
<matplotlib.image.AxesImage at 0xaca22d0>
As we can see, it is as fast as the previous version and it doesn't eat our RAM.
print str(elapsed1/elapsed2) + ' times faster'
1.00003873254 times faster
We have also compressed it.
print 'Size in memory: %s' % image.nbytes
print 'Size in disk: %s' % image.cbytes
print 'Compress ratio: %f' % (float(image.nbytes)/float(image.cbytes))
Size in memory: 1572864 Size in disk: 581160 Compress ratio: 2.706422
What if we want to generate a really big image? Let's try to generate a 30000x20000 Mandelbrot's fractal. This will take a while.
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
658.535179853
print 'Size in memory: %s' % image.nbytes
print 'Size in disk: %s' % image.cbytes
print 'Compress ratio: %f' % (float(image.nbytes)/float(image.cbytes))
Size in memory: 600000000 Size in disk: 4351807 Compress ratio: 137.873761
So we have just generated a 30000x20000 Mandelbrot's fractal that takes 572.2046 MiB of RAM, compressed and gotten a 4.2 MiB without using more than 100 MiB of RAM. That was pretty awesome.