thecodingidiot.com

The InfiniteThe Set

The Set

Implement the Mandelbrot escape-time function and connect it to the pixel buffer. By the end of this page the set is visible — black for points in the set, white for points that escape.


The formula

A complex number has two components: a real part and an imaginary part. Written as c=re+imic = \text{re} + \text{im} \cdot i, where re and im are ordinary decimal numbers and i is the imaginary unit. The Mandelbrot iteration applies this formula repeatedly, starting from zero:

z0=0zn+1=zn2+cz_0 = 0 \qquad z_{n+1} = z_n^2 + c

From the formula to code. To implement z2z^2 in C you need to expand it algebraically. i in the formula is not the loop counter — it is the imaginary unit, with one defining property: i2=1i^2 = -1. Squaring a complex number works like any algebraic expansion:

(re+imi)2=re2+2reimi+im2i2(\text{re} + \text{im} \cdot i)^2 = \text{re}^2 + 2 \cdot \text{re} \cdot \text{im} \cdot i + \text{im}^2 \cdot i^2

Substitute i2=1i^2 = -1:

=re2im2+2reimi= \text{re}^2 - \text{im}^2 + 2 \cdot \text{re} \cdot \text{im} \cdot i

The real part of z2z^2 is re2im2\text{re}^2 - \text{im}^2. The imaginary part is 2reim2 \cdot \text{re} \cdot \text{im}. Adding c to each gives the update rule directly — those two lines in the code are the formula.

Why double, not int or long. The complex plane has non-integer coordinates: points like (0.5,0.27)(-0.5,\,0.27) or (0.743,0.1)(-0.743,\,0.1). int and long store only whole numbers and cannot represent anything between zero and one. double is a 64-bit floating-point type — it stores fractional values and distinguishes coordinates to roughly 15 significant digits. float is 32-bit (about 7 digits of precision): accurate near the default view, but at deep zoom levels adjacent pixels differ in the fourteenth decimal place, and float rounds them to the same value, producing visible smearing. double is the standard choice for fractal renderers.

In C, both components update each step using the previous values:

double new_re = re * re - im * im + c_re;
double new_im = 2.0 * re * im + c_im;
re = new_re;
im = new_im;

This fragment declares variables inline for clarity. In the final mandelbrot.c, all variables are declared at the top of the function scope — C99 style throughout the project. The final implementation also drops new_im entirely: im is overwritten directly, and only new_re needs saving (the mutation order section explains why).

The escape condition is z>2|z| > 2. For a complex number, z|z| is its distance from the origin. A complex number re + im·i is a point on a plane — re along the horizontal axis, im along the vertical. Its distance from the origin is the hypotenuse of a right triangle, by the Pythagorean theorem[1]:

           * z = (re, im)
          /|
         / |
    |z| /  | im
       /   |
      /    |
     *-----*
  (0,0)    re

z=re2+im2|z| = \sqrt{\text{re}^2 + \text{im}^2}

Checking z>2|z| > 2 would require calling sqrt() every iteration. That is unnecessary: sqrt is monotonically increasing, so x>2\sqrt{x} > 2 is true exactly when x>4x > 4. Squaring both sides of the condition eliminates the square root entirely:

z>2    re2+im2>4|z| > 2 \implies \text{re}^2 + \text{im}^2 > 4

re * re + im * im > 4.0

mandelbrot.h

#ifndef MANDELBROT_H
#define MANDELBROT_H
 
int mandelbrot(double c_re, double c_im, int max_iter);
int burning_ship(double c_re, double c_im, int max_iter);
 
#endif

No SDL2 include. The escape-time functions are pure computation — no platform dependency. This is what makes the unit tests possible.


mandelbrot.c

#include "mandelbrot.h"
 
int mandelbrot(double c_re, double c_im, int max_iter)
{
    double re;
    double im;
    double new_re;
    int    i;
 
    re = 0.0;   /* z starts at zero: z0 = 0 + 0i */
    im = 0.0;
    i = 0;
    while (i < max_iter) {
        if (re * re + im * im > 4.0)       /* |z|^2 > 4, i.e. |z| > 2 */
            return (i);
        new_re = re * re - im * im + c_re; /* real part of z^2 + c */
        im = 2.0 * re * im + c_im;         /* imag part of z^2 + c — uses old re */
        re = new_re;
        i++;
    }
    return (max_iter); /* never escaped: point is in the set */
}

Parameters. c_re and c_im are the real and imaginary components of the complex number c — the pixel being tested. max_iter is the iteration ceiling: if z has not escaped after this many steps, the point is treated as inside the set.

Initialisation. re and im represent z, which starts at zero (z0=0z_0 = 0). c is fixed for the whole loop — it is the pixel coordinate passed in by the caller. Each iteration computes the next z from the current one.

The escape check. re * re + im * im > 4.0 is the optimised form of z>2|z| > 2 derived above: comparing the squared magnitude against 4 avoids a sqrt call on every iteration. When the condition is met, the function returns i — the number of steps it took to escape. That count becomes the pixel's colour.

The update lines. These two lines are the algebra from the formula section translated directly into C:

re² - im² + c_re   ←   real part of z² + c
2·re·im  + c_im    ←   imaginary part of z² + c

Mutation order. Both update lines use the current values of re and im — the values from this iteration, not the next. The problem is that im is computed second. If re were overwritten first, the im line would silently use the already-updated re, producing a different (wrong) formula. new_re holds the correct new value while im is still computed from the old re; only then is re updated. This single-variable save is the entire solution.

Return value. If the loop completes without the escape condition ever firing, the function returns max_iter. The caller uses i == max_iter to mean "this point is in the set" and colours it black.


First render

Update render_frame in render.c to call mandelbrot and map the result to a black-or-white pixel:

#include "mandelbrot.h"
 
/* inside render_frame, replacing the placeholder: */
i = mandelbrot(re, im, MAX_ITER);
if (i == MAX_ITER)
    colour = 0xFF000000; /* black — in the set */
else
    colour = 0xFFFFFFFF; /* white — escaped */
s->pixels[y * WIDTH + x] = colour;

The view functions (pixel_to_re, pixel_to_im) are still stubs returning 0. The entire screen maps to (0.0, 0.0), which is inside the set, so the window renders solid black. That is expected — the viewport is the next page.

To confirm the function works before the viewport is ready, temporarily hard-code a few test coordinates in render_frame:

/* temporary test: map pixel (x, y) to a fixed complex range */
re = -2.5 + (double)x / WIDTH * 3.5;
im = -1.25 + (double)y / HEIGHT * 2.5;
i = mandelbrot(re, im, MAX_ITER);
make re && ./infinite

The Mandelbrot set appears: a black cardioid with black discs attached, surrounded by white. It is not interactive yet — the viewport functions are stubs — but the iteration is working. Remove the temporary mapping once confirmed; the real viewport implementation is next.


What you just computed

Every pixel just ran up to 100 iterations of the complex squaring loop. At 800 × 600, that is up to 48 million floating-point operations for a single frame. On one CPU core.

That is why the GPU was built.

Zoom is not wired yet. When it is, deep zooms will require higher MAX_ITER values to reveal detail — and you will feel the cost directly.

Footnotes

  1. Pythagorean theorem - Wikipedia