A square that evaluates negative

The two essays before this one end on verdicts that two builds disagree about, and in both cases there is a right answer: the determinant is 1, the matrix is definite or it is not, and an exact computation settles it. This one has no right answer to reach for, and that is the point of it.

Two traces of noise. The polynomial is (x − 1)⁶ expanded, so its value is a sixth power and cannot be negative; within 10⁻³ of the root it is under 10⁻¹⁸, which is below the rounding of the terms that produced it. One build returns a negative value at 179 of the 401 points and the other at 196, and the two disagree about the sign at 98 — a quarter of the window.

Why this is a different failure #

In one multiply the compiler removed the fused form is exactly right and the other is exactly wrong, because the algorithm was designed around the instruction. In a matrix that is definite on one machine neither form is reliably right, but one of them is right on each matrix and an exact computation says which.

Here neither trace is nearer the truth in any useful sense. The true value is 10⁻¹⁸ and both evaluations return quantities around 10⁻¹⁵, so both are wrong by three orders of magnitude, and asking which is closer to zero is asking which random number is smaller. The largest difference between the two evaluations across the window is 1.44·10⁻¹⁵ — six roundings, exactly what the bound allows — and both are entirely made of them.

What the picture shows is therefore not an accuracy comparison. It is the shape of the region in which a computed sign belongs to the build, and that region has a width, which is computable.

Where the window’s edge is #

The window is where |p(x)| falls below the rounding of the terms that make it up. Horner’s rule evaluating a degree-d polynomial commits about d roundings, each of size u times the largest partial value it holds, so the noise floor is roughly

u · max over the recurrence of |sₖ| ≈ u · Σ|cᵢ| |x|ⁱ,

which near x = 1 is u times the sum of the coefficients’ magnitudes: 1 + 6 + 15 + 20 + 15 + 6 + 1 = 64, giving about 1.4·10⁻¹⁴. The true value reaches that at |x − 1| = (1.4·10⁻¹⁴)^(1/6) ≈ 4.9·10⁻³.

Which is what the sweep measures.

The disagreement count rises by a factor of fifty across four decades of window and then stops. At windows of 10⁻¹, 10⁻², 10⁻³, 10⁻⁴ and 10⁻⁵ the two builds disagree at 2, 29, 98, 106 and 97 of 401 points, and the counts of negative values run 7/6, 53/55, 179/196, 184/196 and 186/205.

The saturation is the part worth having, and it confirms where the edge is. Between 10⁻¹ and 10⁻³ the disagreement grows fiftyfold; between 10⁻³ and 10⁻⁵, across two further decades, it does not grow at all — 98, 106, 97, which is one number with sampling noise on it. That is exactly what a noise floor predicts: once the true value is entirely below the rounding, going closer to the root cannot make the sign more arbitrary than arbitrary. The predicted edge at 4.9·10⁻³ sits precisely where the count stops climbing, so the estimate above is not a rough order of magnitude — it locates the transition to within a factor of two, and it was computed from the coefficients alone.

That number — Σ|cᵢ| divided by |p(x)| — is exactly the summation condition number of the vector that hides it, applied to the terms of a polynomial. The essay’s census puts a residual at 10¹⁷ and this evaluation is the same object: an expression whose terms are enormous relative to their sum.

What a root finder does with it #

The consequence is not that a value is wrong. It is that a branch is taken differently, and the branch is inside every algorithm that locates something by sign.

A bisection asks for the sign at the midpoint and discards half the interval. Inside the window that sign is a property of the build, so two builds discard different halves, converge to different points, and return roots that differ by the width of the window. Their residuals are the same size — both are 10⁻¹⁵, which is what convergence looks like here — and neither is wrong by the standards a root finder applies to itself.

This is the same mechanism the roots are not the coefficients measures for a root’s conditioning: a double root moves like √ε, a root of multiplicity m like ε^(1/m), so half the digits are lost at a double root and five sixths at a sixth-order one, with nothing ill conditioned anywhere. The window here is that fact seen from the other side. The conditioning says the root cannot be located better than 10⁻³ from data at the rounding level; the window says that inside that 10⁻³ the arithmetic’s own decisions are what a bisection is following.

The two agree, and they are different statements. The conditioning is a property of the problem and it is the same for both builds. The window is a property of the arithmetic, and inside it the two builds take different paths through a computation that both complete successfully.

Counting the disagreements is the measurement #

The number worth quoting from the hero figure is not 179 or 196 — those are counts of a wrong sign, which is interesting once. It is 98, the number of points at which the two builds disagree with each other, because that is the quantity a person maintaining two machines actually meets.

It is worth seeing why it is neither 0 nor 401. It is not zero because the two evaluations differ by up to six roundings and the values are the size of one, so a disagreement is likely wherever the noise happens to straddle zero. It is not 401 because the two traces are strongly correlated: they share every input and differ only in how each of six multiply-adds was rounded, so most of the time they land on the same side of zero by accident of sharing most of their arithmetic.

A quarter is therefore the honest expectation, and it is stable across the window: measured, the disagreeing fraction is 0.244 at a reach of 10⁻³, 0.264 at 10⁻⁴ and 0.242 at 10⁻⁵, and it collapses to 0.005 at 10⁻¹ because the true value has come back above the noise. That stability is what makes the width of the window — rather than the count — the thing to compute in advance.

The general shape, since three essays have now produced it #

This is the third instance of one pattern and it is worth writing down before the field moves on to reductions again.

A computation produces a real number. The number is a difference of quantities much larger than itself, so its leading digits are rounding. A consumer compares it with something — zero, a tolerance, another number — and turns it into a discrete outcome. Two conforming builds compute the rounding differently, so the outcome differs, and every diagnostic anybody applies afterwards reports that both builds behaved correctly.

The three instances differ only in what is being compared and whether there is an adjudicator:

The third row is the one that closes the sequence, because it is where the site’s usual move stops working. This collection’s habit is to answer a disagreement by computing the quantity exactly and saying who was right. Here the exact value is 10⁻¹⁸, it is not zero, and it is positive — so the “right” answer is that every negative value in the figure is wrong, which is true and useless, because no build was ever going to return it and no algorithm that consumes the sign can act on it.

What the exact value does establish is where the window is, which is the thing to compute.

Three things that do not fix it #

Higher precision. It moves the edge: quadruple precision would put the noise floor at 10⁻³² and the window at |x − 1| ≈ 5·10⁻⁶. The window is narrower and it is still there, and inside it the same thing happens. This is the shape the whole field keeps producing — precision moves a boundary and does not remove one.

A tolerance. Declaring |p| < 10⁻¹⁴ to be zero converts a sign into three outcomes, positive, negative and cannot tell, which is genuinely better — but a bisection cannot use cannot tell without an answer for what to do next, and the honest answer is to stop and report an interval. Most implementations do not have that outcome and the ones that do return an interval of width 10⁻³, which is the conditioning again.

Turning off the contraction. It makes the trace reproducible between builds that share the flag, and it does not make it correct, and it makes it the less accurate of the two forms elsewhere. It is the option that most looks like a fix and buys the least.

What does work #

Two things, and both are changes of formulation rather than of arithmetic.

Do not expand the polynomial. (x − 1)⁶ evaluated as a sixth power of (x − 1) is exact to a rounding at every point, has no cancellation anywhere, and is never negative. The whole failure is a consequence of writing a polynomial in a basis in which its value is a difference of large terms — which is what a monomial expansion does to any polynomial with a root away from zero, and what the roots are not the coefficients is about when it measures a Chebyshev basis against a monomial one.

Or evaluate with an error-free transformation. Compensated Horner — the same twoProduct and twoSum machinery the reproducible summation policies use — evaluates a polynomial as though in twice the precision, at about four times the cost, and pushes the noise floor to u² instead of u.

That does not remove the window, and how far it narrows it is worth running rather than quoting. Sweeping the negative fraction over the reach:

   reach        plain     compensated
   10⁻²         0.126        0.000
   3.2·10⁻³     0.414        0.000
   10⁻³         0.471        0.002
   10⁻⁴         0.433        0.015
   3.2·10⁻⁵     0.466        0.039
   10⁻⁵         0.469        0.162
   3.2·10⁻⁶     0.454        0.298

The compensated window’s edge is about 10⁻⁵, where the negative fraction first becomes substantial. The essay’s own noise-floor formula gives that when it is applied: (u²·Σ|cᵢ|)^(1/6) is 9.6·10⁻⁶, against 4.4·10⁻³ for the uncompensated (u·Σ|cᵢ|)^(1/6) and a measured 5·10⁻³. Both predictions land, and the narrowing is u^(−1/6) = 457 in the abscissa.

There is one more thing in that table and it undoes the sentence this section wanted to end on. The coefficients here — 1, −6, 15, −20, 15, −6, 1 — are small integers, exactly representable, so the data carries no error at all and there is no conditioning limit for the arithmetic to be measured against. The root’s ε^(1/6) sensitivity is a statement about perturbing coefficients, and nothing has perturbed them. So the remaining window after compensation is not the conditioning taking over from the arithmetic; it is the arithmetic, four hundred times narrower and otherwise unchanged in kind.

Which is the honest version of what compensation buys. It moves the boundary by u^(−1/6) and it does not change what is on either side of it, and the sixth root is what makes a factor of 10¹⁶ in the noise floor into a factor of four hundred in the window. A fractional power is as unkind to a repair as it is to a perturbation, and this is the same 1/n exponent that makes the perturbation analysis in the eigenvalues that are not there so flat, working against the reader rather than for them.

What the two numbers together are worth to a caller #

Putting the confirmed claim beside the corrected one gives the section a shape it did not have.

The width of the window is predictable to within a factor of two from three things a caller already has: the coefficients, the precision, and the multiplicity of the root. That is the number to compute, and it is one line — (u·Σ|cᵢ|·|x|ⁱ ⁄ |leading behaviour|)^(1/m), with the sixth root here being what makes it so much wider than the noise floor suggests.

The fraction of points inside it where two builds disagree is a quarter, and it does not move. That is the number nobody can compute in advance and nobody needs to: it is not far from a half, it is not close to zero, and any policy that survives a quarter survives whatever the real figure is on a different polynomial.

So the practical instruction is entirely about the first. Compute the window; if the query point is inside it, the sign is not information and no amount of care in the evaluation changes that. If it is outside, both builds agree and the question was never interesting. There is no middle case in which being careful about the arithmetic decides anything — which is the uncomfortable part, because being careful about the arithmetic is the reflex this whole subject trains.

Both are the same recommendation this site keeps arriving at from different directions: when a computed quantity is a difference of large terms, the repair is upstream of the arithmetic.

The window is not a defect of Horner’s rule #

Worth heading off, because the obvious response to the figure is to blame the algorithm.

Horner’s rule is the best of the standard evaluations by every measure the subject applies. It uses d multiply-adds where the naive sum of powers uses about 2d operations, its backward error is optimal — the computed value is the exact value of a polynomial whose coefficients differ from the given ones by a relative γ₂d — and it is what every library uses. The failure here is not a failure of Horner’s rule; it is the failure the backward error statement predicts.

A polynomial whose coefficients are perturbed by 10⁻¹⁵ relative has roots that move by (10⁻¹⁵)^(1/6) ≈ 3·10⁻³ for a sixth-order root, and inside that neighbourhood the perturbed polynomial takes both signs. Horner returns the exact value of one such polynomial, and the fused build returns the exact value of a different one. Both are correct, both are answers to slightly different questions, and the questions differ by less than the data can distinguish.

That is the site’s founding identity applied to an evaluation rather than to a solve, and it puts this essay’s window exactly where three errors and one number would predict: the conditioning of the root is enormous, the backward error of the evaluation is at the rounding level, and their product is the width of the region in which nothing can be said.

Where else the sign of a nearly-zero quantity is read #

Naming them, because a reader will recognise more of these than they expect.

A discriminant. b² − 4ac decides whether a quadratic’s roots are real, and near a double root it is the difference of two nearly equal products — exactly the determinant of one multiply the compiler removed with different letters.

A pivot. When symmetry is not enough chooses between a 1 × 1 and a 2 × 2 block by comparing computed quantities whose difference decides the branch.

A convergence test on a difference of iterates. |xₖ₊₁ − xₖ| is a difference of nearly equal things by construction, and its comparison with a tolerance is a tolerance that reads its own residual.

A determinant used to decide orientation, a curvature used to decide a step, a sign used to decide a reflection. Every one is a real number computed from cancelling terms and consumed as a bit.

The field’s general statement covers all of them: a stopping test is a race prices the third, and the pattern is that the arithmetic’s disagreement is small, the consumer is discrete, and the composition is a branch two machines take differently.

What to report instead of a sign #

The recommendation the three FMA essays share, stated once here because this is where it is unavoidable.

A routine that computes a quantity and consumes its sign should return the quantity and a scale, not the sign. For the polynomial: return p(x) and Σ|cᵢ||x|ⁱ, from which a caller computes the noise floor and decides whether the sign means anything. For a Cholesky: return the pivot and the matrix norm. For a determinant: return the value and the sum of the magnitudes of its two products.

In every case the second number costs one accumulation the loop was nearly doing anyway, and it converts an unanswerable question — is this sign real — into an arithmetic one. That is the same move rank is a decision makes for a rank: the honest output is the spectrum with the threshold drawn on it, and the count above the threshold is the caller’s.

The reason it is not standard is that a sign has an obvious return type and a scale does not have an obvious name. sgn is one bit and (value, scale) is a struct, and a routine that returns a struct where the mathematics says sign looks like it is hedging. It is not hedging: the sign genuinely is not determined by the inputs at this precision, and returning one is a claim the arithmetic cannot support.