Cayley Never Saw This

Newton's method applied to z⁵ − 1 = 0 over the complex plane. Each pixel's color shows which of the five roots the iteration converged to; brightness reflects convergence speed. Tap or click to zoom in. Double-tap or double-click to reset.

In 1669, Isaac Newton described a method for finding roots of polynomials. Start with a guess. Evaluate the function and its derivative. Step toward zero. Repeat. The method converges fast, quadratically: each iteration roughly doubles the number of correct decimal places. It became one of the most widely used algorithms in numerical computing, and most of its users never thought much about what it looks like.

In 1879, Arthur Cayley wrote a short paper asking a simple question: given Newton's method applied to a complex polynomial, which starting points converge to which root? For a quadratic, it turns out to be clean. The complex plane splits in two along a straight line, and each half flows to one root without ambiguity. Cayley extended the analysis to cubics, a polynomial with three roots, and found the problem "seems to present considerable difficulty." He published the quadratic case and left the rest open. He had stumbled into a fractal boundary and had no tools to describe it.

The difficulty Cayley noticed: a third basin appears near the boundary between any two basins of attraction. And near that boundary, the first two reappear. The structure repeats at every scale. In the 1970s and 1980s, as raster displays arrived and Benoit Mandelbrot was cataloguing fractal geometry, John Hubbard worked out the topology of Newton basins and showed that the boundary between any two basins is simultaneously the boundary of every other basin. Every boundary point is a limit point of all regions at once. The fractal is not just complicated, it is, in a precise sense, irreducible.

This visualization uses the polynomial z⁵ − 1 = 0, whose five roots sit at the vertices of a regular pentagon on the unit circle. For each pixel, Newton's iteration runs until it converges within 10⁻⁵ of one of the roots. The pixel gets that root's color, dimmed proportional to how many iterations it needed. Fast convergence produces bright color; slow convergence, near the fractal boundary, fades toward black. The rendering walks row by row, so you can watch the image build.

The core of the computation is complex arithmetic: building up powers of z and then applying one division per step:

// z^2, z^3, z^4, z^5 from z = (re, im)
const br = re*re - im*im,    bi = 2*re*im;           // z²
const cr = br*re - bi*im,    ci = br*im + bi*re;     // z³
const z4r = cr*re - ci*im,   z4i = cr*im + ci*re;    // z⁴
const z5r = z4r*re - z4i*im, z5i = z4r*im + z4i*re; // z⁵

// Newton step: z -= (z⁵ - 1) / (5z⁴)
const nr = z5r - 1, ni = z5i;
const dr = 5*z4r,   di = 5*z4i;
const denom = dr*dr + di*di;
re -= (nr*dr + ni*di) / denom;
im -= (ni*dr - nr*di) / denom;

Twelve multiplications and a division per step. That's the entire dynamical system. Cayley could write it down in 1879. What he couldn't do was run it a million times and look at the result. The algorithm was 200 years old before anyone saw what it actually does.