Quaternion Julia with Houdini VEX

Check it on ArtStation : Quaternion Julia.

Quaternion Julia is a type of 3D fractals created using 4D vectors called quaternions. These 4D vectors have a complex structure with one real value and three imaginary values. For more info about it, check this nice article from Paul Bourke.

For the project setup, I used the same one as in the Mandelbox with VEX project but with different code in the VolumeWrangler.

The main idea here is to use quaternions to change each voxel's position during the iterations. We will use a for loop to compute the terms of the complex series where is the the quaternion computed at the iteration n, and is a quaternion identifying the Julia set.

The implemetation with VEX is quite simple because the quaternion multiplication is already implemetend in VEX with qmultiply. Note that the first value in a quaternion is the real part and the three others are the imaginary parts (each of them corresponds to an axis in the 3D space) so the the formula will be like this :

vector4 q = set(d, @P.x, @P.y, @P.z);
q = qmultiply(q, q) + qc;

As we are calculating a series for several iterations, we have to stop when the point is outside the Julia set or when we judge that we iterated enough. In this case, the loop will look like this :

for(int i = 0; i < maxiter; i++)
{
    q = qmultiply(q, q) + qc;
    if(length(q) >= L) {
        break;
    }
}

We are working with voxels inside of a volume wrangler node, so we need to decide what happens to the voxel when it goes outside the Julia set. The trick here is to initialize the @densitiy attribute to 1.0 and then change it to 0.0 if it goes outside. The loop will look like this :

@density = 1.0;
for(int i = 0; i < maxiter; i++)
{
    q = qmultiply(q, q) + qc;
    if(length(q) >= L) {
        @density = 0.0;
        break;
    }
}

Note that maxiter, d, L and qc are not defined in the code above, we can define them using channels. The whole code of the volume wrangler will look like this :

vector4 qc = chp("qc");
int maxiter = chi("MaxIter");
float L = chf("L");
float d = chf("d");

vector4 q = set(d, @P.x, @P.y, @P.z);

@density = 1.0;
for(int i = 0; i < maxiter; i++)
{
    q = qmultiply(q, q) + qc;
    if(length(q) >= L) {
        @density = 0.0;
        break;
    }
}

I suggest you to use these values for maxiter, d andL :

maxiter = 12
d = 0.0
L = 4.0

And then you can play around with qc to find the fractal that you want. For instance, the fractal in the left picture at the beginning of the page is identified by this quaternion :

qc = (-0.162,0.163,0.560,-0.599)

ObeidaZakzak/QuaternionJulia.md

Quaternion Julia with Houdini VEX