Derivative Extrapolation

Robert P. Munafo, 2023 Jul 8.

"Derivative extrapolation", also "differential approximation", is a perturbation method using the iterates of a single reference point, along with the iterated derivative, a linear model for extrapolating nearly points, and a "rebasing" method to remove roundoff error when the iteration comes near the origin.

It is effectively a degree 1 approximation of the recurrence relation, rather than degree 3 as in cubic deltas. The method is called "derivative extrapolation" because it extrapolates iterate values based on a derivative of the reference point's iterates; alternatively it is called "differential estimation" because the differential (difference between a point's iterate and a nearby point's iterate) is being estimated (approximated).

Given the actual iterate values of a "reference point" (whose c value must be near the c of the iteration to be approximated), the iterates of any "nearby" point of interest can be approximated by calculating successive values of the difference between its iterates and the reference point's iterates. The definition of "nearby" is expressed in terms of a view area (typically a rectangle) and the use of a period-finding technique to locate the lowest-period nucleus within that view area (even if its mu-atom is much smaller than a pixel). Period-finding methods are described in: the period article, section "Finding the Period of a mu-Atom", and the articles Jordan curve method and Synchronous-Orbit Algorithm both which have illustrations.

In comparison to cubic deltas, the use of fewer terms means that iteration cannot proceed through any step that puts the reference point very close to the origin; however this is not a problem if that condition is checked on every iteration. When needed a rebasing calculation is performed, in which a new high precision key point value is created based on the previous (good) values; then iteration continues as usual.

In turn, using a linear model makes it easier to perform the condition test with accurate results (no "false negatives"). Such a test was introduced by "Zhuoran" on Fractal Forums in 2021 Dec.

The practical benefit is far fewer "glitches", making it no longer necessary to use higher-degree extrapolation models to get good images. For this reason DE has mostly supplanted cubic deltas. It is used in the more recent versions of Kalles Fraktaler and contemporaneous deep zooming programs.

The Model

Given a reference point C that is known to be periodic (a nucleus or Misiurewicz point) and relatively close to the pixels we want to draw, its iterates (preperiod plus limit cycle) are calculated at full precision, then each iterate is rounded to a lower-precision Z_n that can be at lower precision (typically single or double precision).

For every other pixel, define its parameter in terms of C and an offset: c=C+Δc. Its iterates are defined similarly:

z_n=Z_n+Δz_n for all n (1)

We wish to calculate the Δz_n, which can be done at lower precision because c is close to C and the iterates won't diverge much until the point is close to escape.

The initial value for any iteration is z=0, and so z₀ = Z₀+Δz₀ = 0 amd Δz₀ = 0. Iteration proceeds using the standard recurrence relation:

z_n+1 = z_n² + c (2)

Use (1) to substitute z_n on the right, expand and regroup:

z_n+1 = (Z_n+Δz_n)² + C+Δc
= Z_n² + 2 Z_n Δz_n + Δz_n² + C + Δc
= Z_n² + C + 2 Z_n Δz_n + Δz_n² + Δc

Now we use (2) to replace Z_n²+C with Z_n+1:

z_n+1 = Z_n+1 + 2 Z_n Δz_n + Δz_n² + Δc

and use (1) to replace z_n+1 with Z_n+1+Δz_n+1 on the left:

Z_n+1+Δz_n+1 = Z_n+1 + 2 Z_n Δz_n + Δz_n² + Δc

So we can subtract Z_n+1 from both sides to get:

Δz_n+1 = 2 Z_n Δz_n + Δz_n² + Δc (3)

Therefore the code for the main iteration step is simply:

/* Calculate new dz */ dz := 2 * Z[n] * dz + dz * dz + c;

This matches the two published implementations (computer code). Zhuoran's code has:

dz = 2 * dz * Reference[RefIteration] + dz * dz + dc;

Claude Heiland-Allen's code is refactored slightly to save one multiplication (here with my comments added):

Z = Z_ptr[m]; // Z_ptr is the precomputed array of reference iterates z = (2 * Z + z) * z + c; // "z" is their name for Δz

Relation to Derivative Iteration

Now consider the derivative calculation for constant c, as when making an image of a Julia set. The derivative (with respect to z_n of the recurrence relation is:

d/dz_n (z_n²+c) = 2 z_n

d (z_n²+c) = 2 z_n dz_n

d (z_n+1) = 2 z_n dz_n

This is the first part of (3); the rest of (3) is just the right side of (2). Thus we are augmenting the ordinary iteration formula by also performing the derivative iteration, allowing us to calculate estimates of the differences Δz_n between our desired iterates z_n and the reference iterates Z_n.

revisions: 20230620 first version; 20230622 add model equations and derivations of iteration; 20230703 move the BLA text to its own heading; 20230708 rewrite mose of the text

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 09.

s.27