Friday, January 15, 2010

Profiling adventures and cython - setting the stage

Summary This post is the first in a series dedicated to
examining the use of profiling and, eventually, using cython, as a
means to improve greatly the speed of an application. The intended
audience is for programmers who have never done any profiling and/or
never used cython before. Note that we will not make use of cython
until the third post in this series.

Preamble

Python is a great multi-purpose language which is really fun to use. However, it is sometimes too slow for some applications. Since I only program for fun, I had never really faced a situation where I found Python's speed to be truly a limiting factor - at least, not until a few weeks ago when I did some exploration of a four-colouring grid problem I talked about. I started exploring ways to speed things up using only Python and trying to come up with different algorithms, but every one I tried was just too slow. So, I decided it was time to take the plunge and do something different. After considering various alternatives, like using shedskin or
attempting to write a C extension (which I found too daunting since I don't know C), I decided to try to use cython.

cython, for those that don't know it, is a Python look-alike language that claims

to make writing C extensions for the Python languageas easy as Python itself.

After looking at a few examples on the web, I concluded that such a rather bold statement might very well be true and that it was worthwhile trying it out on a more complex example. Furthermore, I thought it might be of interest to record what I've done in a series of blog posts, as a more detailed example than what I had found so far on the web. As I was wondering if an esoteric problem like the four-colouring grid challenge mentioned previously was a good candidate to use as an example, by sheer serendipity, I came accross a link on reddit by a new programmer about his simple Mandelbrot viewer.

Who does not like fractals? ... Furthermore, I have never written a fractal viewer. This seemed like a good time to write one. So, my goal at the end of this series of posts, is to have a "nice" (for some definition of "nice") fractal viewer that is fast enough for explorations of the Mandelbrot set. In addition, in order to make it easy for anyone having Python on their system to follow along and try their own variation, I decided to stick by the following constraints:

With the exception of cython, I will only use modules found in the
standard library. This means using Tkinter for the GUI.

The program should work using Python 2.5+ (including Python 3).

So, without further ado, and based on the example found on the reddit link I mentioned, here's a very basic fractal viewer that can be used as a starting point.

At this point, perhaps a few comments about the program might be useful

I have tried to write the code in the most straightforward and pythonic way, with no thought given to making calculations fast. It should be remembered that this is just a starting point: first we make it work, then, if needed, we make it fast.

The function mandel() is the simplest translation of the Mandelbrot fractal iteration into Python code that I could come up with. The fact that Python has a built-in complex type makes it very easy to implement the standard Mandelbrot set algorithm.

I have taken the maximum number of iterations inside mandel() to be 20, the same value used in the post I mentioned before. According to the very simple method used to time the application, it takes about 2 seconds on my computer to draw a simple picture. This is annoying slow. Furthermore, by looking at the resulting picture, and trying out with different number of iterations in mandel(), it is clear that 20 iterations is not sufficient to adaquately represent the Mandelbrot set; this is especially noticeable when exploring smaller regions of the complex plane. A more realistic value might be to take 100 if not 1000 iterations which takes too long to be practical.

Tkinter's canvas does not have a method to set the colour of individual pixels. We can simulate such a method by drawing a line (for which there is a primitive method) of length 1.

The screen vertical coordinates ("y") increase in values from the top towards the bottom, in opposite direction to the usual vertical coordinates in the complex plane. While the picture produced is vertically symmetric about the x-axis, I nonetheless wrote the code so that the inversion of direction was properly handled.

This basic application is not really useful as a tool for exploring the Mandelbrot set, as the region of the complex plane it displays is fixed. However, it is useful to start with something simple like this as a first prototype. Once we know it is working we can move on to a better second version. So, let's write a fancier fractal viewer following the outline below:

class Viewer(object):
'''Base class viewer to display fractals'''
# The viewer should be able to enlarge ("zoom in") various regions
# of the complex plane. I will implement this
# using keyboard shortcuts.
#
self.parent.bind("+", self.zoom_in)
self.parent.bind("-", self.zoom_out)
def zoom_in(self, event):
def zoom_out(self, event):
def change_scale(self, scale):
# Since one might want to "zoom in" quickly in some regions,
# and then be able to do finer scale adjustments,
# I will use keyboard shortcuts to enable switching back
# and forth between two possible zooming mode.
# A better application might give the user more control
# over the zooming scale.
self.parent.bind("n", self.normal_zoom)
self.parent.bind("b", self.bigger_zoom)
def normal_zoom(self, event, scale=1):
def bigger_zoom(self, event):
# Set the maximum number of iterations via a keyboard-triggered event
self.parent.bind("i", self.set_max_iter)
def set_max_iter(self, event):
# Like what is done with google maps and other
# similar applications, we should be able to move the image
# to look at various regions of interest in the complex plane.
# I will implement this using mouse controls.
self.parent.bind("<button-1>", self.mouse_down)
self.parent.bind("<button1-motion>", self.mouse_motion)
self.parent.bind("<button1-buttonrelease>", self.mouse_up)
def mouse_down(self, event):
def mouse_motion(self, event):
def mouse_up(self, event):
# Presuming that "nice pictures" will be eventually produced,
# and that it might be desired to reproduce them,
# I will include some information about the region of the
# complex plane currently displayed.
def info(self):
'''information about fractal location'''

Furthermore, while I plan to use proper profiling tools, I will nonetheless display some basic timing information as part of the GUI as a quick evaluation
of the speed of the application.

Finally, since I expect that both the function mandel() and the drawing method draw_fractal to be the speed-limiting factors, I will leave them out of the fractal viewer and work on them separately. If it turns out that the profiling information obtained indicates otherwise, I will revisit this hypothesis.

Here is a second prototype for my fractal viewer, having the features described above.

# mandel1a.py
def mandel(c, max_iterations=20):
'''determines if a point is in the Mandelbrot set based on deciding if,
after a maximum allowed number of iterations, the absolute value of
the resulting number is greater or equal to 2.'''
z = 0
for iter in range(0, max_iterations):
z = z**2 + c
if abs(z) >= 2:
return False
return abs(z) < 2

And, finally, I implement the missing functions for the viewer in a new main application.