10/16/22
This one is just for fun- a brief look at my first introduction to programming a visualization of a math phenomenon, and problems encountered along the way.
Calling the Mandelbrot set beautiful or mesmerizing is as cliche as it gets in math. My copy of James Gleick's "Chaos", which is older than me, has a full color section dedicated to images of the Mandelbrot set more intricate than you'll see here. Nonetheless, when I first read about the mysterious images generated by a few simple rules, I set out to see it for myself.
At the start of this little project a couple years ago, my programming experience was mostly relegated to taking given datasets and breaking them down into more useful information- MATLAB to process experimental data and rediscover physics equations, and Python to summarize statistics of financial time series. I hadn't yet generated anything from scratch.
The setup for generating the set is really simple. Take a complex number, throw into this equation:
...and iterate a couple times. If the tester complex number you start with blows up to infinity after iterating, then that number is not in the Mandelbrot set. If the tester number stays bound after iterating, then that number is in the set. C is your test number and z starts at zero, so for example, starting with a test number of 2:
Obviously, iterating the equation with a test value of 2 results in the z values running away toward infinity. However, let's start with a test value of -1:
Iterating the equation starting with -1 leaves the z values range bound, no matter how many times you iterate. So, -1 is in the Mandelbrot set. Now we test all of the values on the complex plane and see what happens to the z values. Here was the first image my code produced, after about 30 minutes of runtime:
I was ecstatic. As crude as it looks, it was unmistakably an image of the Mandelbrot set. My reference images, again from the "Chaos" book, accompanied by the text, "The Mandelbrot Set emerges. In Benoit Mandelbrot's first crude computer printouts, a rough structure appeared, gaining more detail as the quality of the computation improved," looked like this:
Zooming in and changing the range of test values:
I attempt to "gain more detail" by increasing the number and granularity of the test values:
At this point I ran into the issue of time. The two images above took around 12 hours each to generate. My computer was left running in the evening and I'd wake up to a new image. If the point I tried to zoom into was in an area of the complex plane devoid of any points in the Mandelbrot set, I'd wake up to a totally blank image... not a great use of time. I set out to try and optimize the code so I could explore faster.
I separated out each individual calculation in the code and timed them. I tried using shortcuts based on guilty knowledge- feeding the algorithm a list of points I already knew were in the Mandelbrot set based on previous runs, so the computer wouldn't have to iterate over them again. Despite decreasing the number of calculations that the code would theoretically have to do, the actual execution time wasn't budging much. I completely missed the actual issue - something I hadn't even timed at the very bottom of the code, because I was sure the inefficiencies lied in the math.
The images above, behind the scenes, are scatter plots. If a test value passes the Mandelbrot test, it is plotted on the complex plane as a dot. Upon searching the internet for solutions to speed up the code, I came across a post of someone else complaining about how their large dataset scatter plots were extremely slow. The response on the post was more or less, "Why the hell are you using a scatter plot to visualize your datasets? We know that Python scatter plots are extremely slow and inappropriate for the amount of data you want to show." Whoops...
The final product I'm after is an image file. Why not just set up the code to output images, rather than scatter plots to save as images? Mapping pixels to a grid of complex numbers to test takes 1 line of code. The PIL Python library can create images from arrays. With some adjustment, the code can now produce images of the Mandelbrot set in seconds or minutes at (custom) high resolution, where before it took hours. Taking a step back to realize that my efforts to optimize where misplaced was the difference it took. Here are some of my favorites:
You can run the code for yourself in a new Jupyter Notebooks file. Enjoy:
Parameters:
origin = point on the complex plane where the top left pixel of the image will be located
side = side length, basically a zoom function
resolution = pixel density, final image will be resolution x resolution pixels
threshold = number of iterations to perform on each pixel/test value. If z values are still contained with absolute value less than 2 after threshold iterations are performed, test value is deemed to be included in the Mandelbrot set
Total run time will be proportional to resolution x resolution x threshold
Code:
import numpy as np
from PIL import Image
#SETTINGS
#weird lump
origin = [-0.14525, 0.65135]
side = 0.0006
#top and bottom swirls
#origin = [-0.147, 0.652]
#side = 0.003
#centered
#origin = [-2.2, 1.5]
#side = 3
#seahorse valley
#origin = [-0.7498, 0.1020]
#side = .007
#upper right off main cardiod
#origin = [-0.0, 0.73]
#side = 0.2
resolution = 1000
threshold = 500
#FUNCTIONS
#pixel to coordinate conversion
def coord(a, b):
return [origin[0] + a * side/resolution, origin[1] - b * side/resolution];
#step function
def step(z,c):
return z**2 + c
#testing a single value multiple times
def tester(c):
n = 0 + 0j #initial z value
steps = threshold #max number of times to iterate
i = 0
while i < steps:
if abs(n) > 2: #broken past 2
#print('broke out in '+ str(i) + ' steps')
return False
break
n = step(n,c) #iterate
i += 1
if i == steps:
#add to the mandelbrot set
#print('added to the set')
return True
#EXECUTION
#set up empty array for image
data = np.zeros((resolution, resolution, 3), dtype=np.uint8)
i = 0
j = 0
while j < data.shape[0]:
while i < data.shape[1]:
# print (i, j, coord(i,j))
if tester(complex(coord(j,i)[0], coord(j,i)[1])):
data[i,j] = [207, 181, 59]
i += 1
i = 0
j += 1
Image.fromarray(data)
The current setup will display something like this: