# Cubic Deltas

Robert P. Munafo, 2023 Jun 17.

The cubic delta approximation technique, pioneered by Kevin Martin (also called series approximation of degree 3), is a perturbation method using a single reference point and third-order polynomial (i.e. cubic) extrapolation.

A single reference point x is iterated at full precision, with its
iterates called x_{n}. The variable naming is different, so first we
will restate the Mandelbrot calculation (as seen in the iteration
article):

x_{0} = complex number to be iterated

x_{n+1} = x_{n}^{2}+x_{0} (eq.1)

Any other desired point y_{0} is iterated by first noting its
difference from x_{0}, i.e. ∂=y_{0}-x_{0} or equivalently
y_{0}=x_{0}+∂. To iterate a
different point y_{n} involves calculating the difference of each
y_{n} from x_{n}. These two are equivalent:

∂_{n} = y_{n} - x_{n} (and the same for n+1)

y_{n} = x_{n} + ∂_{n} (and the same for n+1)

As for x, we iterate y by the same formula as (eq.1) above:

y_{n+1} = y_{n}^{2}+y_{0}

Replace the y_{n} with x_{n}+∂_{n}, and replace y_{0} with
x_{0}+∂:

y_{n+1} = (x_{n}+∂_{n})^{2}+x_{0}+∂

Expand the binomial:

y_{n+1} = x_{n}^{2} + 2x_{n}∂_{n} + ∂_{n}^{2} + x_{0} + ∂

Because x_{n+1} is x_{n}^{2}+x_{0}, we can combine those two:

y_{n+1} = x_{n+1} + 2x_{n}∂_{n} + ∂_{n}^{2} + ∂

Move the x_{n+1} to the other side:

y_{n+1} - x_{n+1} = 2x_{n}∂_{n} + ∂_{n}^{2} + ∂

Because ∂_{n+1}=y_{n+1}-x_{n+1}, we now have an expression for
∂_{n+1}:

∂_{n+1} = 2x_{n}∂_{n} + ∂_{n}^{2} + ∂ (eq.2)

### The Cubic Model

To approximate the ∂_{n} values without calculating both x_{n} and
y_{n}, we use a model that extrapolates what will happen as
∂_{n} is iterated, based on the iterations of x_{n}. The
particular model we use is a cubic polynomial with coefficients A,
B, and C. These coefficients also change with each iteration so
they have subscripts:

∂_{n} = A_{n}∂ + B_{n}∂^{2} + C_{n}∂^{3} + o(∂^{4}) (eq.3)

The "o(∂^{4})" part means that there is an additional quantity that is
assumed to be of a magnitude comparable to ∂^{4}.

The initial ∂, which we can also think of as "∂_{0}", has the coefficients:

A_{0} = 1 B_{0} = 0 C_{0} = 0

Now we can derive a way to iterate the values of these coefficients.
Substitute (eq.3) into (eq.2) omitting the o(∂^{4}) part, expand, and
regroup terms according to the power of ∂:

∂_{n+1} = 2x_{n}∂_{n} + ∂_{n}^{2} + ∂

= 2x_{n}(A_{n}∂ + B_{n}∂^{2} + C_{n}∂^{3})
+ (A_{n}∂ + B_{n}∂^{2} + C_{n}∂^{3})^{2} + ∂

= 2x_{n}A_{n}∂ + 2x_{n}B_{n}∂^{2} + 2x_{n}C_{n}∂^{3}
+ (A_{n}∂)^{2} + (B_{n}∂^{2})^{2} + (C_{n}∂^{3})^{2}
+ 2A_{n}∂B_{n}∂^{2} + 2A_{n}∂C_{n}∂^{3}
+ 2B_{n}∂^{2}C_{n}∂^{3}
+ ∂

= (2x_{n}A_{n} + 1)∂
+ (2x_{n}B_{n} + A_{n}^{2})∂^{2}
+ (2x_{n}C_{n} + 2A_{n}B_{n})∂^{3}
+ o(∂^{4}) + o(∂^{5}) + o(∂^{6})

Now we see that:

A_{n+1} = 2x_{n}A_{n} + 1

B_{n+1} = 2x_{n}B_{n} + A_{n}^{2}

C_{n+1} = 2x_{n}C_{n} + 2A_{n}B_{n}

Notice that we need to multiply the high-precision values x_{n} by
the lower-precision values A_{n} etc. but this multiplication can be
done at the lower precision by truncating the iterates x_{n} (which
need only be done once and the values reused for every y_{n} pixel).

The above was developed and implemented in Java by Kevin Martin in the program "SuperFractalThing", announced to fractalforums in Apr 2013.

### Later History

This method was reimplemented by Botond Kosa, later improved by Karl Runmo and implemented in his program Kalles Fraktaler. (That program also applies these methods to multibrot sets and many other recurrence relations.) The cubic model breaks down when iterating points near islands whose period is less than the number of iterations; yielding incorrect dwell values for pixels. This was clearly visible in images, and called "glitches" in the fractal community. For some time, glitches were fixed manually by designating additional reference points in the image to be iterated at full precision.

Eventually it because clear that automatic detection of some kind was
needed; and in 2014 Pauldelbrot posted to fractalforums an automatic
detection method. It was based on second-order perturbation theory
(approximating the iterated ∂_{n} values as ∂ plus ε_{n}), which
he developed in his program "Nanoscope". It effectively
locates embedded Julia sets and related nodes within the image
so that the pixels near them can be iterated separately based on their
own reference point, chosen to be near the related island.

Sometime around this point the bilinear approximation technique was adopted in place of cubic deltas; see that article for more.

From the Mandelbrot Set Glossary and Encyclopedia, by Robert Munafo, (c) 1987-2024.

Mu-ency main page — index — recent changes — DEMZ

This page was written in the "embarrassingly readable" markup language RHTF, and was last updated on 2023 Jul 03. s.27