# What is the algorithm hiding beneath the complexity in this paper?

So, I am a computer scientist (at least, I'm working to become one..) and I asked a question on here concerning some mathematics behind the Mandelbrot set. A reply I recieved pointed me to this paper. But as I said, I'm a computer scientist. There is a lot of terminology that I don't understand in that paper.

My question is twofold:

• What is the algorithm outlined in this paper?
• What should I do when confronted with these types of answers?

The paper concerns the dynamics that arise when we iterate a complex polynomial that has a parabolic fixed point. The dynamics near a parabolic fixed point are very slow and, as result, the Julia set is hard to compute. The paper shows that, in spite of the challenges, it possible to compute the image of such a set in a reasonable amount of time - specifically, in polynomial time.

A parabolic fixed point is a fixed point $z_0$ (i.e. $f(z_0)=z_0$) such that $|f'(z_0)|=1$ with $f'(z_0)^n=1$ for some integer $n$. Since $$f(z_0+h) \approx f(z_0) + f'(z_0)h = z_0 + f'(z_0)h,$$ a point near $z_0$ doesn't really move much towards or away $z_0$. The point is neither attractive nor repulsive and the dynamics are quite slow near such a point. If, by contrast, $|f'(z_0)| = r < 1$, then $$|f(z_0+h) - f(z_0)| \approx r|h|,$$ so that the point $z_0$ is attractive and points nearby move quickly towards $z_0$ under iteration.

To understand the challenges and the solution, it helps to look at a concrete example, like the Julia set of $f(z)=z+z^5$, which looks like so:

There is a standard approach to draw the filled Julia set of a polynomial. The filled Julia set is, by definition, the set of all complex numbers $z_0$ with the property that iteration from $z_0$ leads to a bounded orbit. The Julia set itself is the boundary of the filled Julia set. The Julia set is of interest because it is invariant under the action of the polynomial and the dynamics restricted to the Julia set are chaotic.

Using this definition, the standard approach to determine if a point $z_0$ is in the filled Julia set of a polynomial is to iterate $f$ from $z_0$ until one of two things happens:

1. The absolute value of the iterate exceeds some value past which we know the orbit diverges to $\infty$ ($\sqrt{2}$ works for this example), or
2. We exceed some pre-specified bailout (maybe 100 or 1000 or 100000).

In the first case, we are certain that $z_0$ lies outside the filled Julia set - we shade the point white. In the second case, we strongly suspect that $z_0$ is in the filled Julia set, but we're not certain as we might just not have iterated enough. We shade the point black.

If we perform this procedure for a dense grid of points in a rectangular domain of the complex plane, then we generate an approximate image of the filled Julia set - the higher the bailout and the more dense the set of points, the better the approximation. This simple strategy is very easy to implement. Here's a $500\times500$ image using a bailout of 2000 on the iteration count:

There is a standard way to detect the boundary between the two regions to yield the actual Julia set. Unfortunately, it is not particularly impressive at this resolution:

So, how do we improve the situation? One approach might be to just increase the iteration count. Unfortunately, this just doesn't work. The reason is that $z_0=0$ is a parabolic fixed point and the dynamics of $f$ near there are very slow. How many iterates, for example, must we apply to move $z_0=1/32$ to $1/16$? The answer is nearly $2000000$! And this is just one point that's just a little bit close to zero. We need points much closer to zero to get a reasonable resolution.

The solution proposed by Braverman has two main portions

1. Use a more clever classification than just "escapes or not". The classification is based on the Leau-Fatou flower theorem.
2. Expand $f^{2^n}(z)$ to obtain polynomials representing high order iterates. These polynomials are used only near the fixed point and only a relatively small number of terms of the polynomial are actually required.

We can describe the more clever classification scheme in the context of this example easily enough. Suppose we iterate $f(z)=z+z^5$ only 100 times. It turns out that most iterates have either exceeded 2 in absolute value, or they have landed in one of the four prominent leafs attached to the origin. We shade points that escaped white and the remaining points according which of the quadrants the point has landed in. The result looks like so:

If we then apply the boundary scanning technique to this, we get the first image above.

It should be pointed out that many of the ideas in Braverman's paper have been in practical use for some time. The point of the paper is not mainly to present a good computational technique but to prove that the algorithm can run in polynomial time. As such, it's a much harder read than it needs to be if your primary interest is simply how to generate Julia sets of such functions.

• Thanks for a great answer. :) I have two main questions. First, what do you do with $f^{2^n}(z)$? Also, what is the technique used to get the boundary? Feb 4 '15 at 0:40
• @StevenFontaine No problem, I'm happy to do it. I happened to read that paper a couple of years ago. The boundary detection technique is pretty standard in image processing. Photoshop can do it, as can any good image processing library. I specifically used convolution with the third kernel described on Wikipedia's image processing with a kernel page. Feb 4 '15 at 2:34
• The $f^{2^n}$ thing is quite a bit more involved or I would have put it in the answer in the first place. The rough idea is that you have a polynomial that represents $f^{2^n}$ for surprisingly large $n$. Imagine, for example, being able to compute $f^{1024}(z_k)$ with one polynomial application. Then you can make much faster progress with your iterations. Importantly, you really only need a few terms of the polynomials as $0.1^{20}$, for example, is already smaller than machine precision. Feb 4 '15 at 2:38
• @Mark McClure : How do you expand $f^{2^n}(z)$ ? Can ypu show some programs/code that use such technique ?
• @Adam Yes, I could get some Mathematica code to illustrate that part as well. It's not super easy and I've got a busy weekend / life right at the moment. One quick comment - when I say "expand", I mean that we should expand $f^{2^n}$ in a Taylor series about the fixed point. I certainly don't mean that we should expand the polynomial in full. We really only need a few terms since, if $z$ is close to the fixed point $z_0$, then $(z-z_0)^k$ will be quite small. Mathematica's Series command works quite well for this. I've never been able to get Maxima's taylor command to work, though. Nov 17 '17 at 18:26