Smooth coloring of Mandelbrot

I got half nerd-sniped (by my dad, as it quite often happens) after my blog post on my fractal renderer – “but… your continuous coloring, there… how does it work?”. I had implemented it 8 months ago without really asking myself the question, so after answering “eh, some magic”, I finally got curious for real, and I started digging a bit more into it.

Starting from the Wikipedia article footnotes, I ended up on a couple of articles: Renormalizing the Mandelbrot Escape and Smooth Escape Iteration Counts, as well as smooth iteration count for generalized Mandelbrot sets, that mostly allowed me to piece together some understanding of it. I will admit that my own understanding is still a bit messy/sloppy – my aim here is to try to give some understanding somewhere between “apply this function and you get smooth coloring” and “spectral analysis of eigenvalues in Hilbert spaces” (it’s APPARENTLY related; I’m pretty sure I’ve known at SOME point what this group of words meant, but right now it’s barely a fuzzy concept. Math: if you don’t use it you lose it ๐Ÿ˜ฆ ).

Mandelbrot set and discrete coloring

Before talking about continuous coloring, let me set a few concepts and notations so that we’re on the same page.

We iterate, from a point c = x + y\mathrm{i} (where (x,y) are coordinates on the plane), the sequence z = z^2 + c, where z is a starting constant, typically 0, and we look at the points for which this sequence diverges to infinity or does not. The points for which the sequence does not diverge to infinity (i.e., the sequence is bounded) are part of the Mandelbrot set.

Counting to infinity is, however, a bit long. The first thing we can do is – if we know early that the sequence goes to infinity, we can stop counting: the sequence diverges. We happen to know that, for z = z^2 + c, the sequence diverges as soon as its magnitude |z| > R (or equivalently, x^2 + y^2 > R). (Note that the minimum value for this to work is R = 2.) So we can iterate on every point and, as soon as its magnitude reaches 2, conclude that the point is not in the Mandelbrot set. This process gives another information: the iteration at which this condition is met – which is an indication about how fast the sequence diverges.

That only solves a part of the problem – because by definition, the points for which the sequence is bounded will never stop the process in that model. So the second thing we do is to consider that if a sequence has not diverged after a given, fixed number of iterations, then probably it won’t.

Finally, for every possible number of iterations, we define a color (in a way that eventually makes the result pretty/interesting); we also define a color for points whose sequence does not diverge (and that are considered in the set); and we color each point according to that.

Now the thing with this approach is that the number of iterations tend to gather in “bands” – as in the following picture.

Discrete coloring bands

Smooth coloring

The bands come from the fact that two neighboring pixels “tend to” have the same number of iterations, and hence the same color. If the range of possible values is smooth, and if two neighboring points end up not being associated to the exact same value, and if the palette maps them to different colors, then the coloring should be smooth – even with a fairly reduced palette.

The idea is that we want to represent whether a sequence “barely” reached the escape radius or whether it could have done so almost on the previous iteration, so far it is from the threshold. Ideally, we represent that by a function that goes from 0 to 1, we get a nice smooth range, and that yields a nice smooth coloring.

Now many sources go from this statement to “and now, like a rabbit out of the hat, take \log(\log(z)) and we’re done”. Some of them mumble something about “deriving it from Douady-Hubbard potential“, but it’s not really satisfying either. I gave a quick look to that, and I may revisit the topic eventually, but right now the maths go wayyyy over my head ๐Ÿ˜ฆ

Still, the main question I had was “well, if we’re only interested in mapping ‘the space between the bands’, why not just drop a linear operator on that and be done?”. My reasoning was as follows: we know that at iteration i, we’re above the radius of convergence R, and we were not in the previous iteration, that means that we can bound it. If z_{i-1} and z_i are the values of z at escape iteration i, we have, by definition of i, |z_i| > R and |z_{i-1}| \leq R. To bound z_i from the above, we write |z_i| = |z_{i-1}^2 + c| (by definition of z), which is less than |z_{i-1}|^2 + |c| (by triangular inequality), which is in turn less than R^2 + |c|.

Hence, our “escape margin” |z| - R is between 0 and r^2 + |c|, and we want to map that to something between 0 and 1, so let’s just divide those two things, define my fractional iteration count as \displaystyle i + 1 - \frac{|z| - R}{r^2 + |c|} and we’re done, right?

Hmmm. It DOES look smooth-er. It’s not… smooth. I have two guesses here. First: the triangular inequality is a bit too rough for my estimate. Second: it may be that I’m generally speaking missing something somewhere. Maybe a bit of both.

Now that’s where the “magic” happens. I’ve read enough on the topic to have some inkling that there’s actual, proper math behind that makes it work, but I’m not able to understand, let alone explain it in simple terms ๐Ÿ˜‰

The fractional iteration count that makes things work is \displaystyle i + 1 - \frac{\ln \ln |z|}{\ln 2}. The best intuition/simplified explanation I’ve seen for this comes from Renormalizing the Mandelbrot Escape, where they have the following quote:

The following simplified and very insightful explanation is provided by Earl L. Hinrichs:
“Ignoring the +C in the usual formula, the orbit point grows by Z := Z^2. So we have the progression Z, Z^2, Z^4, Z^8, Z^16, etc. Iteration counting amounts to assigning an integer to these values. Ignoring a multiplier, the log of this sequence is: 1, 2, 4, 8, 16. The log-log is 0, 1, 2, 3, 4 which matches the integer value from the iteration counting. So the appropriate generalization of the discrete iteration counting to a continuous function would be the double log.”

Smooth coloring

So, there. I’m going to leave the topic there for now – that was quite fun to dig into, but I do need more advanced math if I want to dig deeper. It might happen – who knows – but not now ๐Ÿ™‚

New year Marzipan

Around a year ago, I started playing with a small fractal generator that I called Marzipan (because, well, Mandelbrot – basically marzipan, right). Okay, in all fairness, I re-used that name from a previous small fractal generator that I had coded in Javascript several years ago ๐Ÿ™‚

It’s primarily (for now) a renderer of colored Mandelbrot sets. The Mandelbrot set is fairly well-known:

Mandelbrot set

This image is a made on a 2D plane whose upper-left coordinate is at (-2, -1) and lower-end coordinate is at (1, 1). For each point (x, y) of this plane, we consider the complex number z = x + yi, we do some operations (called iterations) repeatedly on z, and we see what is the result. If z grows larger and larger (if the sequence of operations diverges), the point z is not in the Mandelbrot set, and we color it grey. If it does not, or if we give up before seeing that it diverges, the point z in in the Mandelbrot set, and we color it black. The longer we wait before deciding that a point is in the set, the more precise the boundary is (because we have “more chances” that it diverges if it is going to diverge).

And it’s technically possible to zoom as much as you want on the border of the set and to get “interesting” results at infinite amount of zoom – every little part of the border has an infinite amount of detail.

Now one of the reason why Mandelbrot (and other similar computation) renderings are quite popular is that they’re colorful and pretty. In particular, I have a fair amount of screen wallpaper coming from Blatte’s Backgrounds, that contain some of this kind of images. The connection between my set above and the colorful images I linked is not a priori obvious.

Enter: coloring algorithms. The idea is to try to represent graphically how the points outside of the set diverge. The first algorithm I implemented was a escape time algorithm: if the point diverges after 5 iterations, color it in a given color, after 6, in an other color, after 50, in yet another color, and so on. And the fastest way to generate a color palette is to just generate it randomly (and to affect one color to each possible number of iterations), which can yield… fairly ugly images ๐Ÿ˜‰

Random coloring of the escape time

A variation of that approach is to affect to each number of iteration a color that is “proportional” (for instance, more or less saturated) to the number of points that actually reach that number of iterations.

Histogram coloring of the escape time

Then the next idea is to go from “discrete coloring” to “continuous coloring” – right now, we have bands of colors that may have some aesthetic quality, but that may not work for what one has in mind in the end. To achieve that, we add a “continuous” component to our iteration computation (so that it’s not integer anymore) and we map it to the color palette.

Continuous coloring of the escape time

The other coloring algorithm that I started exploring today is coloring by orbit traps. Instead of considering when the iteration escape, we look at how it escapes, and we try to represent that. The first idea is to take an arbitrary, “interesting” point (from an aesthetic point of view, and mostly found by trial and error at this stage), and to look at how close the iterations of the escaping points come to this fixed point). The colored values are the distances to this point. (And the image at the beginning of this post is a (low) zoom from this one ๐Ÿ™‚ ) Note: on this one I also tweaked the palette computation to get more contrast. That was fun ๐Ÿ™‚

Coloring of the distance to an orbit point (-0.25 + 0.5i)

Generally speaking, this project is also for me a nice visual sandbox to play around – on top of practicing my C++, my goal is to generate pretty images, but that typically requires a fair amount of “quality of life” updates:

  • a very basic set of command-line options so that I could generate images without hard-coding all the values
  • quicker than I would have thought: a minimal Qt UI that allows me to zoom and increase/decrease the number of iterations – and right now I kind of feel the need to expand that UI so that it’s… actually useable (being able to change parameters on the UI, re-scaling the window to fit the ratio of the rendered input, that sort of things)
  • yesterday, I sped up the rendering by… well, adding threads ๐Ÿ˜› (via the QtConcurrent library).

Generally speaking, it’s a fun project – and it’s actually something I can go back to quite quickly (once I go over the shock of “urgh C++” – I actually DO like C++, but I did AoC in Go, and it’s a fairly different language) and implement a thing or two here or there, which is nice. For instance, a few months ago, I went “given my current code, can I add support for Julia sets within 10 minutes before going to bed?” and the answer was yes:

Julia set for z = -0.4+0.6i

Advent of Code 2018

In the past few years, December has been for me “the month of Advent of Code”. Advent of Code is an advent calendar with puzzles that can mostly be solved by programming: the input is a problem description and a user-specific input (as far as I know, there’s a set of a “few” pre-validated/pre-tested inputs), and you have to compute a number or a short string, and provide that as a result. If you have the correct result, you unlock a star and the second part of the puzzle, which allows to unlock a second star. Over the course of the 25 first days of December, you consequently get to collect 50 stars.

I like the puzzles, I like the fluff around it (typically a more-or-less-sensical story around the Christmas theme), I like the duration and effort involved, I like that I learn something every year – all in all, I’m a fan ๐Ÿ™‚ There’s a “global leaderboard” that tracks all the users who are in the first 100 to solve each puzzle; although I’m on a “favorable” timezone for the puzzles (they come up at 6am for me), I have no hope of ever scoring points on the leaderboard; my claim to fate this year is that I was in the first 1000 players for 11 stars out of 50 ๐Ÿ™‚

This year, as last year, I solved the whole thing in Go; my Go is still quite bad, probably not-a-all Go-like, and I tend to try to avoid spending too much time second-guessing and refactoring my code in that kind of contexts (it’s good training against analysis paralysis ๐Ÿ˜‰ ), and it’s available on my GitHub if you’re not afraid of bad code and spoilers ๐Ÿ™‚ (Oh, and also: truly horrifying code organization.)

After this general introduction, I want to talk a bit about my experience of the individual problems, so beware: spoilers ahead!

Continue reading “Advent of Code 2018”