hi fpsunflower,

return 0.5f * R * (float) Math.log(R) / dR;

Is the leading term 0.5 for all powers? Is the log supposed to be the natural log or the log in some other base? Does anyone have the reference to the explanation of how this formula actually works?

The log is the natural logarithm (base e¡2.71...), yes.

The formula comes, more or less, from this sequence of steps: (assume we have the regular quadratic case)

1. We first construct a uniformization function ph(c) (also known as Bottcher coordinates) with will map the exterior od the M set to the unit circle. It's like a deformation function, a smooth deformation in fact (with no wholes or spikes, or as mathematicians call it, a holomorphica map). Since z->z^2 + c and as z grows c has less importance, in the limit (far from the origin) z->z^2. So the uniformization map is:

phi(c) = Zn ^ (1/2^n)

where n is the number of iterations (ideally infinite). Far from the origin the map is almost identity. You can see a picture of phi(c) in the first image here:

http://iquilezles.org/trastero/fieldlines/2. The Hubbard-Douady potential is the next idea. It's defined as G(c) = log(phi(c)) with natural logarithm. This function is useful for demonstrating the connecteness of the M set (and is related to the theory of external rays).

G(c) = log |phi(c)| = log( Zn ^ (1/2^n) ).

You can think of G as a magneitc filed emanating from the M set, but increasing as you move away from it. It is zero inside M. As we said pchi(c) = c when far from M, so G = log|c| as you move away from M. You can see how it looks here:

http://iquilezles.org/trastero/potential.

G = log|Zn^(1/2^n))| = log |Zn| / 2^n

As with any scalar field, you can compute the gradient of it, or in other words, the maximun slope of the field, if it was interpreted as a terrain. The gradient is related to the derivatives, so you must know how to take derivatives of logarithms. Basically,

grad(G) = |dZn| / (|Zn| * 2^n)

Now, you can approximate a smooth function with a linear version of it (a first order Taylor series) with the gradient (first ortder derivative), like you can say that G(C+eps) = G(C) + eps·|grad(G)|. If we take C inside the M set and we want to compute the distance to it in a close point (C+eps), we can isolate eps int he previous formula and get

eps = distance to M = (G(C+eps)-G(C) )/grad(G(C)). But G(C) is zero because it is inside M, so distance = G/drad(G). Replacing now G and grad G with what we got before, and canceling out the 2^n terms, we get:

dist = |Z|·log |Z| / |dZ|

Now, I'm missing the 2 factor, and I'm not sure why. In fact, apparently the real formula is dist = 2·|Z|·log |Z| / |dZ| and then by something called the "Koebe 1/4 theorem" we know the distance must be not more than a quarted of the previous quantity, so:

dist >= 0.5 · |Z|·log |Z| / |dZ|

Apparently the 1/4 factor holds for 3D too. The 2 that I'm missing in the numerator factor Im not sure where it comes from, probably from one of those 2^n thingies, so I would say it is related to the power of the M set indeed. In that case, the distance estimator would be

DE = (power/4) · |Z|·log |Z| / |dZ|

but Im not sure. Perhapse somebody here can confirm?