Visualising the Julia Set

I previously did a project visualising the Mandelbrot set using Python and really enjoyed it so I decided to do a similar project: visualising the Julia set. This time I used C++ and OpenGL so that the program would be faster and the images would be of higher quality.

Definition

The Julia set is very similar to the Mandelbrot set. It is defined as the set of all complex values $z = x + iy$ in the complex plane where the recursion relation $$ z_{n+1} = z_n^2 + c $$ is bounded as $n \to \infty $, and $ z_0 = x + iy $. That is, a value $z$ is in the Julia set if $ | z_n | \leq 2 $ for all $n > 0$.

The difference between the Julia set and the Mandelbrot set lies in the initial condition $z_0$ and the constant $c$. In the Mandelbrot set, we have $z_0 = 0$. For the Julia set, we have $z_0 = x + iy$ (the point we wish to determine is in the set). Furthermore, in the Mandelbrot set, $c$ takes on each coordinate value. However, in the Julia set, $c$ is a constant values.

Code

I recommend checking out the OpenGl tutorial and playing around with some examples before diving into the code for computing the fractals, since a lot of the information relating to OpenGL is quite new. I had never used OpenGL before and I found the tutorial quite helpful.

Since the main.cpp file relates to setting up the GLFW window, almost all of the code from the tutorial examples can be reused.

main.cpp

#include <glad/glad.h>
#include <GLFW/glfw3.h>
#include <glm/glm.hpp>
#include <glm/gtc/matrix_transform.hpp>
#include <glm/gtc/type_ptr.hpp>

#include <iostream>
#include <shader.h>


// Global variables for glfw window object
const int screen_width = 1080; 
const int screen_height = 1080; 


// Global variables for moving around the window 
float center_x {-0.5f}; 
float center_y {-0.25f};
float zoom {1.0f};


// Function to process input: checks if certain keys are pressed
// This function lets us move around on the screen
void processInput(GLFWwindow *window) 
{
    // If we press ESC, close window 
    if (glfwGetKey(window, GLFW_KEY_ESCAPE) == GLFW_PRESS)
    {
        glfwSetWindowShouldClose(window, true);
    }

    // Up
    if (glfwGetKey(window, GLFW_KEY_UP) == GLFW_PRESS)
    {
        center_y = center_y + 0.05f * zoom;

        if (center_y > 1.0f)
        {
            center_y = 1.0f;
        }
    }

    // Down 
    if (glfwGetKey(window, GLFW_KEY_DOWN) == GLFW_PRESS)
    {
        center_y = center_y - 0.05f * zoom;

        if (center_y < -1.0f)
        {
            center_y = -1.0f;
        }
    }

    // Left
    if (glfwGetKey(window, GLFW_KEY_LEFT) == GLFW_PRESS)
    {
        center_x = center_x - 0.05f * zoom;

        if (center_x > 1.0f)
        {
            center_x = 1.0f;
        }
    }

    // Right
    if (glfwGetKey(window, GLFW_KEY_RIGHT) == GLFW_PRESS)
    {
        center_x = center_x + 0.05f * zoom;

        if (center_x < -1.0f)
        {
            center_x = -1.0f;
        }
    }

    // Zoom Out 
    if (glfwGetKey(window, GLFW_KEY_LEFT_SHIFT) == GLFW_PRESS)
    {
        zoom = zoom * 1.02f;

        if (zoom > 1.0f) 
        {
            zoom = 1.0f; 
        }
    }

    // Zoom In
    if (glfwGetKey(window, GLFW_KEY_LEFT_CONTROL) == GLFW_PRESS)
    {
        zoom = zoom * 0.98f;

        if (zoom < 0.00001f) 
        {
            zoom = 0.00001f;
        }
    }
}

// Function for each time the window size changes
void framebuffer_size_callback(GLFWwindow *window, int height, int width) 
{
    // When the window size changes this function gets called back

    // viewport matches the new window dimensions
    glViewport(0,0,width, height);
}


float vertices[] = {
         1.0f,  1.0f, 0.0f, // top right
         1.0f, -1.0f, 0.0f, // bottom right
        -1.0f, -1.0f, 0.0f, // bottom left
        -1.0f,  1.0f, 0.0f  // top left
    };

unsigned int indices[] = {
    0, 1, 3, // first triangle
    1, 2, 3  // 2nd triangle
};


// Main 

int main()
{

    // initialize glfw and configure
    glfwInit(); 
    glfwWindowHint(GLFW_CONTEXT_VERSION_MAJOR, 3);
    glfwWindowHint(GLFW_CONTEXT_VERSION_MINOR, 3);
    glfwWindowHint(GLFW_OPENGL_PROFILE, GLFW_OPENGL_CORE_PROFILE);

    // Using Apple so I need this line
    glfwWindowHint(GLFW_OPENGL_FORWARD_COMPAT, GL_TRUE);

    // glfw window object
    GLFWwindow* window = glfwCreateWindow(screen_width, screen_height, "LearnOpenGL", NULL, NULL); 
    if (window == NULL)
    {
        std::cout << "Failed to create GLFW window" << std::endl; 
        glfwTerminate(); 
        return -1; 
    }

    glfwMakeContextCurrent(window);
    glfwSetFramebufferSizeCallback(window, framebuffer_size_callback);

    // Check GLAD loaded correctly 
    if (!gladLoadGLLoader((GLADloadproc)glfwGetProcAddress))
    {
        std::cout << "Failed to load GLAD" << std::endl; 
        return -1; 
    }



    // Vertex Objects
    // ---------------------------------------------------------------------------
    unsigned int VAO, VBO, EBO; 
    glGenVertexArrays(1, &VAO); // Create Vertex Array Object
    glGenBuffers(1, &VBO); // Create Vertex Buffer Object
    glGenBuffers(1, &EBO); // Create Element Buffer Object

    glBindVertexArray(VAO); // Bind VAO
    
    glBindBuffer(GL_ARRAY_BUFFER, VBO); // Bind the buffer
    glBufferData(GL_ARRAY_BUFFER, sizeof(vertices), vertices, GL_STATIC_DRAW); // Copy vertices into buffer memory

    glBindBuffer(GL_ELEMENT_ARRAY_BUFFER, EBO); // Bind Element Buffer
    glBufferData(GL_ELEMENT_ARRAY_BUFFER, sizeof(indices), indices, GL_STATIC_DRAW); // Copy indices into buffer memory

    glVertexAttribPointer(0,3, GL_FLOAT, GL_FALSE, 3 * sizeof(float), (void*)0);
    glEnableVertexAttribArray(0);


    // Shaders
    // ---------------------------------------------------------------------------
    Shader myshader("../include/shader/shader.vs", "../include/shader/shader.fs");

    glEnable(GL_DEPTH_TEST);

    // loop to keep the window open
    while (!glfwWindowShouldClose(window))
    {
        // input
        processInput(window);

        glClearColor(0.2f, 0.0f, 0.2f, 1.0f);
        glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);

        myshader.use(); 

        // Create uniforms to be used in fragment shader
        myshader.set_uniform_float("zoom", zoom); 
        myshader.set_uniform_float("center_x", center_x); 
        myshader.set_uniform_float("center_y", center_y); 

        glBindVertexArray(VAO);
        glDrawElements(GL_TRIANGLES, 6, GL_UNSIGNED_INT, 0);
        //glBindVertexArray(0);
        
        glfwSwapBuffers(window); 
        glfwPollEvents(); 
    }

    // Delete VAO, VBO, EBO
    glDeleteVertexArrays(1, &VAO); 
    glDeleteBuffers(1, &VBO); 
    glDeleteBuffers(1, &EBO);

    // close window
    glfwTerminate(); 

    return 0; 
}


An important aspect of OpenGL are shaders, which are small programs that rest on the GPU. Shaders handle all of the computations required, and save a lot of computation time since they are GPU programs. Two important shader programs are the vertex and fragment shaders ( in fact they're required by OpenGL if we want to render anything).

The vertex shader handles the processing of individual vertices and is responsible for projecting our coordinates onto a 2-dimensional space and normalizing the coordinates so that they're visible in the OpenGL window screen.

Here's the code for my vertex shader file. It's really simple and straight forward, since the input data requires no processing (it's already in normalized coordinates).

shader.vs

#version 330 core

// pass in positions from VBO
layout (location = 0) in vec3 pos; 

void main()
{
    // use gl_Position built-in variable
    gl_Position = vec4(pos.xyz, 1.0);
}


The only special part of the vertex shader is the gl_Position variable. This is a built-in OpenGL variable, and essentially tells the GPU the location of our vertex.

The fragment shader file is much more complicated. This file is where all of the computations (which determine if a point is in the Jula set) occur, and is executed over all of the fragments/pixels in the space.

The file has two helper functions. The first is calculate_iterations(), which performs the recursive relation outlined at the start of this page on a point. It returns the number of iterations the recursive relation was used before becoming unbounded, or returns the value of the MAX_ITERS macro. If the MAX_ITERS value is returned, this means that the point is in the Julia set.

The second helper function is julia_colours(). This calls calculate_iterations and returns a vector which corresponds to the colouring we will see on the screen. If the point is in the set, assigns it a black colour. If not, it assigns a different colour.

The main() function simply calls the julia_colours() function and assigns it to the frag_color output variable.

shader.fs

#version 330 core

in vec4 gl_FragCoord; 
out vec4 frag_color; 

#define MAX_ITERS 500

// Pass in our previous global variables using OpenGL uniforms
uniform float center_x;
uniform float center_y;
uniform float zoom;


// Helper function to calculate the number of iterations we reach
int calculate_iterations()
{
    // Initialise z_0 as our coordinate
    float zx = ((gl_FragCoord.x / 1080.0f - 0.5f) * zoom + center_x) * 4.0;
    float zy = ((gl_FragCoord.y / 1080.0f - 0.5f) * zoom + center_y) * 4.0;

    // Initialise our constant c
    // We can change these values to get different shapes
    float cx = 0.356;
    float cy = 0.356;

    // Iterate to see how many iterations before we become unbounded
    int iters = 0; 
    while (iters < MAX_ITERS)
    {
        float zx_temp = zx;
        zx = (zx * zx - zy * zy) + cx; 
        zy = (2.0 * zx_temp * zy) + cy;

        // Check if we become unbounded
        if ( zx * zx + zy * zy > 4.0)
        {
            break;
        }
        ++iters; 
    }

    return iters; 
}

// Function to determine the colours of our Julia set
vec4 julia_colours()
{
    // Call our helper function to compute the iterations
    int iters = calculate_iterations(); 

    // If the point is in the Julia set, give it a black colour
    if (iters == MAX_ITERS)
    {
        gl_FragDepth = 0.0f;
        return vec4(0.0f, 0.0f, 0.0f, 1.0f);
    }

    // otherwise we colour it
    return vec4(0.0f, float(iters) / MAX_ITERS, 0.0f, 1.0f);
}


// Main 
void main()
{
    frag_color = julia_colours();
}


It's common to create a Shader class to deal with all writing and compiling of OpenGL shader programs. It makes our life easier, and cleans up the code in the main.cpp file. The class structure is relatively simple. It consists of a constructor which requires the relative paths to the vertex and fragment shaders as input. The constructor calls two private member functions, std::string Shader read_shader_file(const char* shader_path) and void Shader::make_shader(unsigned int program_ID, const char* shader_path, GLenum Gl_SHADER) which read the shader files and compile the OpenGL shader object, respectively. In order to pass our uniforms to the shader files, I also defined member functiosn void Shader::set_uniform_float(const std::string &name, float val) const and void Shader::set_uniform_vec4(const std::string &name, glm::vec4 vec) const. Finally, I defined a void Shader::use() method to call in main.cpp which allows us to make use of the shader files.

shader.h

#ifndef shader_h
#define shader_h


#include <glad/glad.h>
#include <glm/glm.hpp>
#include <string>
#include <iostream>
#include <fstream>
#include <sstream>


class Shader
{
public: 
    // Data Members
    // ---------------------------------------------
    
    // Shader Program ID
    unsigned int ID; 

    // Member Functions
    // ---------------------------------------------

    // Constructor & Destructor
    Shader(const char* VertexPath, const char* FragPath); 
    ~Shader(); 

    // Method to activate / use our shaders
    void use(); 

    // Methods to set our OpenGL uniforms
    void set_uniform_float(const std::string &name, float val) const; 
    void set_uniform_vec4(const std::string &name, glm::vec4 vec) const; 

private: 
    // Methods relating to creating the shader ob
    void make_shader(unsigned int program_id, const char *shader_path, GLenum GL_shader); 
    std::string read_shader_file(const char* shader_path);
};

#endif


All of the source code for the class methods defined in the header files are written in the shader.cpp implementation file, below.

shader.cpp

#include <glad/glad.h>
#include <shader.h>

// Public Member Functions
// ---------------------------------------------

// Constructor
Shader::Shader(const char* VertexPath, const char* FragPath)
{
    ID = glCreateProgram(); 

    // Create our Vertex and Fragment Shader
    make_shader(ID, VertexPath, GL_VERTEX_SHADER);
    make_shader(ID, FragPath, GL_FRAGMENT_SHADER);

    // Link shaders to our Shader Program
    glLinkProgram(ID);

    // Check Program linked correctly
    int success; 
    char infolog[512];

    glGetProgramiv(ID, GL_LINK_STATUS, &success);
    if(!success)
    {
        glGetProgramInfoLog(ID, 512, NULL, infolog);
        std::cout << "Shader Program Failed to Link\n" << infolog << std::endl; 
    }
}


// Destructor
Shader::~Shader()
{
    if (ID != 0)
    {
        glDeleteProgram(ID);
    }

    ID = 0;
}

// Use the shader
void Shader::use()
{
    glUseProgram(ID); 
}

// Set OpenGL uniforms as a float
void Shader::set_uniform_float(const std::string &name, float val) const
{
    glUniform1f(glGetUniformLocation(ID, name.c_str()), val);
}

// Set OpenGL uniform as a vector
void Shader::set_uniform_vec4(const std::string &name, glm::vec4 vec) const
{
    glUniform4f(glGetUniformLocation(ID, name.c_str()), vec.x, vec.y, vec.z, vec.w);
}


// Private Member Functions
// ---------------------------------------------

// Function to read the a shader file 
std::string Shader::read_shader_file(const char* shader_path)
{
    std::ifstream f(shader_path); 

    // Check if the file was correctly opened
    if(! f.is_open()){
        std::cout << "Failed to open file: " << shader_path << std::endl;
    }

    // Retrieve the stream of characters from file associated buffer object
    std::stringstream shader_stream; 
    shader_stream << f.rdbuf();  
    f.close(); 

    // Turn our stringstream object into a string
    std::string shader_string = shader_stream.str().c_str();
    return shader_string;
    
}



void Shader::make_shader(unsigned int program_id, const char* shader_path, GLenum GL_SHADER)
{
    
    // Read our shader file and extract as a string
    std::string shader_string = read_shader_file(shader_path);

    // Pointer to the shader string as a C string
    unsigned int shader; 
    const char* shader_code = shader_string.c_str();
    
    shader = glCreateShader(GL_SHADER);
    glShaderSource(shader, 1, &shader_code, NULL );
    glCompileShader(shader); 


    // Check the shader configured correctly
    int success; 
    char infolog[512];
    glGetShaderiv(shader, GL_COMPILE_STATUS, &success);
    if(!success)
    {
        glGetShaderInfoLog(shader, 512, NULL, infolog); 
        std::cout << "Shader Failed: " << shader_path << std::endl; 
        std::cout << infolog << std::endl;
    }

    // Attach Shader to Shader Program
    glAttachShader(program_id, shader);
}


Pretty straightforward!

Results

I used a couple different values for $c$. First, I used $c = (0.356, i0.356) $ and this is the result I got,

If I use $c = (0.265, i0) $, I get this shape,