X-Git-Url: http://nsz.repo.hu/git/?p=musl;a=blobdiff_plain;f=src%2Fmath%2Ffmal.c;h=87e30fcf65394f71b34027081feb9cab0f013ca4;hp=cbaf46ebec926516f3abd6afa7fd76212ff6a4aa;hb=1e5eb73545ca6cfe8b918798835aaf6e07af5beb;hpb=d09a83f613c1d06442ed920ec55a0e5eedacb422;ds=sidebyside diff --git a/src/math/fmal.c b/src/math/fmal.c index cbaf46eb..87e30fcf 100644 --- a/src/math/fmal.c +++ b/src/math/fmal.c @@ -115,7 +115,7 @@ static inline long double add_and_denormalize(long double a, long double b, int if (bits_lost != 1 ^ (int)(u.bits.manl & 1)) sum.hi = nextafterl(sum.hi, INFINITY * sum.lo); } - return (ldexp(sum.hi, scale)); + return scalbnl(sum.hi, scale); } /* @@ -162,6 +162,7 @@ static inline struct dd dd_mul(long double a, long double b) */ long double fmal(long double x, long double y, long double z) { + #pragma STDC FENV_ACCESS ON long double xs, ys, zs, adj; struct dd xy, r; int oround; @@ -228,7 +229,7 @@ long double fmal(long double x, long double y, long double z) } } if (spread <= LDBL_MANT_DIG * 2) - zs = ldexpl(zs, -spread); + zs = scalbnl(zs, -spread); else zs = copysignl(LDBL_MIN, zs); @@ -254,23 +255,25 @@ long double fmal(long double x, long double y, long double z) */ fesetround(oround); volatile long double vzs = zs; /* XXX gcc CSE bug workaround */ - return (xy.hi + vzs + ldexpl(xy.lo, spread)); + return xy.hi + vzs + scalbnl(xy.lo, spread); } if (oround != FE_TONEAREST) { /* * There is no need to worry about double rounding in directed * rounding modes. + * TODO: underflow is not raised correctly, example in downward rounding: + * fmal(0x1.0000000001p-16000L, 0x1.0000000001p-400L, -0x1p-16440L) */ fesetround(oround); adj = r.lo + xy.lo; - return (ldexpl(r.hi + adj, spread)); + return scalbnl(r.hi + adj, spread); } adj = add_adjusted(r.lo, xy.lo); if (spread + ilogbl(r.hi) > -16383) - return (ldexpl(r.hi + adj, spread)); + return scalbnl(r.hi + adj, spread); else - return (add_and_denormalize(r.hi, adj, spread)); + return add_and_denormalize(r.hi, adj, spread); } #endif