# Mandelbrot II

# python 3 version of code for making your own mandelbrot fractal
#
# blog and book at http://makeyourownmandelbrot.blogspot.co.uk
from pylab import *
# set the location and size of the atlas rectangle
xvalues = linspace(-0.22, -0.21, 1000)
yvalues = linspace(-0.70, -0.69, 1000)

# size of these lists of x and y values
xlen = len(xvalues)
ylen = len(yvalues)
# mandelbrot function, takes the fixed parameter c and the maximum number of iterations maxiter, as inputs
def mandel(c, maxiter):
# starting value of complex z is 0+0i before iterations update it
z = complex(0,0)

# start iterating and stop when it's done maxiter times
for iteration in range(maxiter):

# the main function which generates the output value of z from the input values using the formula (z^2) + c
z = (z*z) + c

# check if the (pythagorean) magnitude of the output complex number z is bigger than 4, and if so stop iterating as we've diverged already
if abs(z) > 4:
break
pass
pass

# return the number of iterations we actually did, not the final value of z, as this tells us how quickly the values diverged past the magnitude threshold of 4
return iteration
# create an array of the right size to represent the atlas, we use the number of items in xvalues and yvalues

atlas = zeros((xlen,ylen))

# go through each point in this atlas array and test to see how many iterations are needed to diverge (or reach the maximum iterations when not diverging)
for ix in range(xlen):
for iy in range(ylen):

# at this point in the array, work out what the actual real and imaginary parts of x are by looking it up in the xvalue and yvalue lists
cx = xvalues[ix]
cy = yvalues[iy]
c = complex(cx, cy)

# now we know what c is for this place in the atlas, apply the mandel() function to return the number of iterations it took to diverge
# we use 40 maximum iterations to stop and accept the function didn't diverge
atlas[ix,iy] = mandel(c,120)

pass
pass
# plot the array atlas as an image, with its values represented as colours, peculiarity of python that we have to transpose the array
imshow(atlas.T, interpolation="nearest", cmap="jet")
gcf()