The arithmetic underneath

The direction the error leans

The size of one rounding error is set by the precision. How ten thousand of them combine is set by something else entirely — the rounding mode — and the fitted exponents are 0.47 for round-to-nearest and 1.01 for round-toward-infinity, on identical data at identical precision.

Worth reading first: What a float can hold · The order they are added in.

Every figure on this site so far takes the same step and none of them shows it. A computed value falls between two representable numbers; one of them is returned. Which one?

The answer this site has assumed for four phases is the nearer, with ties to even, which is what IEEE 754 does by default and what roundBits did unconditionally. It is one of four answers the standard defines, and there is a fifth in silicon on several current accelerators that the standard does not define at all.

The choice barely changes the size of a single error — every mode is within a factor of two of every other on one operation. What it changes is how the errors of many operations combine, and that is a difference of a whole power of n.

Relative error of a sum against the number of terms, at 24 significand bitsFive curves of relative error against the number of terms, both axes logarithmic. Three rise with unit slope; two rise with about half that.10²10³10⁴10⁻⁸10⁻⁷10⁻⁶10⁻⁵10⁻⁴terms summedrelative error+∞, −∞ 1.01zero 1.00stochastic 0.50nearest 0.47√n against n, fittednearest, fitted exponent0.47stochastic, fitted exponent0.5toward +∞, fitted exponent1twelve seeds averaged at each sizethe slope is the bias, not the precision
Fig. 1 Five modes, identical data, identical precision, and two families of curve. Two rise like the square root of the number of terms and three rise linearly. Drag the precision and every curve moves down together without changing slope — until fourteen bits, where the sum stagnates and the distinction disappears.

Both halves of that are claims about a slope holding across a range, so both are read off the range.

Relative error of a sum against the number of terms, at 14 significand bitsFive curves of relative error against the number of terms, both axes logarithmic. Three rise with unit slope; two rise with about half that.10²10³10⁴10⁻⁵10⁻⁴10⁻³10⁻²10⁻¹terms summedrelative error+∞, −∞ 1.03zero 0.97stochastic 0.50nearest 0.51√n against n, fittednearest, fitted exponent0.51stochastic, fitted exponent0.5toward +∞, fitted exponent1twelve seeds averaged at each sizethe slope is the bias, not the precision
Fig. 2 Fourteen bits, the lowest precision at which the separation survives. Fitted exponents: nearest 0.51, stochastic 0.50, toward-zero 0.97, downward 0.97, upward 1.03.

Half a slope and a whole one, with nothing between them. Six bits further up the slider they are still the same two numbers:

Relative error of a sum against the number of terms, at 20 significand bitsFive curves of relative error against the number of terms, both axes logarithmic. Three rise with unit slope; two rise with about half that.10²10³10⁴10⁻⁷10⁻⁶10⁻⁵10⁻⁴10⁻³terms summedrelative error+∞, −∞ 1.00zero 1.00stochastic 0.60nearest 0.58√n against n, fittednearest, fitted exponent0.58stochastic, fitted exponent0.6toward +∞, fitted exponent1twelve seeds averaged at each sizethe slope is the bias, not the precision
Fig. 3 Twenty bits: 0.58, 0.60, 1.00, 1.00, 1.00.

Six bits of precision, and the five exponents are the same five numbers. Eight more and they still are — which is what makes the fourteen-bit floor below the slider a property of the sum rather than of the fit.

Relative error of a sum against the number of terms, at 28 significand bitsFive curves of relative error against the number of terms, both axes logarithmic. Three rise with unit slope; two rise with about half that.10²10³10⁴10⁻⁹10⁻⁸10⁻⁷10⁻⁶10⁻⁵terms summedrelative error+∞, −∞ 1.01zero 1.00stochastic 0.54nearest 0.46√n against n, fittednearest, fitted exponent0.46stochastic, fitted exponent0.54toward +∞, fitted exponent1twelve seeds averaged at each sizethe slope is the bias, not the precision
Fig. 4 Twenty-eight bits: 0.46, 0.54, 1.00, 1.00, 1.01.
Relative error of a sum against the number of terms, at 11 significand bitsFive curves of relative error against the number of terms, both axes logarithmic. Three rise with unit slope; two rise with about half that.10²10³10⁴10⁻⁴10⁻³10⁻²10⁻¹110¹10²terms summedrelative error+∞, −∞ 1.49zero 0.87stochastic 0.59nearest 1.09√n against n, fittednearest, fitted exponent1.1stochastic, fitted exponent0.59toward +∞, fitted exponent1.5twelve seeds averaged at each sizethe slope is the bias, not the precision
Fig. 5 And eleven bits, below the floor: nearest fits 1.09 and upward 1.49 — the unbiased mode looking biased and the biased one looking worse than linear.
bits nearest stochastic toward zero downward upward
11 1.09 0.59 0.87 0.87 1.49
12 0.96 0.58 0.93 0.93 1.18
13 0.54 0.52 0.96 0.96 1.04
14 0.51 0.50 0.97 0.97 1.03
15 0.54 0.61 0.99 0.99 1.01
17 0.52 0.52 1.00 1.00 1.01
20 0.58 0.60 1.00 1.00 1.00
24 0.47 0.50 1.00 1.00 1.01
28 0.46 0.54 1.00 1.00 1.01

Two exponents and six precisions. From fourteen bits up the two unbiased modes sit at 0.46 to 0.61 and the three biased ones at 0.97 to 1.03 — a half and a one, with no overlap and no drift. The distinction is not a feature of double precision or of any particular width; it is the difference between errors that cancel and errors that do not, and it survives a factor of sixteen in the unit roundoff.

And the first three rows are the instrument, which the figure says out loud. At eleven bits nearest fits 1.09 and upward fits 1.49 — the unbiased mode looking biased and the biased one looking worse than linear — because the sum has stagnated and the fit is measuring stagnation. The readout prints below fourteen bits the sum stagnates and the separation is gone at exactly those three stops and not at the others, so the caveat is generated rather than remembered.

One column is a duplicate of another at every precision. Toward-zero and downward read 0.87/0.87, 0.93/0.93, 0.96/0.96, 0.97/0.97, 0.99/0.99 and 1.00/1.00 all the way down — identical to two decimals at nine stops. That is not a coincidence and it is worth saying: every term of this sum is positive, so rounding toward zero is rounding downward, and the figure draws five modes of which four are distinct. A reader counting rounding modes on this picture would count one too many.

The five answers

mode what it returns
nearest round to the nearer, ties to even. The default everywhere.
zero truncate. What a cast to an integer does, and what a great deal of fixed-point hardware does when nobody chose.
up toward +∞.
down toward −∞.
stochastic round up with probability equal to the distance up. Not in IEEE 754; a hardware feature of several recent accelerators and a software trick before that.

The first four are functions of the value. The fifth is a distribution over the same two answers, with a mean equal to the value itself, and it is the next essay’s subject as much as this one’s.

Round-to-nearest and round-toward-zero are the two a programmer meets by accident. The first because it is the default; the second because a truncation is what a cast does, and because integer and fixed-point arithmetic truncate unless told otherwise.

Rounding 1.42 to 4 significand bits, five waysA number line of the representable values in one octave, with the exact value marked between two of them and a row for each rounding mode showing which neighbour it returns.11.131.251.381.51.631.751.882exact 1.42nearest1.375toward zero1.375toward +∞1.5toward −∞1.375stochastic36% upstochastic, measuredposition between neighbours0.36share of draws rounded up0.369 representable values between 1 and 2four functions and one distribution
Fig. 6 One value, four functions and a distribution, on the number line of what four significand bits can hold. The bar at the bottom is the stochastic mode’s split, measured over two thousand draws against the value’s position between its neighbours. Drag the value and watch the three deterministic rows jump between neighbours while the split moves continuously.

Bias, and why it decides the exponent

For round-to-nearest, an inexact result is as likely to be pushed up as down over any reasonable set of inputs, so the errors of n operations partly cancel. Their sum is a random walk and grows like √n.

For a directed mode, every single rounding goes the same way. Nothing cancels; the errors add. The sum grows like n.

That is the whole mechanism and the measurement is a fit. Sum n like-signed values at 24 significand bits, average over twelve seeds so that a random walk is being measured rather than one sample of one, and fit the exponent in error ∝ n^p across three decades of n:

mode fitted exponent
nearest 0.47
stochastic 0.50
toward zero 1.00
toward −∞ 1.00
toward +∞ 1.01

Two families, cleanly separated, on data that differs in nothing. At n = 2,000 the gap between round-to-nearest and round-toward-infinity is a factor of 172.

The data is like-signed on purpose. A sum of values that change sign lets the values cancel, which hides what the errors are doing — and cancellation in the values is the order they are added in’s subject rather than this one’s; the question here is only whether the rounding errors cancel, so the data is arranged not to.

Truncation is biased towards zero, not downwards

Round-toward-zero deserves separating from the other two directed modes, because its bias has a different shape and it is the one most likely to arrive uninvited.

Rounding up is biased upwards on the number line. Rounding down is biased downwards. Truncation is biased towards zero, which means it is biased downwards for positive values and upwards for negative ones — so on a sum of like-signed terms it is a systematic underestimate of the magnitude, and on a sum of mixed signs its errors partly cancel again for a reason that has nothing to do with being unbiased.

Measured at n = 2,000 the signed relative errors are: rounding up +3.98·10⁻⁵, rounding down −4.16·10⁻⁵, truncation −4.16·10⁻⁵ on this all-positive data — identical to rounding down, as it must be — against round-to-nearest’s 2.31·10⁻⁷.

The exponent is the same 1.00 either way. What differs is that the sign of a truncation’s bias depends on the data, so a programmer who checks one test case with positive inputs and finds an underestimate learns nothing about what negative inputs will do.

The one guarantee this site can offer

The directed pair is not only a hazard. Computed together they are the whole of interval arithmetic, and they give something no error estimate on this site can:

Round the entire computation up, and again round the entire computation down. The true answer lies between the two results. Not probably. Not usually.

Every intermediate value in the first run is at least the true one and every intermediate in the second is at most, so the ordering propagates. There is no tolerance in the statement, no constant, and no hypothesis about the data. It is the only claim of that kind anywhere on this site.

And it is loose. The bracket widens by roughly an ulp per operation, so its width grows like nu while the actual error grows like √n·u, and the gap between the guarantee and the thing guaranteed opens by √n.

That last sentence is the claim, and the table it was first written under does not support it. One sum at each size, at 24 bits, gives ratios of 6,827 · 29 · 19 · 173 · 89 · 910 · 182 at n = 10 · 30 · 100 · 300 · 1,000 · 3,000 · 10,000. Those numbers move by two orders in both directions and are not monotone anywhere.

The reason is that the two quantities are not the same kind of quantity. The width is nearly deterministic — it is an ulp per operation, whatever the data. The round-to-nearest error is the random walk this essay’s first half is about, so the ratio inherits the whole of the walk’s variance. The n = 10 entry is the extreme case: that seed’s error is 6.7·10⁻¹¹ against a median of 4.5·10⁻⁸, lucky by a factor of 670, and it produces a ratio in the thousands at exactly the size where the ratio is smallest.

Twenty-four seeds per size, medians:

terms 10 30 100 300 1,000 3,000 10,000
median ratio 11.2 17.9 50.0 78.2 115 264 496
one seed 6,827 29 19 173 89 910 182

Monotone, and it fits: the width’s exponent is 0.988 — an ulp per operation, which is n — the error’s is 0.442, which is the walk, and the difference is 0.546. The gap opens like √n, exactly as the argument says. The argument was right and the evidence originally offered for it could not have shown it either way, which is the failure mode this site keeps meeting from the other side: a number that happens to agree with a correct claim is not evidence for it. It is the same discipline one matrix is an anecdote applies to a single test case, arriving here about a single seed rather than a single matrix.

The containment never fails — that is asserted at every size and at every seed — and by ten thousand terms the interval is five hundred times wider than the error it is bounding. That is the shape of the trade: an interval method is true at every size and useful only while the width stays below what the application can accept, and the width grows faster than the error does.

The remedy is not more precision — and precision is not the only thing a format has, since the other half of a format is the exponent range and it decides a different set of failures. Both quantities here scale with the unit roundoff, so buying bits moves both curves down together and does not close the gap between them — which the figure’s slider shows directly.

The width of a guaranteed bracket, and the error it guarantees, at 24 bitsTwo rising curves on logarithmic axes. The upper one, the bracket width, rises faster than the lower one, the actual error.10¹10²10³10⁴10⁻¹¹10⁻¹⁰10⁻⁹10⁻⁸10⁻⁷10⁻⁶10⁻⁵10⁻⁴10⁻³terms summedrelative widthbracket widthactual errora bound with no tolerance in itbracket contains the truth7width at 10,000 terms4.4·10⁻⁴error it bounds2.4·10⁻⁶every sum computed twice, up and downtrue at every size, useful at none of the large ones
Fig. 7 The guarantee and the error it guarantees, on the same axes. Both rise, the upper one faster, and the slider moves them down together without changing the distance between them.

What one rounding costs, and why it is the wrong question

It is worth being precise about the quantity the modes do share, because the essay’s claim is a claim about the one they do not.

For a single operation, round-to-nearest has an error of at most half an ulp and every other mode here has an error of at most a whole ulp. So the worst-case error of one rounding differs by exactly a factor of two across the whole set, and every textbook bound on a computation — the nu of a dot product, the backward error of an elimination, the ‖A − QR‖/‖A‖ this site prints on every factorisation figure — inherits that factor of two and nothing else.

Those bounds are therefore almost mode-independent, and they are worst-case bounds: they assume every rounding goes the wrong way, which is exactly what a directed mode does and what round-to-nearest does not. So a directed mode attains the classical bound and round-to-nearest sits far below it — which is why a small residual is not a small error can be read two ways depending on which mode produced the residual, and the distance between them is the √n against n this essay measures.

That reframes the usual complaint that error bounds are pessimistic. The bounds are not pessimistic about arithmetic in general — they are tight for a computation that rounds in one direction, and loose for one that rounds to nearest, by a factor of √n. A reader who has been told that “in practice errors grow like √n rather than like n” has been told a fact about the rounding mode, and it is usually presented as a fact about floating point.

Changing the mode is not a local decision

The measurements above treat the rounding mode as a parameter, because in this site’s arithmetic it is one: arith(bits, { mode }) returns a context, every routine takes a context, and two contexts can coexist in the same program. That is a simulation choice and it is not what the hardware offers.

On a real processor the rounding mode lives in a control register. It is a property of the thread, not of the operation, so setting it changes the behaviour of everything running on that thread — including library code that was written and tested under the default. A library that computes a special function by a polynomial tuned for round-to-nearest, or that relies on a particular splitting identity being exact, may silently return something else when the register changes underneath it. Compilers add a second problem: an optimiser is entitled to reorder or constant-fold arithmetic unless told the rounding mode may vary, and under a directed mode the reordering is no longer harmless.

So interval arithmetic in production is not two runs with a register flipped between them. It is either a library that sets and restores the mode around each operation, at a substantial cost per operation, or one that emulates the directed results from round-to-nearest ones — which is arithmetic on top of arithmetic and slower again.

That is worth knowing before reading the guarantee in the previous section as a cheap one. The bracket costs at least twice the work and usually much more, and its width grows like n. The reason it is used anyway, where it is used, is that no other technique on this site can say the answer is in here without a hypothesis.

The mode nobody chose

There is one more reason to care, and it is the one this phase actually arrived through.

An eight-bit format is not usually the arithmetic. It is the storage, and the operation that matters is the quantisation: taking a value in double or single precision and choosing the nearest representable eight-bit one. That is a rounding, it happens once per stored value rather than once per arithmetic operation, and its mode is decided by whatever the quantisation routine does.

A routine that truncates — because it was written with a cast, or because the fixed-point hardware it targets truncates — introduces a bias towards zero in every stored value at once. The values of a tensor are then all slightly smaller in magnitude than they should be, systematically, and the sums computed from them inherit an error growing like n rather than √n. At 24 bits that is the factor of 172 measured above; at four significand bits the individual errors are 2⁻⁴ rather than 2⁻²⁴, and whatever the factor is, it multiplies something twenty decades larger.

Nothing about that is exotic and nothing warns. The quantised values are plausible, the norms are close to right, and the defect is a systematic shrinkage that a spot check on a few entries will not find. It is the same shape as the truncation bias above and it arrives through a different door — not from an arithmetic loop but from a conversion nobody thought of as arithmetic.

The block-format essay quantises 640 values at a time and rounds to nearest throughout, deliberately, so that the errors it measures are the format’s and not the rounding’s.

Where the distinction stops holding

The separation between the two families has a floor, and the floor is worth knowing because a reader who takes the exponents away as constants will apply them where they are false.

At fourteen significand bits and above, the fits are 0.5 and 1.0 as advertised. Below thirteen they collapse together. At eleven bits a sum of ten thousand terms of size one has a running total near ten thousand, whose ulp at that precision is 8 — so the addends stop arriving entirely, and what is being measured is no longer an accumulation of rounding errors but the size of what was never added.

Round-to-nearest fits at 1.09 there and rounding up at 1.49, both above the linear rate. The distinction this essay is about has not weakened; it has been replaced by a different failure, and the fitted exponent is reporting that failure rather than the one it was set up to measure.

That failure is stagnation, and it is where the fifth mode earns its keep.

What is asserted, and what refuses

Every mode does what its name says, on eight values spanning positive, negative and tiny: rounding up is never below the value, rounding down never above, the bracket is at most one ulp wide, and truncation agrees with rounding-down for positive values and with rounding-up for negative ones — four properties, each of which fails for at least one other mode, so the set distinguishes rather than merely accepting.

The exponents separate: below 0.8 for the two unbiased modes and above 0.85 for the three biased ones, at every precision from fourteen bits up.

And the collapse below thirteen bits is asserted too, in the other direction: round-to-nearest’s exponent must be above 0.9 there and rounding up’s above it again. A slider that quietly stopped checking at one end would be drawing the wrong claim there.

The bias has a direction, asserted by sign rather than by magnitude, on all three biased modes.

The bracket contains the truth at every size on the sweep, and is more than ten times wider than the error it bounds at the largest.

The claim that precision alone decides how an error accumulates is refused. Two fits at the same 24 bits, in two modes, are fed to assertClose inside rejects() — 0.47 against 1.01, and the assertion must fail.

And the bracket used backwards is refused. Rounding down to get an upper bound and up to get a lower one produces an interval containing nothing, and the assertion that it brackets the truth is fed to rejects().

What is left

Stochastic rounding, which is the next essay: unbiased like round-to-nearest, and the only mode that escapes the stagnation this one ends on.

Interval arithmetic as a method. What is measured here is a two-run bracket on a sum. A real interval library propagates intervals through every operation, which is tighter in some places and looser in others, and is a subject rather than a comparison.

The rounding of a division or a square root, which this site’s simulated arithmetic performs by computing in double and rounding once — correct for addition and multiplication of narrow operands and not exactly right for the other two, as assertMatchesHardware measures rather than assumes. A directed mode makes that double-rounding gap slightly more visible, and nothing here quantifies it.

Round-half-away-from-zero, the fifth IEEE mode, which differs from round-to-nearest only on exact ties. Ties are a measure-zero event on random data and are common on the decimal fractions a person types, which is why the two modes disagree constantly in a spreadsheet and never in this essay’s experiments. Measuring it would need data chosen to land on ties, which is a different experiment rather than another column in this one.

Where a rounding bias is squared, repeatedly

A rounding mode’s bias matters most where the same operation is applied over and over to its own output. The clearest instance on this site is the squaring phase of a matrix exponential, where a result is multiplied by itself s times and each squaring doubles whatever error it was carrying.

What links here

Computed from the collection, not written here: the essays that point at this one.

Reads more easily once this is understood

Essays that name this one as worth reading first.

Shares its objects with

Essays that name at least two of the same things, and that neither author linked.

Named objects

A flat tag is an object no other essay names yet.

Directed roundingError accumulationIEEE 754Interval arithmeticRounding modesSummationTruncationUnit roundoff