biComplexMod5
preset values

Hypercomplex fractals

3D/4D fractals explored with stochastic IIM (Inverse Iterative Method)

Abuot Inverse Iteration Method:

Below I try to explain, in a very elementary way, the IIM (Inverse Iteration Method), and how I use it in glChAoS.p: this is important to know how to use the parameters

The classic Julia set is repulsive with respect to the polynomial fc(z) = z^2 + c, which means that the orbits of points near the set move farther away from it. With the more known DE (Distance Estimator) method, we explore, point by point a region of C space, and assign a color to the speed with which the point flee from orbit. Conversely, the Julia set is "attractive" with respect to the inverse function fc^-1 = sqrt(z-c), because the orbits of points far away from the Julia set come closer to it With the IIM, we choose a initial point on the C space, and we iterate until we reach the orbit.
There are 2 things to observe:
1) the fc^-1 = sqrt(z-c) formula produce two results fx^-1 = +sqrt(z-c) and fy^-1 = -sqrt(z-c) (sqrt property)
So if z0 is initial point { fx^-1(z0), fy^-1(z0) } are two points that will be closer to the Julia set.
From these, with next iteration, we have { fx^-1(fx^-1(z0)), fx^-1(fy^-1(z0)), fy^-1(fx^-1(z0)), fy^-1(fy^-1(z0)) }...
... and so on.
After k iteration, we have 2^(k-1) results: and a typical binary tree, of depth = k-1, to explore.
glChAoS.P uses a stochastic approach obtained by choosing one of the two roots, random, at each stage of the recursion up to a determined depth (selectable), and when we reach it, restart the iteration from starting point.
Statistically, sooner or later, we explore all the branches of the binary tree.
Formally, this can be considered a typical example for an iterated function system (IFS).
2) not all Julia sets generated by a single starting point can be visualized by the inverse iteration method, and a single point may not be enough to define the full object
In glChAoS.P, when the depth has been reached, the re-iteration restart from starting point plus a "possible" random displacement, if are defined "min" and "max" values.
For example, if starting point is p(0,0,0) and min = -1.0 and max = 1.0, the new starting point will be randomly generated in a sphere with center p(0,0,0) and radius = 1. In this way there are more possibilities to define the full set, but are also created some "garbage" starting points on the space: you can visually find a balance.

 

 

3D/4D hypercomplex fractals:

The point clouds is generated from sequences of numbers pn{xn,yn,zn} ⇒ pnR3, nN, where n0→∞ denotes the step of the iteration process starting from a initial p0{x0,y0,z0} point.
In the cloud each next point is function of the previous one:

\[ \eqalign { x_{i+1} = \xi(x_i, y_i, z_i) & \\ y_{i+1} = \phi(x_i, y_i, z_i) & \\ z_{i+1} = \psi(x_i, y_i, z_i) & } \qquad \Bigg\{ \eqalign { & x, y, z \in R \\ & [0, i, n_{\rightarrow\infty}[ \text{   } \Rightarrow i,n \in N \\ } \]

In this section we will see also 4D attractors:

\[ \eqalign { x_{i+1} & = \xi(x_i, y_i, z_i, w_i) & \\ y_{i+1} & = \phi(x_i, y_i, z_i, w_i) & \\ z_{i+1} & = \psi(x_i, y_i, z_i, w_i) & \\ w_{i+1} & = \zeta(x_i, y_i, z_i, w_i) & } \qquad \Bigg\{ \eqalign { & x, y, z, w \in R \\ & [0, i, n_{\rightarrow\infty}[ \text{   } \Rightarrow i,n \in N \\ } \]

where wi coordinate concur to attractor evolution, but in the rendering is plotted the point pi{xi,yi,zi,wi=0}


In the computational code:


 

In the ATTRACTORS window of glChAoS.P:


Colors are indicative of point speed: distance between pi and pi+1


 
You can to start wglChAoS.P with a specific attractor directly from  explore  button. Select lowResources for low resources devices (e.g. mobile devices)
 

Resolution: X        render in new window
 
  touchScreen         Full Rendering Mode

 

 

\[ \eqalign { x_{i+1} & = b \centerdot \eqalign { \frac{(x_i-y_i^2-d)\centerdot(x_i^2+y_i^2)-x_i^3z_i^2}{c \centerdot y_i} }\\ y_{i+1} & = b\\ z_{i+1} & = \frac{-z_i}{\sqrt{2(\rho - d \centerdot (x_i^2+y_i^2)/c)}} } \quad \left\{ \eqalign { \quad \rho & = \sqrt{x_i^2+y_i^2+z_i^2} \\ a & = random(\pm) \sqrt{y_i^2+x_i^2 \centerdot ( 2+ z_i^2/(x_i^2+y_i^2)) - 2 \rho x_i } \\ b & = random(\pm) \sqrt{ rho-x_i-a}/2 \\ c & = y_i^4+x_i^2 \centerdot (y_i^2-z_i^2) \\ d & = a \centerdot (x_i \centerdot (\rho+x_i) + y_i^2) } \right. \]

 

explore
        const float xQ = p.x * p.x, yQ = p.y * p.y, zQ = p.z * p.z;
        float r = sqrtf(xQ+yQ+zQ);
        float a = sign1 * sqrtf(yQ + xQ*(2.f+zQ/(xQ+yQ)) - 2.f*p.x*r);
        float b = sign2 * sqrtf(r - p.x - a) * .5f;
        float c = yQ*yQ + xQ * (yQ - zQ);
        float d = a * (p.x * (r+p.x) + yQ);
        pNew = vec3(b * ((p.x*yQ-d) * (xQ+yQ) - p.x*xQ*zQ) / (p.y*c),
                  b, -p.z/sqrtf(2 * (r - d *(xQ+yQ)/c)));
    

\begin{align} x_{i+1} & = \\ y_{i+1} & = \\ z_{i+1} & = \\ \end{align}

 

explore
    float r = length(p);
    int absOrder = abs(degreeN);
    int k1 = (absOrder - (p.z < 0 ? 0 : 1)), k2 = (3*absOrder + (p.z < 0 ? 4 : 2));
    int dk = ((absOrder%2)==0 && (kPhi << 2)>k1 && (kPhi << 2) < k2) ? (sign(p.z) * ( (absOrder % 4) ? 1 : -1)) : 0;
    float theta  = (atan2f(p.y,p.x) + (2 * kTheta + dk) * T_PI) / float(degreeN);
    float phi    = (asinf(p.z/r)    + (2 * kPhi   - dk) * T_PI) / float(degreeN);
    float cosphi = cosf(phi);
    pNew = vec4(powf(r, 1.0f/float(degreeN)) * vec3(cosf(theta)*cosphi,sinf(theta)*cosphi,sinf(phi)), 0.f);

\[ \eqalign { x_{i+1} & = 0.5 \centerdot (c3_r+c4_r)\\ y_{i+1} & = 0.5 \centerdot (c3_i+c4_i)\\ z_{i+1} & = 0.5 \centerdot (c4_i-c3_i)\\ w_{i+1} & = 0.5 \centerdot (c3_r-c4_r)} \quad \Rightarrow \quad \left\{ \eqalign { \quad & c1\{x_i, y_i\} \\ & c2\{z_i, w_i\} \\ & c3 = random(\pm) \sqrt{c1-c2} \\ & c4 = random(\pm) \sqrt{c1+c2} \\ & c1, c2, c3, c4 \in C } \right. \]

 

explore
    int rnd = fastRand32::xorShift();
    float sign1 = (rnd&1) ? 1.f : -1.f, sign2 = (rnd&2) ? 1.f : -1.f;
    vec4 p(pt - ((vec4 &)*kVal.data()+kRnd));

    std::complex z1(p.x, p.y), z2(-p.w, p.z);
    std::complex w1 = sign1 * sqrt(z1 - z2), w2 = sign2 *sqrt(z1 + z2);
    pNew = vec4(w1.real()+w2.real(), w1.imag()+w2.imag(), w2.imag()-w1.imag(), w1.real()-w2.real())*.5f;

\[ \eqalign { x_{i+1} & = \sigma \centerdot a\\ y_{i+1} & = \sigma \centerdot b \centerdot y_i\\ z_{i+1} & = \sigma \centerdot b \centerdot z_i\\ w_{i+1} & = \sigma \centerdot b \centerdot w_i} \quad \Rightarrow \quad \left\{ \eqalign { \quad \sigma & = random (\pm) \\ \rho & = \sqrt{x_i^2+y_i^2+z_i^2+w_i^2} \\ a & = \sqrt{(x_i+\rho)/2} \\ b & = (\rho-x_i) \centerdot a / (y_i^2+z_i^2+w_i^2) } \right. \]

 

explore
    float xQ = p.x * p.x, yQ = p.y * p.y, zQ = p.z * p.z, wQ = p.w * p.w;
    float r = sqrtf(xQ + yQ + zQ + wQ);
    float a = sqrtf((p.x+r)*.5);
    float b = (r-p.x) * a / (yQ + zQ + wQ);
    pNew = sign * vec4(a, b*p.y, b*p.z, b*p.w);

Previous formulas are a stochastic adaptation of Paul Nylander's Wolfram Mathematica code.

 

 

 

Modified biComplex Julia formula

 

The following formulas were born by chance from BiComplex Julia original formula:
I had not implemented yet 4D iteration, so I have set w parameter to mantain w0 start value, with no increment in iterated formula.
What appeared was pleasant, so I included them in the program... rather, I went to look for others:

 

 

 

And finally I added an element to search for them: biComplex Julia Explorer: