Mandelbulb math
Not just for Mandelbulbs!
Discoveries and problems
Domain colouring in 3D
Tricomplex / face-centred complex construct
Caveats
Nomenclature
Background
Shape metrics
Metric view of shape inversion
Constructing custom metrics
Doyle spirals
Mapping plane to line and back
Accumulation-iteration scheme
Wrap-around iteration
Practical considerations
Eulerian and Hamiltonian cycles
[2019-06-28] The Mandelbulb is a 3D analogue of regular Mandelbrot sets of polynomials of the form zn + c. This requires a 3D extension of complex numbers, in particular the complex product.
The basic idea of Daniel White and Paul Nylander is ingenious in its simplicity. It begins with the polar coordinate form of complex multiplication:
The natural 3D analogue of polar coordinates are spherical coordinates. They are likewise defined by the distance from the origin, but now there are two angles to specify the direction. The analogous multiplication then becomes(r1, φ1)(r2, φ2) = (r1*r2, φ1 + φ2)
That's all there is to it! Division follows as straightforwardly as with complex numbers, so rational functions are easily defined. The same goes for roots, with similar multiplicity caveats as in the complex plane. Addition and subtraction are defined in the "Cartesian" sense using the raw xyz coordinates.(r1, φ1, θ1)(r2, φ2, θ2) = (r1*r2, φ1 + φ2, θ1 + θ2)
So what's new here? For some reason, everyone just wants Mandelbrot sets, and the new math was only created to 3Dify those, by generalizing the complex power of a single number. White's page even argues how boring other 3D fractals are in comparison to the unending variety found in Mandelbrot sets. There is some merit to this, considering that a single Mandelbrot set "contains" all possible Julia sets of the respective function.
Still, such a narrow focus misses an enormously wider potential of new mathematics and art. It's a bit like inventing fluid mechanics and only using it for bullets. Or using Fourier transforms to solve differential equations and nothing else; as Cambridge physicist Michael Hobson put it, this is akin to using Mona Lisa as a wrapper for fish and chips.
I started to explore Julia bulbs in late 2018, and I soon realized the math would also work with rational functions and roots. I also experimented with potential problems and found the following two issues:
Spherical coordinates in general have a number of caveats, and these explain a part of the problems with Mandelbulb math. For example, crossing a pole makes φ jump by π, so the mapping between Cartesian and spherical coordinates is not continuous. Gimbal lock is another issue around the poles: with θ at 0 or π, φ is undefined, so moving away from the pole in a controlled direction needs some thought.
The periodicity and non-injectivity of trigonometric functions also accounts for some of the issues. Generally, things are fine and well-defined when working within the spherical coords, but back-and-forth conversions cause issues. The Cartesian→spherical mapping is not always unique.
For pure mathematicians, Mandelbulb math might not be very fruitful. But for mathematical artists, it can be a wonderful tool. There are of course other ways to define higher-dimensional algebras; one I have tested in passing is Bicomplex numbers.
Domain colouring is a relatively traditional technique for visualizing complex numbers, for example using
hue = φ/2π saturation = r value = 1which corresponds to the complex plane mapped onto the colour wheel. There is an infinite variety of other mappings, so going 3D is not going to make it any simpler. For example, the extra angle θ can be used to vary any of the HSV colour coordinates.
Simpler alternatives can take the RGB components straight from the xyz coordinates, for example
(r, g, b) = (1 + x, 1 + y, 1 + z)a variant of which was used to find this two-faced quadratic Mandelbulb.
[2019-07-11] On July 3rd, 2019 I came up with another way to extend complex math into 3D: do a given complex transformation in each of the xy, yz and zx planes, and add the results together. Applying this to the classic Mandelbrot set of z2 + c produced this wonderful rose-like form. Getting such a nice surprise felt like a good sign, so I continue to explore the method further. I named the form Tricomplex Rose, partly influenced by the prominent 3-fold symmetry, thus naming the entire approach.
Given f: R2→R2 and z ∈ R3, w = tricomplex(f)(z) can be compactly defined by using GLSL, as it allows simple access to each of the xy, yz and zx planes:
vec3 w = vec3(0.0);
w.xy += f(z.xy);
w.yz += f(z.yz);
w.zx += f(z.zx);
There is nothing particularly complex-analytic here, any 2D function
can be used. For example,
this 2D IFS
was turned
into 3D using
the Tricomplex approach. While it preserved some of the main
features, it gained some interesting deformations as it went from a
perfect circle into a not-so-perfect sphere.
In general, triplicating a transformation this way does not preserve its algebraic properties. For example, inverting transformations becomes harder, as the 3 intermediate results are folded into each other.
Unlike Mandelbulb math, the Tricomplex approach is not an algebraic system. It is merely a way to construct 3D transformations out of 2D ones, intended for artistic rather than scientific purposes.
All 2D functions can be triplicated, but some are more tricomplex than others. For example, the complex cube does not make well-contained Julia sets this way. It has more symmetry than complex square to begin with, and more terms cancel each other out in the sum.
Using "Tricomplex" for this scheme is somewhat a misnomer. Bicomplex numbers are a full 4D system formed of two 2D complex planes, and tricomplex numbers are the next step in the same series of Multicomplex numbers. So my 3D scheme pales in comparison to actual triply complex systems. It is also not strictly connected to complex-analytic (holomorphic) transformations. But it is not really a number system, either, as noted above.
Nevertheless, this is a 3D system particularly useful with complex transformations. Many 2D functions can be easily generalized to any dimensionality, but complex numbers are very much a 2D-only system. So any way of extending some complex ideas into 3D is welcome.
I might also call this "face-centred complex" in analogy with cubic crystal systems. This gives a more concrete visualization, by considering the faces of the cube as complex planes. The fcc lattice also has the same 3-fold symmetry about the body diagonals. The only downside with this name is its length.
This was almost an accidental discovery, a by-product of sorts. I was looking for a 3D fractal with a simple Cartesian definition, for easier calculation of surface normals via differentiation. The results, my first fractal art with exact reflective lighting, are the Julia sets of the Tricomplex Rose.
[2020-06-29] Shape inversion is a generalization of circle/sphere inversion. I started to play with it in late 2018, having learned about it in Gregg Helt's Bridges 2018 paper. As soon as I started coding my own hexagonal inversion test, I stumbled upon a new, metric approach to shape inversions. This turned out computationally efficient and conceptually satisfying compared to the regular shape inversion, while producing exactly the same results. I have used it in a multitude of OpenGL demos, and by the latter half of 2019 I had enough trust in it (and myself) to write a paper on it.
The metric approach is particularly powerful with the optimization presented in Kento Nakamura's Bridges 2016 paper for rendering inversion fractals. In essence, points are only inverted if they are inside the inverting shape. We can use a single value of the shape metric, first to check if we are inside the shape, and then to perform the inversion. For a system of mixed shapes, we only need to change the metric function for each shape, the rest of the logic can stay the same.
See also my Bridges 2020 short film and my Inversion fractals playlist on Youtube. These are mostly about regular sphere inversions, but some shape inversions are also featured.
As a by-product of the metric approach, I also got a method for constructing metrics from shapes. This is mainly useful for applications other than inversion, and it is also outlined in the paper.
| Some of my Voronoi diagram art using nontrivial metrics on social media |
|---|
| ig/BrvKDtcFg6_ = fb/274819503198185 |
| ig/BrfXiJwgHHU |
| ig/BjhNMIlF2dB |
[2020-07-19] A few days ago, I learned about Doyle spirals from the wonderful artworks of Tis Veugen. I recall seeing such circle packings earlier, but this time it caught my attention as I had just been working with spirals and circle inversions. There were some recipes floating around, but they seemed cumbersome for my GPU demos. So I came up with my own approach:
The lattice must be scaled and optionally rotated for a clean wrap-around of 2*pi along the imaginary (or y) axis. This can be done in countably infinite ways, matching a given pair of integer parameters (p, q). Steady movement in a suitable direction will cause the spiral to rotate and/or expand; this can be directed so that one set of spirals is apparently still while zooming in or out.
This approach is a continuous approximation to a discrete problem — a bit like approximating a sum by an integral. In practice, I've found it visually indistinguishable from the real thing. With very small parameters such as (3, 5) the original discs become deformed, but for example (6, 11) looks fine. This is basically because the complex exponential is analytic/holomorphic. Such functions map small discs to small discs; on an infinitesimal scale, they only do shift and uniform scaling with rotation.
Math demos involving infinite lattices are often practical to do with OpenGL fragment shaders, which operate in a reverse fashion. So for my first public version of the 2D spiral, I started with the complex logarithm and continued with the shift/scaling/rotation.
To make an "S" form of two interconnected spirals, add stage 3: a Möbius transformation that maps 0 to a and infinity to b — these points a and b are the "eyes" of the swirls. This stage has no risk of deformation, the Möbius always maintains circles as circles (though some of these may become straight lines of infinite radius). There are many suitable functions for a given pair of points a and b; for my demos I just used the first/simplest that came to mind: f(z) = (a + b*z)/(1 + z). Other choices will simply change the rotational offset about these points.
This 3D dual spiral was done with ray marching, likewise in the reverse process with a fragment shader. This required some 3D extensions to the complex functions; I did not have to aim for a general 3D complex field, which is theoretically impossible anyway, but some aspects were important to bring to 3D. In particular, uniform scaling so that the spheres don't become flat or oblong ellipsoids. As a nice example of this idea, the complex multiplication is equivalent to scaling and rotation. So I did the rotation in the xy plane, and scaled uniformly in all of xyz.
[2020-09-15] Starting with this Collatz conjecture demo, I sometimes need bijective mappings between integer plane coordinates (x, y) and a single integer n. This is a simpler problem than the real-number version of space-filling curves, and there are rather simple solutions.
The Hilbert curve is a known as a space-filling curve, but its discrete approximations lend themselves nicely to the integer problem. Wikipedia provides C code for both directions of the mapping, and in fact I adapted the xy-to-n function from there to GLSL for the above demo. The other direction is used in this Pac-Man demo where I adapted the example code into Python instead.
A rather different application of the Hilbert curve is found in this photo rasterization demo. There we do not use the actual n for each xy, but we need to know which way the line is turning at each point. This is done in GLSL by comparing the n-values of neighbouring cells and looking for a unit change.
In the earlier Collatz conjecture demo, I also used the square spiral, also known as the Ulam spiral. This can be simply described with a diagram:
37 36 35 34 33 32 31
38 17 16 15 14 13 30
39 18 5 4 3 12 29
40 19 6 1 2 11 28
41 20 7 8 9 10 27
42 21 22 23 24 25 26
43 44 45 46 47 48 49
Another of my demos utilizing this is
Tribonacci
Traces. I also like using the square spiral as a more general
shape, such as in
this photo
rasterization demo.
Having used these two in a number of demos, I wanted to try something else, while sticking with the square lattice. For a visualization of the number-theoretical Möbius function, I turned to boustrophedon, as in the following diagram:
1 2 3 4 5 6 7 8
16 15 14 13 12 11 10 9
17 18 19 20 21 22 23 24
32 31 30 29 28 27 26 25
33 34 35 36 37 38 39 40
48 47 46 45 44 43 42 41
49 50 51 52 53 54 55 56
64 63 62 61 60 59 58 57
Now imagine that people
siht ekil etirw ot desu
and I'm not even doing
-tel eht sa ylreporp ti
ters are not mirrored :)
While the see-saw motion might seem messy for written text, it is
great for integer math, as the next n is always one step away
— there is no jump from the end of one line to the beginning
of the next. One can imagine how this would justify its use in
writing as well.
For the Möbius function demo, I also tested the square spiral, but it was not nearly as nice as boustrophedon. In the Ulam spiral, the square-spiral formation reveals diagonal lines of primes, and the same diagonal concentration is even stronger with the Möbius function that deals with the prime factorizations. The Hilbert curve is commonly used to break up such patterns, but in this case it went too far in the jumbled end. So it is good to have alternatives for these things.
[2021-02-13] In August 2020 I got the idea of making music from the Fibonacci word (IG, FB, YT). The first voice would just play the 0s and 1s of the word itself, and I wanted to derive other voices from the same sequence.
Iteration is a typical way to generate structure in math art. A sequence can be regarded as a function, with f(n) being the nth element. But with f(n) ∈ {0, 1}, and in particular f(0) = 0 and f(1) = 1, iterating f(f(f(n))) etc. would not get very far. So I used the sum of the function values as the argument instead. Denoting v(i, t) as the ith voice at time t, this becomes
v(1, t) = f(t)In this picture (IG, FB) the first voice runs along the bottom row, and the higher rows build on it using this scheme. The music version has a delay before adding each new voice, with the values accumulating from that point on.
v(i, t) = f(v(i-1, 0) + v(i-1, 1) + ... + v(i-1, t)) for i > 1
This scheme lends itself to any sequence or function with a limited, discrete range. So far I've followed the Fibonacci word ones with:
Tribonacci word music (IG, FB, YT)The stills are also available at the gallery.
Möbius function music (IG, FB, YT)
Tribonacci and Tetranacci word stills (IG, FB)
Möbius function stills (IG, FB)
Tetranacci word video (IG, FB)
Pentanacci word video (IG, FB)
Möbius function video (IG, FB)
Fibonacci sequence mod 2 video (IG, FB)
Octanacci sequence mod 3 music (YT, IG, FB)
Digits of the Euler-Mascheroni constant base 3 music (YT, FB, IG)
Abacaba sequence music (YT, FB, IG)
As the Möbius function has a balanced range {-1, 0, 1} I use mod 3 to get {0, 1, 2} instead. This makes the sum increasing, but also makes the values readily usable as indices for programming.
Simple, repeating sequences make more regular structures, but they can also be interesting. For example, using 0, 1, 0, 1, 0, 1, 0, 1, ... in this scheme makes the Sierpinski triangle; in this case the period matches the range. If the two don't match, even a very short period can get interesting with this scheme, e.g. the Fibonacci sequence mod 2 = 1, 1, 0, 1, 1, 0, ... in the video linked above.
For other ways to make 2D art out of number sequences, see also Mapping plane to line and back above.
The argument of the first voice function can also be regarded as a sum, namely that of time steps. This gives a more unified feel to the scheme:
v(1, t) = f(1 + 1 + ... + 1)We might then ask where these 1s come from. Why not from the same function as all the other terms of sums? That is indeed a possibility. We need something to run the first round (t = 0) of iterations, but then we can use the last voice value to seed the first voice for t = 1:
v(1, t) = f(1 + f(N, 0) + ... + f(N, t-1))So after bootstrapping the system with a single value, it then wraps around without any further input. I have tested this, but I haven't found the results any more interesting than the regular version. I also find it nice that you can read the original sequence clearly from the first voice/row.
v(i, t) = f(v(i-1, 0) + v(i-1, 1) + ... + v(i-1, t)) for 1 < i ≤ N
In the cases of the *nacci words and Möbius mod 3, sum(f(n)) < n in the long run. This is good news for programming, because we can prepare the list of values in advance, and it need not be longer than the number of time steps. But for the wrap-around case it means that things will get slower and slower, which is another reason to avoid it.
[2024-02-07] I started making demos of graph-theoretical cycles in June 2023, first Eulerian cycles using Hierholzer's algorithm. Soon I got into Hamiltonian cycles too, first using regular backtracking with recursive functions, and later adding Warnsdorf's rule.
In January 2024 I discovered a simple, fast and intuitive method for generating Hamiltonian cycles. This breakthrough took cycle-finding from exponential hardness into low polynomial times: a graph of a few thousand vertices can be solved in a second using Python, although it is not a universal solution for all graphs. I then learned of other, more professional graph theorists working on similarly fast methods, one of which is referenced in my paper.