diff options
| author | Szabolcs Nagy <nsz@port70.net> | 2018-04-01 20:02:01 +0000 |
|---|---|---|
| committer | Rich Felker <dalias@aerifal.cx> | 2018-04-02 18:57:41 -0400 |
| commit | 282b1cd26649d69de038111f5876853df6ddc345 (patch) | |
| tree | 4287eb355348f3c47d0b2075ecc514c4f3361052 /src/math | |
| parent | 729fef0a9358e2f6f1cd8c75a1a0f7ee48b08c95 (diff) | |
| download | musl-282b1cd26649d69de038111f5876853df6ddc345.tar.gz musl-282b1cd26649d69de038111f5876853df6ddc345.tar.bz2 musl-282b1cd26649d69de038111f5876853df6ddc345.tar.xz musl-282b1cd26649d69de038111f5876853df6ddc345.zip | |
fix fmaf wrong result
if double precision r=x*y+z is not a half way case between two single
precision floats or it is an exact result then fmaf returns (float)r.
however the exactness check was wrong when |x*y| < |z| and could cause
incorrectly rounded result in nearest rounding mode when r is a half
way case.
fmaf(-0x1.26524ep-54, -0x1.cb7868p+11, 0x1.d10f5ep-29)
was incorrectly rounded up to 0x1.d117ap-29 instead of 0x1.d1179ep-29.
(exact result is 0x1.d1179efffffffecp-29, r is 0x1.d1179fp-29)
Diffstat (limited to 'src/math')
| -rw-r--r-- | src/math/fmaf.c | 2 |
1 files changed, 1 insertions, 1 deletions
diff --git a/src/math/fmaf.c b/src/math/fmaf.c index aa57feb6..80f5cd8a 100644 --- a/src/math/fmaf.c +++ b/src/math/fmaf.c @@ -50,7 +50,7 @@ float fmaf(float x, float y, float z) /* Common case: The double precision result is fine. */ if ((u.i & 0x1fffffff) != 0x10000000 || /* not a halfway case */ e == 0x7ff || /* NaN */ - result - xy == z || /* exact */ + (result - xy == z && result - z == xy) || /* exact */ fegetround() != FE_TONEAREST) /* not round-to-nearest */ { /* |