The Julia Set

The Mandelbrot Set’s Less Famous Sibling

Y. Natsume
Cantor’s Paradise

--

Fractals are one of my fondest topics in mathematics. Out of the many intersections between art and science, fractals to me are perhaps the most elegant. I have always been fascinated by fractals, and the first piece of code I ever wrote was to create the Sierpinski triangle fractal.

The Julia set with c = 0.285+0.01i.

The Julia set is one of the most beautiful fractals to have been discovered. Named after the esteemed French mathematician Gaston Julia, the Julia set expresses a complex and rich set of dynamics, giving rise to a wide range of breathtaking visuals. Also as we shall see later on, the Julia set is related to the much more famous Mandelbrot set!

In this article, we will explore how to use Python to create the Julia set!

Complex Numbers

The Julia set lies on the complex plane, which is a space populated by complex numbers. Each point you see in the image above corresponds to a specific complex number on the complex plane with value: z = x + yi,
where i = √-1 is the imaginary unit.

Complex numbers have very interesting dynamics which is what gives the Julia set its complex dynamics, but for now let us take a look at real numbers, the numbers that we encounter in everyday life.

If you take the square of the real number 10, you get 100. Squaring 100 gives 10,000, and squaring 10,000 we get 100,000,000. If we continue squaring the result, we will get an even larger number in the next result, and an even larger number in the next. Doing this enough times, we eventually tend to infinity, and we can say that the result is not bounded.

On the other hand, if you take the square of the real number 0.1, you get 0.01. Taking the square of 0.01 leads to 0.0001, and as with above if we keep taking the square of the result, we get 0.00000001 in the next number, a microscopically small number in the next and an even smaller number in the next. Doing this enough times we eventually tend to zero, and we can say that the result is bounded.

It turns out that complex numbers also behave similarly. If we take the square of the complex number 1 + i, we get 2i. If we square that we get -4. Repeating this we end up with 16, 256, 65536, 4294967296 and so on until we tend to infinity. This is similar to what happens if we iteratively square 10.

Likewise, if we take the square of 0.5 + 0.5i, we get 0.5i, -0.25, 0.0625, 0.00390625 and so on until we tend to 0. This is similar to what happens if we iteratively square 0.1.

Iterative Operations on the Complex Plane

Now instead of doing this iterative squaring operation on a single complex number, we can do this for an grid of complex numbers z lying on the complex plane. The complex plane is a two dimensional plane, with the horizontal axis corresponding to real numbers and the vertical axis corresponding to imaginary numbers. Any complex number can be plotted on the complex plane as shown in the figure below.

The complex plane.

In addition to the squaring operation on the grid of complex numbers z, we also add an arbitrary complex number c after the squaring operation to get the operation z := z² + c.

After multiple iterations of z := z² + c, some complex numbers in z will blow up to infinity, while some will remain bounded. During the iterative process we track how many iterations it takes for a complex number in z to blow up to infinity. It is this number of iterations measurement for each complex number in the grid that contains the image of the fractal in the image above!

Creating the Julia Set with Python

The Python code to create the Julia set is given below. The function julia_set takes the arguments c an arbitrary complex number, num_iter the number of iterations to perform, N the number of grid points on each axis in the grid, and X0 the limits of the grid.

import numpy as np
import matplotlib.pyplot as plt
def julia_set(c = -0.835 - 0.2321 * 1j, num_iter = 50,
N = 1000, X0 = np.array([-2, 2, -2, 2])):
"""
This function creates the Julia set.
Inputs
c: (complex) arbitrary complex number to add to the grid z.
num_iter: (int) number of iterations to perform.
N: (int) number of grid points on each axis of z.
X0: (array) limits of the grid z.
Outputs
x: (array) values of the real x-axis used in the grid.
y: (array) values of the imaginary y-axis used in the grid.
F: (array) the complex grid containing the Julia set.
"""
# Limits of the complex grid.
x0 = X0[0]
x1 = X0[1]
y0 = X0[2]
y1 = X0[3]
# Set up the complex grid. Each element in the grid
# is a complex number x + yi.
x, y = np.meshgrid(np.linspace(x0, x1, N),
np.linspace(y0, y1, N) * 1j)
z = x + y
# F keeps track of which grid points are bounded
# even after many iterations of z := z**2 + c.
F = np.zeros([N, N])
# Iterate through the operation z := z**2 + c.
for j in range(num_iter):
z = z ** 2 + c
index = np.abs(z) < np.inf
F[index] = F[index] + 1
return np.linspace(x0, x1, N), np.linspace(y0, y1, N), F

During each step in the for loop, after z = z ** 2 + c we check which points in z are still smaller than np.inf. The number of iterations taken by each point to blow up to infinity is recorded in F. By plotting F as a 2 dimensional image, the Julia set finally reveals itself. For example, the code below results in the image shown at the top of this post!

The light areas in the image correspond to points in the complex plane z that blow up to infinity very quickly after only several iterations, while the dark areas correspond to points that are bounded even after many iterations!

x, y, F = julia_set(c = 0.285 + 0.01 * 1j, num_iter = 200, 
N = 1000, X0 = np.array([-1.5, 1.5, -1.5, 1.5]))
plt.figure(figsize = (10, 10))
plt.pcolormesh(x, y, F, cmap = "binary")
plt.axis('equal')
plt.axis('off')
plt.show()

As mentioned earlier, the Julia set has a rich and complex set of dynamics. This can be explored by changing the value of c! For example if we use the value c = -0.835–0.2321i, we get the visualization below, which is completely different from the one we saw earlier (which used c = 0.285+0.01i)!

The Julia set with c = -0.835–0.2321i.

The Mandelbrot Set

In fact, we can take things one step further, and define c to be the entire complex grid, by using c = x + y instead using the value of the input argument in the code above. In this case, the resulting fractal is so famous that it has its own name — the Mandelbrot set! I’ll leave it up to you to modify the Python code above to be able to generate the image below!

The Mandelbrot set.

Closing Remarks

This ends our exploration in using Python to create the Julia set, where we showed that a simple iterative squaring operation of complex numbers can lead to complexity and beauty. Can you think of other simple iterative operations which lead to the same result? Thank you for reading!

--

--