The Julia set

About the code

The code in this article was written using Code::Blocks and SDL 2.
You can read here a guide to install this software.
Although it is based on SDL, I don't use its functions directly. I have written a small library with a few basic
functions to ease the understanding and the portability to another language.
You can read more about this lib here.

The code in this article uses the complex number class we added in the Mandelbrot set.

The Mandelbrot's cousin

The Julia set is a cousin of the Mandelbrot set.
They are so close that generally when mathematicians want to demonstrate a property of the Mandelbrot set, they
first prove it in the Julia set.
Remember the formula we used for the Mandelbrot set. 
		z0 = 0
		zn+1 = zn2 + cpix
Where everything is a complex number and cpix is the complex number associated with the current pixel.
In fact, as z1 is always equal to cpix, because 02 = 0, we could also use this formula which would give
us the same set. 
		z0 = cpix
		zn+1 = zn2 + cpix
The Julia set use nearly the same formula, except that at each iteration instead of adding cpix, we add
a constant complex number we chose at the beginning of the drawing. 
		z0 = cpix
		zn+1 = zn2 + c
So it's easy to modify one of the programs we used for the Mandelbrot set. 
		#include <stdio.h>
		#include <stdlib.h>
		#include <math.h>
		#include "main.h"
		#include "Graphics.h"
		#include "System.h"
		#include "math/Complex.h"

		#define SCREEN_WIDTH    640
		#define SCREEN_HEIGHT   480
		#define MAX_ITERATIONS  256

		void drawJulia(int offsetx, int offsety, Color* palette)
		{
			for (int y = 0; y < SCREEN_HEIGHT; y += 8)
			{
				for (int x = 0; x < SCREEN_WIDTH; x += 8)
				{
					// convert the pixel pos to a complex number
					Complex z;
					z.re = (x + offsetx - SCREEN_WIDTH  / 2.0) * 4.0 / SCREEN_WIDTH;
					z.im = (y + offsety - SCREEN_HEIGHT / 2.0) * 4.0 / SCREEN_WIDTH;

					// julia parameter
					Complex c(-0.70176, -0.3842);

					int iteration = 0;

					while (z.squareMag() <= 16 &&
					       iteration < MAX_ITERATIONS)
					{
						// compute the new value
						z = z * z + c;
						iteration++;
					}

					// draw the pixel
					Color   col;

					if (iteration == MAX_ITERATIONS)
					{
						col = Color(0, 0, 0);
					}
					else
					{
						col = palette[iteration];
					}

					gfx.setPixel(x + offsetx, y + offsety, col);
				}
			}
		}

		int main(int argc, char* argv[])
		{
			// init the window
			gfx.init("Julia", SCREEN_WIDTH, SCREEN_HEIGHT);
			gfx.init2D();
			gfx.clearScreen(Color(0, 0, 0, SDL_ALPHA_OPAQUE));

			static Uint8 pattern[64] =
			{
				 0, 32,  8, 40,  2, 34, 10, 42,
				48, 16, 56, 24, 50, 18, 58, 26,
				12, 44,  4, 36, 14, 46,  6, 38,
				60, 28, 52, 20, 62, 30, 54, 22,
				 3, 35, 11, 43,  1, 33,  9, 41,
				51, 19, 59, 27, 49, 17, 57, 25,
				15, 47,  7, 39, 13, 45,  5, 37,
				63, 31, 55, 23, 61, 29, 53, 21
			};

			static Color    palette[256];

			// set up the palette
			for (int i = 0; i < 256; ++i)
			{
				float angle = 2.0 * M_PI * (float)i / 256.0;
				float offset = 2.0 * M_PI / 3.0;
				angle += 40.0 *  M_PI / 128.0;
				angle *= 1.3;

				palette[i].r = 128 + 128 * sin(angle);
				palette[i].g = 128 + 128 * sin(angle * 1.1 - offset);
				palette[i].b = 128 + 128 * sin(angle * 1.5 + offset);
			}

			for (int i = 0; i < 64; i++)
			{
				int j = 0;
				while (pattern[j] != i)
					j++;

				drawJulia(j % 8, j / 8, palette);

				// display and process the events
				gfx.render();
				sys.wait(50);
				sys.processEvents();
				if (sys.isQuitRequested() == true)
					exit(EXIT_SUCCESS);
			}

			// wait until we quit
			while (sys.isQuitRequested() == false)
			{
				sys.processEvents();
				sys.wait(50);
			}

			gfx.quit();

			return EXIT_SUCCESS;
		}

		Download source code
		Download executable for Windows
This program draws the Julia set corresponding to the number -0.70176 -0.3842 * i


The Julia family

As the Julia set corresponds to a given complex number, there is not only one set, but a whole family of sets.
To explore them, you could change the value of the constant in the previous program an recompile it.
But I wanted to do something easier, so we will use the mouse coordinates to change this parameter.

First we will need to listen to the SDL_MOUSEMOTION event in CSystem::processEvents() 
		void CSystem::processEvents()
		{
			SDL_Event event;

			if (sdl.PollEvent(&event) != 0)
			{
				if (event.type == SDL_QUIT)
				{
					mQuitRequest = true;
				}
				else if (event.type == SDL_MOUSEMOTION)
				{
					mMouseX = event.motion.x;
					mMouseY = event.motion.y;
				}
			}
		}
Then, we will modify our code to set the constant complex number according to the current mouse position, and
re-start the drawing every time the mouse moves. 
		int main(int argc, char* argv[])
		{
			// init the window
			[...]

			// set up the palette
			[...]

			int lastMouseX = 0;
			int lastMouseY = 0;

			while (sys.isQuitRequested() == false)
			{
				lastMouseX = sys.mMouseX;
				lastMouseY = sys.mMouseY;

				// julia parameter
				Complex c;
				c.re = (sys.mMouseX - SCREEN_WIDTH  / 2.0) * 4.0 / SCREEN_WIDTH;
				c.im = (sys.mMouseY - SCREEN_HEIGHT / 2.0) * 4.0 / SCREEN_WIDTH;

				for (int i = 0; i < 64; i++)
				{
					int j = 0;
					while (pattern[j] != i)
						j++;

					drawJulia(j % 8, j / 8, c, palette);

					gfx.render();
					sys.processEvents();

					if (lastMouseX != sys.mMouseX ||
						lastMouseY != sys.mMouseY ||
						sys.isQuitRequested() == true)
					{
						break;
					}
				}

				sys.processEvents();
			}

			gfx.quit();

			return EXIT_SUCCESS;
		}

		Download source code
		Download executable for Windows
Now if you move the mouse you should see the set change in real time and take various shapes.

The link with Mandelbrot

The Julia sets have a deep link with the Mandelbrot set.
The shapes they take when the complex constant varies are linked to the corresponding point in the Mandelbrot.
If you take a point far away from the Mandelbrot, you will see only a few islands of points.
As you move towards the frontier of the set, the shapes will get more complex, the islands will connect together
and spirals will appear.
And you will only see black areas in the Julia set if the corresponding point lies inside the Mandelbrot set.

To explore that visually, in this 3rd program, we will draw the Mandelbrot set on left of the window and the
Julia set on the right.


If you move your mouse over the Mandelbrot, the corresponding Julia set will appear on the right.
		#include <stdio.h>
		#include <stdlib.h>
		#include <math.h>
		#include "main.h"
		#include "Graphics.h"
		#include "System.h"
		#include "math/Complex.h"

		#define SCREEN_WIDTH    1280
		#define SCREEN_HEIGHT   480
		#define MAX_ITERATIONS  256

		// pixel coordinates to complex number
		Complex pixToComp(int x, int y)
		{
			Complex c;
			double w = SCREEN_WIDTH / 2.0;
			double h = SCREEN_HEIGHT;
			c.re = (x - w / 2.0) * 4.0 / w;
			c.im = (y - h / 2.0) * 4.0 / w;

			return c;
		}

		void drawJulia(int offsetx, int offsety, Complex c, Color* palette)
		{
			for (int y = 0; y < SCREEN_HEIGHT; y += 8)
			{
				for (int x = 0; x < SCREEN_WIDTH / 2; x += 8)
				{
					// convert the pixel pos to a complex number
					Complex z = pixToComp(x + offsetx, y + offsety);

					int iteration = 0;

					while (z.squareMag() <= 16 &&
					       iteration < MAX_ITERATIONS)
					{
						// compute the new value
						z = z * z + c;
						iteration++;
					}

					// draw the pixel
					Color   col;

					if (iteration == MAX_ITERATIONS)
					{
						col = Color(0, 0, 0);
					}
					else
					{
						col = palette[iteration];
					}

					gfx.setPixel(SCREEN_WIDTH / 2 + x + offsetx, y + offsety, col);
				}
			}
		}

		void drawMandelbrot(Color* palette)
		{
			for (int y = 0; y < SCREEN_HEIGHT; y++)
			{
				for (int x = 0; x < SCREEN_WIDTH / 2; x++)
				{
					// convert the pixel pos to a complex number
					Complex c = pixToComp(x, y);
					c.re -= 0.7;

					Complex z;
					int iteration = 0;

					while (z.squareMag() <= 16 &&
					       iteration < MAX_ITERATIONS)
					{
						// compute the new value
						z = z * z + c;
						iteration++;
					}

					// draw the pixel
					Color   col;

					if (iteration == MAX_ITERATIONS)
						col = Color(0, 0, 0);
					else
						col = palette[iteration];

					gfx.setPixel(x, y, col);
				}
			}
		}

		int main(int argc, char* argv[])
		{
			// init the window
			gfx.init("Julia 3", SCREEN_WIDTH, SCREEN_HEIGHT);
			gfx.init2D();
			gfx.clearScreen(Color(0, 0, 0, SDL_ALPHA_OPAQUE));

			static Uint8 pattern[64] =
			{
				 0, 32,  8, 40,  2, 34, 10, 42,
				48, 16, 56, 24, 50, 18, 58, 26,
				12, 44,  4, 36, 14, 46,  6, 38,
				60, 28, 52, 20, 62, 30, 54, 22,
				 3, 35, 11, 43,  1, 33,  9, 41,
				51, 19, 59, 27, 49, 17, 57, 25,
				15, 47,  7, 39, 13, 45,  5, 37,
				63, 31, 55, 23, 61, 29, 53, 21
			};

			static Color    palette[256];

			// set up the palette
			for (int i = 0; i < 256; ++i)
			{
				float angle = 2.0 * M_PI * (float)i / 256.0;
				float offset = 2.0 * M_PI / 3.0;
				angle += 40.0 *  M_PI / 128.0;
				angle *= 1.3;

				palette[i].r = 128 + 128 * sin(angle);
				palette[i].g = 128 + 128 * sin(angle * 1.1 - offset);
				palette[i].b = 128 + 128 * sin(angle * 1.5 + offset);
			}

			int lastMouseX = 0;
			int lastMouseY = 0;

			// julia parameter
			Complex c;

			drawMandelbrot(palette);

			while (sys.isQuitRequested() == false)
			{
				lastMouseX = sys.mMouseX;
				lastMouseY = sys.mMouseY;

				if (sys.mMouseX < SCREEN_WIDTH / 2)
				{
					c = pixToComp(sys.mMouseX, sys.mMouseY);
					c.re -= 0.7;
				}

				for (int i = 0; i < 64; i++)
				{
					int j = 0;
					while (pattern[j] != i)
						j++;

					drawJulia(j % 8, j / 8, c, palette);

					gfx.render();
					sys.processEvents();

					if (lastMouseX != sys.mMouseX ||
						lastMouseY != sys.mMouseY ||
						sys.isQuitRequested() == true)
					{
						break;
					}
				}

				sys.processEvents();
			}

			gfx.quit();

			return EXIT_SUCCESS;
		}

		Download source code
		Download executable for Windows

Links

Video of the second program

Science & Vie Micro n°27 page 83

Video: Julia Set in processing