Implement fmal, some fma bugfixes

1 files changed, 8 insertions, 0 deletions

diff --git a/sysdeps/ieee754/dbl-64/s_fma.c b/sysdeps/ieee754/dbl-64/s_fma.c
index 911682eceb..3b0bfd5ce6 100644
--- a/sysdeps/ieee754/dbl-64/s_fma.c
+++ b/sysdeps/ieee754/dbl-64/s_fma.c
@@ -43,6 +43,12 @@ __fma (double x, double y, double z)
       || __builtin_expect (u.ieee.exponent + v.ieee.exponent
 			   <= IEEE754_DOUBLE_BIAS + DBL_MANT_DIG, 0))
     {
+      /* If z is Inf, but x and y are finite, the result should be
+	 z rather than NaN.  */
+      if (w.ieee.exponent == 0x7ff
+	  && u.ieee.exponent != 0x7ff
+          && v.ieee.exponent != 0x7ff)
+	return (z + x) + y;
       /* If x or y or z is Inf/NaN, or if fma will certainly overflow,
 	 or if x * y is less than half of DBL_DENORM_MIN,
 	 compute as x * y + z.  */
@@ -165,6 +171,8 @@ __fma (double x, double y, double z)
     }
   else
     {
+      if ((u.ieee.mantissa1 & 1) == 0)
+	u.ieee.mantissa1 |= fetestexcept (FE_INEXACT) != 0;
       v.d = a1 + u.d;
       int j = fetestexcept (FE_INEXACT) != 0;
       feupdateenv (&env);