Re: [musl] Wrong rounding in printf when precision is not set to max

[Rich Felker <dalias@libc.org>]

On Tue, Jan 02, 2024 at 11:16:15AM +0100, Florian Weimer wrote:
> * Rich Felker:
> 
> > Note that, even then, there are lots of gotchas. The x87's "single"
> > and "double" modes are not actually IEEE single and double, but
> > nonstandard types with the IEEE significand precision but excess
> > exponent range that's truncated via double-rounding whenever the
> > intermediate is spilled to memory.
> 
> Doesn't setting the control word avoid double rounding?

No. The first rounding is of the exact mathematical result to
x87-single (23-bit significand, 16-bit exponent) or x87-double (52-bit
significand, 16-bit exponent), then there's a second rounding on store
to IEEE single or IEEE double. The result differs on overflow or
underflow. For example:

- a+b-c (including special case a+b-b) can be finite in x87-double
  (even after final store), but infinite in IEEE double.

- a+b-c (including special case a+b-b) can be exact in x87-double, but
  incur rounding in IEEE double due to a+b becoming denormal.

If you can be certain no intermediates will ever overflow or underflow
(and no denormals appear in inputs), then I believe x87 reduced
precision modes can faithfully emulate IEEE arithmetic. But this is a
nontrivial condition to satisfy in general.

Rich