Select Page

Both the Julia set and the Mandelbrot set are among the most famous fractals sets and it is likely,that you have seen some images of those sets already. They are basically endlessly – repeating patterns in complex plane. To see whether a point from the complex plane lies in a fractal set or not, one has to compute an equation over and over again and see whether it converge or diverge. If for the given point the equation does converge, the point lies in the set. If it diverges, then it doesn’t. Both sets use the same formulas:

• Julia Set Equation: z = z2 +c
• Mandelbrot Set Equation: z = z2 + c

The main difference is however that for a Julia set, the c – constant is a scalar value while for the mandelbrot set, it changes based on the area of interest ( Usually X → [-2,2] and Y → [-2 , 2] ) and therefore is generally a matrix (Complex Matrix). Amount of iterations of the above equations should be “high enough” to properly decode the divergence / convergence criterion. I did not estimated the right amount of iterations, although ~100 seemed just about legal to show the sets properly. The value might be however dependent upon the complex plane discretization and/or floating point precision. The second important parameter is the convergence/divergence threshold. IE, if the function value is 100, is the function considered to be diverging or not? The general rule is that some points close to set’s borders will converge/diverge slower that other points, which are far away from the set. As such, this parameter will mostly affect the set’s borders.

The reason why I have created this post is however that computing these sets is extremely costly operation. Lets say we want to visualize the Julia set in FullHD resolution (1920 x 1080 pixels). We will choose 100 iterations. Note that all computations are complex. As such, (A + Bi) x (C + Di) = (AC – BD) + (CB + AD)i. IE 4 floating point multiply operations and 2 floating point additions, or more precisely 4 additions (by incorporation the complex constant). As a result, the computation of the set will require : 1920 x 1080 x 100 x ( 4 + 4) floating point operations. Whats more, this task is perfectly suited for acceleration by GPU architecture:

``````/****************************************************************
* @brief Main Julia Kernel
* @param Plane: Julia Function Values
* @param C: Custom Julia Set Complex Constant
* @param Iterations: Amount of iterations processed for each point
* @param Threshold: Threshold for computing
* @param NumIter: Amount of Iterations per kernel call
* @param Nx X Dimensioon of the Grid
* @param Ny Y Dimension of the Grid
*****************************************************************/
__global__ void JuliaKernel(cufftComplex* Plane, cufftComplex C, int* Iterations, float Threshold, int NumIter, int Nx, int Ny) {
int tix = threadIdx.x + blockIdx.x * blockDim.x;
int tiy = threadIdx.y + blockIdx.y * blockDim.y;
if (tix < Nx && tiy < Ny) {
int tid = tix + tiy * Nx;

cufftComplex LA = Plane[tid];
int IterProcessed = Iterations[tid];

// Chunk Iteration Processing Loop
for (unsigned int i = 0; i < NumIter; i++) {
cufftComplex Res = LA * LA + C;
cufftReal Absolute = abs(Res);

if (Absolute < Threshold) {
LA = Res;
IterProcessed++;
}
}