A double root halves the error instead of squaring it

Where Newton’s method goes instead coloured the plane by which root each starting point reaches, and the essays that followed it asked which starting points reach a root at all. None of them asked how long the trip takes. For the three cube roots of unity the answer is short everywhere except near the fractal boundary: a start within reach of a root gets there in a few steps, because near a root Newton’s method squares its error.

That squaring is not a law of nature. It depends on the root being simple — on the polynomial crossing zero with a non-zero slope. At a root where the curve only touches zero, a double root, the method still converges, but its character changes completely. How fast the staircase arrives noted in a sentence that the rate becomes linear; drawn over the whole plane, that sentence turns into most of the picture.

Squaring, drawn #

For $z^3 - 1$ the map is the familiar three-armed pattern seen from a new angle. Around each root is a pale disc of starting points that arrive within three steps, then a ring that needs four or five, and the slow starts — ten, fifteen, twenty-five steps — are confined to thin filaments along the boundaries between basins, where a start is flung far away before it settles. Averaged over the square, a start needs 6.93 steps to come within $10^{-7}$ of a root.

The pale discs are the squaring made visible. If the error is $e$, the next error is about $C e^2$ for a constant $C$ set by the polynomial, so an error of $10^{-1}$ becomes $10^{-2}$, then $10^{-4}$, then $10^{-8}$. The number of correct digits doubles every step, and a start that is close enough needs almost nothing more. The constant is $C = p''(r)/(2p'(r))$, so a root at which the curve is steep and nearly straight is reached fastest of all, and for $z^3 - 1$ it is the same at all three roots, which is why the three pale discs are the same size. Outside the discs, the count grows by one for each halving of the distance to the disc’s edge — until the boundary is reached, where a start’s first step can send it anywhere.

Halving, drawn #

The polynomial $(z - 1)^2(z + 1) = z^3 - z^2 - z + 1$ has a simple root at $-1$ and a double root at $1$. Its speed map, at the top of the essay, is another picture entirely. Around the simple root at $-1$ there is the usual pale disc. Around the double root at $+1$ there is none: the whole right-hand side of the plane is a single band of sixteen to twenty-five steps, and the average over the square rises from under seven to 18.80.

The reason is a short calculation. Near a root $r$ of multiplicity $m$, the polynomial behaves like $c(z - r)^m$ and its derivative like $mc(z - r)^{m-1}$, so the Newton step is

$$z - \frac{p(z)}{p'(z)} \approx z - \frac{z - r}{m}, \qquad\text{new error} \approx \left(1 - \frac1m\right) \times \text{old error}.$$

For a simple root, $m = 1$, the factor is zero — which is why the next term, the square, takes over. For a double root the factor is $\tfrac12$. Each step halves the error, and halving from $1$ to $10^{-7}$ takes 23 steps however the start was chosen. The flat band on the right of the map is that count, the same almost everywhere, because once linear convergence sets in, the only thing that varies from start to start is how far away it began. Twenty-odd steps is simply the price of one double root.

The error curves show the same thing one start at a time. On a logarithmic scale a squaring error plunges — its exponent doubles each step — while a halving error descends in a straight line of slope $\log_{10}\tfrac12 \approx -0.3$ per step. A triple root is slower still, with factor $\tfrac23$: more than thirty steps to gain five digits. Each measured run settles to the predicted factor, $\tfrac12$ and $\tfrac23$, to two decimal places.

The Newton map’s own slope #

Everything above can be said in the language of iteration, and saying it that way explains why the factor is $(m - 1)/m$ rather than something else. Newton’s method is the iteration of a map, $N(z) = z - p(z)/p'(z)$, and a root of $p$ is a fixed point of $N$. How fast an iteration approaches a fixed point is decided by the map’s derivative there: derivative $\lambda$ means the error is multiplied by $\lambda$ each step, and derivative 0 means the linear term vanishes and the square takes over.

Differentiating, $N'(z) = p(z)\,p''(z)/p'(z)^2$. At a simple root $p(r) = 0$ and $p'(r) \ne 0$, so $N'(r) = 0$: the fixed point is superattracting, the name for the case where the square rules. At a root of multiplicity $m$, numerator and denominator both vanish, and substituting $p \approx c(z - r)^m$ gives $N'(r) = m(m - 1)/m^2 = (m - 1)/m$. So a double root is an attracting fixed point of slope $\tfrac12$ — attracting, but ordinarily so, like the fixed point of $x \mapsto x/2$.

The doubled step changes the map to $N_2(z) = z - 2p/p'$, whose slope at a root of multiplicity $m$ is $1 - 2/m$. At the double root that is 0, superattracting again; at a simple root it is $-1$, a fixed point that neither attracts nor repels, around which the iteration swings back and forth. The map’s slope at each fixed point is the whole story of the speed maps, and the doubled step moves one slope to 0 by moving another to $-1$.

A ceiling on the digits #

The double-root curve does something else: it stops. Near $10^{-8}$ the error stops falling and wanders, and no further steps improve it. This is not the method’s fault, and it is not a failure of convergence. It is rounding.

A computer evaluates $p(z)$ with a relative error of about $10^{-16}$ in each operation. Near a simple root, $p(z)$ is about $c(z - r)$, and it stays distinguishable from rounding noise until $|z - r|$ is about $10^{-16}$ — full precision. Near a double root, $p(z)$ is about $c(z - r)^2$, which drops below the noise when $|z - r|$ is about $\sqrt{10^{-16}} = 10^{-8}$. Closer than that, the computed value of $p$ is noise, the Newton step is noise, and the iteration wanders in a cloud of radius about $10^{-8}$ around the root.

A root of multiplicity $m$ can be located to only about $16/m$ digits in double precision, by any method that evaluates the polynomial, because the polynomial itself cannot see closer. That is a property of the problem, not of Newton’s method: a double root is ill-conditioned, and a perturbation of the coefficients by $10^{-16}$ moves it by $10^{-8}$. The speed maps use a target of $10^{-7}$ precisely so that the double root can be reached at all.

A root that is almost double #

Real polynomials rarely have exact double roots, because exact coincidences do not survive measurement. They have roots that are close — and close roots behave like a double root until the iteration gets near enough to tell them apart. The map for $(z - 1)(z - 1.1)(z + 1)$ shows it. Each of the two roots near 1 has a small pale disc of its own, 0.1 wide, where squaring applies. Outside those discs, the two roots are seen from a distance at which they merge, and the approach is the halving approach of a double root until the gap between them is resolved. The boundary between their basins runs as a vertical wall between them, with a fractal fringe that spreads to the top and bottom of the square.

Between the two close roots sits a root of the derivative — the roots of the slope stay inside the region the roots span, and here it lies just between 1 and 1.1 — and the Newton step, which divides by $p'$, is enormous there. A start near it is thrown far across the plane, which is why the wall between the two basins is fringed with the same fractal filaments that separate the basins of $z^3 - 1$. The closer the two roots, the closer the critical point sits to both, and the narrower the region in which either can be reached quickly.

Whether a polynomial has a repeated root at all is decided by a single number computed from its coefficients, the discriminant, which what the coefficients already know computes for the cubic: it vanishes exactly when two roots coincide, and is small when two are close. So the halving is not an artefact of the exact double root. It is what Newton’s method does at any distance from a cluster of roots that is large compared with the cluster’s width, and it is the typical behaviour for polynomials whose roots crowd together — the situation where two roots run into each other describes from the side of the coefficients.

Doubling the step, and what it breaks #

If the multiplicity $m$ is known, the repair is obvious: the step fell short by a factor of $m$, so take $m$ times the step,

$$z - m\,\frac{p(z)}{p'(z)}.$$

Ernst Schröder proposed it in 1870. At a double root with $m = 2$ the linear factor becomes $1 - m \cdot \tfrac1m = 0$ and squaring returns. The error curve for the doubled step at the double root plunges like the simple root’s, to the rounding floor in four steps.

But a polynomial has more than one root, and the doubled step is wrong for all the others. At the simple root $-1$, where the correct factor is $m = 1$, doubling the step gives a linear factor of $1 - 2 = -1$: the error is reversed rather than reduced, and the iteration oscillates across the root without approaching it. The speed map shows the price. The double root’s basin is fast now, with an average of 5.32 steps among the starts that arrive — but 14.8% of the square never arrives at all, and it is the region that should have been the easy basin of the simple root.

A repair tuned to one root breaks the others, because the method has one step size for the whole plane and the roots need different ones.

A repair that needs no multiplicity #

There is a way out that needs no knowledge of $m$ and works at every root at once. The function $p/p'$ has a simple root wherever $p$ has a root of any multiplicity: near a root of multiplicity $m$ it behaves like $(z - r)/m$, which crosses zero with a non-zero slope. So apply Newton’s method to $u = p/p'$ instead of to $p$. Every root of $p$ is a simple root of $u$, and the squaring returns everywhere — at the price of computing $p''$, since $u' = 1 - p\,p''/p'^2$.

The other way to remove the multiplicity is to divide it out. The greatest common divisor of $p$ and $p'$ contains exactly the repeated factors — $(z - 1)$ for $(z - 1)^2(z + 1)$ — and dividing $p$ by it leaves a polynomial with the same roots, all simple. That gcd is computed by Euclid’s algorithm on polynomials, the same computation the remainders that count the roots uses, where a chain ending in a non-constant polynomial is the signal that a root repeats. The method that counts the roots also removes their multiplicity.

Both repairs work in exact arithmetic. In floating point, where coefficients are measured rather than exact, “the gcd of $p$ and $p'$” becomes a question of how close to a common root two computed polynomials must come before they are declared to share one, and the answer is a tolerance, not a theorem.

What the speed maps cannot show #

The maps are sampled. Each is a grid of 130 by 130 starting points, each coloured by the steps it took, and the fractal filaments along the basin boundaries are finer than any grid. A band drawn as a single colour contains points that take many more steps, and the averages are averages over the grid.

The target is a choice. Every count depends on stopping at $10^{-7}$, a target chosen so that the double root, whose rounding floor is near $10^{-8}$, can be reached. A target of $10^{-12}$ would show the double root’s basin as never arriving, which is true of the arithmetic and false of the mathematics.

The rounding floor is a property of the arithmetic. In exact rational arithmetic, or with more digits, the double root’s error curve would continue its straight descent for as long as the precision allowed; the floor at $10^{-8}$ belongs to the 53 bits of a double-precision number, and a computation in 106 bits would move it to about $10^{-16}$.

And the error curves are one start each. They start from 3 on the real line, and their rates are the limiting rates near the root; far from the root the method does something else — for the triple root it takes a few steps before the $\tfrac23$ factor settles in — and the curves show the transition rather than hide it.

Still open: the best way to find all roots at once #

Newton’s method finds one root from one start, and finding all roots of a polynomial means either many starts — covering rather than avoiding gives a guaranteed set of them — or a method that iterates on all the roots simultaneously. Methods of the second kind, such as the Weierstrass–Durand–Kerner iteration and Aberth’s method, converge quadratically or cubically near simple roots and are what numerical libraries use.

Whether they converge from essentially every starting configuration is not known in general. For the Weierstrass–Durand–Kerner iteration the answer turned out to be no: Reinke, Schleicher and Stoll showed in 2020 that for some cubic polynomials a whole open set of starting configurations never converges — the all-roots analogue of the trapped region an area that never finishes found for Newton’s method. For Aberth’s method, the one most libraries use, no comparable result is known in either direction. A related question is the cost of accuracy near a cluster: how many steps a method needs, as a function of how close the roots are and how many digits are wanted. For Newton’s method the answer is known in outline — the halving phase lasts for about $\log_2$ of the ratio between the distance to the cluster and its width — and for the simultaneous methods it is much less clear. The one-root question has an answer — Newton’s method has open sets of failure, and a finite covering set avoids them — and the all-roots question, which is the one that matters in practice, does not yet.

Speed as a picture #

A basin map says where a starting point goes; a speed map says how fast, and the two tell different stories. The three cube roots have speed maps whose slowness lives on the fractal boundary. A double root has one whose slowness is everywhere in its basin, because the method’s step is the right size only at simple roots, and near a double root it is half the right size at every step.

The speed maps are therefore a second kind of picture of the same dynamics: the basin maps colour a start by where it ends, these by how the map’s slope treats it on the way. At a simple root the slope is zero and the map is fast from anywhere nearby; at a double root the slope is a half and the map is slow from everywhere in the basin; at a cluster of close roots the slope is nearly a half until the cluster is resolved.

When an iteration’s fixed point is not simple, its rate is the first thing to lose — squaring becomes halving, a precision of sixteen digits becomes eight — and the repair has to address the whole map rather than one fixed point, which is why Newton’s method on $p/p'$ succeeds where a doubled step does not.