https://gcc.gnu.org/bugzilla/show_bug.cgi?id=96600

--- Comment #5 from Jakub Jelinek <jakub at gcc dot gnu.org> ---
Actually, I think this is a glibc bug rather than gcc, because the values
computed at compile time look right, rather than those at runtime.
Consider:
int e = 69;
int main() {
  long double a =
-__builtin_ldexpl(1.9446689187403240306919491832695730985733566864714824565497322973045558e+00l,
e);
  long double b =
-__builtin_ldexpl(1.9446689187403240306919491832695730985733566864714824565497322973045558e+00l,
69);
  long double tmp1, tmp2;
  long double aa = __builtin_modfl(a, &tmp1);
  long double bb = __builtin_modfl(b, &tmp2);
  long double c = __builtin_truncl(a);
  long double d = __builtin_truncl(b);
  __builtin_printf ("a %.34La %.16a %.16a\n", a, (double) a, (double) (a -
(double) a));
  __builtin_printf ("b %.34La %.16a %.16a\n", b, (double) b, (double) (b -
(double) b));
  __builtin_printf ("aa %.34La %.16a %.16a\n", aa, (double) aa, (double) (aa -
(double) aa));
  __builtin_printf ("bb %.34La %.16a %.16a\n", bb, (double) bb, (double) (bb -
(double) bb));
  __builtin_printf ("tmp1 %.34La %.16a %.16a\n", tmp1, (double) tmp1, (double)
(tmp1 - (double) tmp1));
  __builtin_printf ("tmp2 %.34La %.16a %.16a\n", tmp2, (double) tmp2, (double)
(tmp2 - (double) tmp2));
  __builtin_printf ("c %.34La %.16a %.16a\n", c, (double) c, (double) (c -
(double) c));
  __builtin_printf ("d %.34La %.16a %.16a\n", d, (double) d, (double) (d -
(double) d));
  return 0;
}

The tmp1, tmp2, c and d values should be the same, i.e. the integral part.
But only tmp2, c and d are the same,
-0x1.f1d5d27f89914ab9200000000000000000p+69
while tmp1 is -0x1.f1d5d27f89914ab9200000000000000000p+69.
c is truncl computed at runtime, and tmp1 is the integral part from modfl
computed at runtime, so they should be the same, but they are not.
Seems glibc for modfl in this case just uses the upper double mantissa bits and
ignores the lower double mantissa bits.
Please file against glibc.

truncl has:
      hi = trunc (xh);
      if (hi != xh)
        {
          /* The high part is not an integer; the low part does not
             affect the result.  */
          xh = hi;
          xl = 0;
        }
      else
        {
          /* The high part is a nonzero integer.  */
          lo = xh > 0 ? floor (xl) : ceil (xl);
          xh = hi;
          xl = lo;
          ldbl_canonicalize_int (&xh, &xl);
        }
i.e. if the double trunc of the highpart double is equal to the highpart double
(this case), it does another floor or ceil.
While the modfl code assumes that if the exponent is > 52 and smaller than 103,
the fractional part is the second double, which is wrong.

The (i1&0x7fffffffffffffffLL) in there is weird because there is
i1 &= 0x000fffffffffffffLL; earlier, so it could just test i1, but more
importantly, the } else { ... code assumes that the low double must have a
particular exponent (it could have any smaller than the higher exponent minus
50something).

Reply via email to