bug in frexpl module

[Bruno Haible]

Hi Paolo,

frexpl of 1.0L should return the exponent 1 and mantissa 0.5L. Glibc does
this. The substitute code in frexpl.c, however, returns the exponent 0 and
mantissa 1.0L, which IMO is wrong.

OK to fix this?

By the way, now that gnulib supports unit tests, how about moving the test
code at the end of this file into a test module (modules/frexpl-tests,
tests/test-frexpl.c)?

Also, assuming that the exponent has at most 20 bits is a bit risky. One
can very well imagine 'long double' floating-point formats with an 31-bit
or 32-bit exponent. Care to change the size of the array to 64 instead of 20?
That should be safe for the near future.


2007-02-24  Bruno Haible  <address@hidden>

        * lib/frexpl.c (frexpl): Correct return values for x = 1.0L.

--- lib/frexpl.c        18 Feb 2007 15:10:28 -0000      1.5
+++ lib/frexpl.c        24 Feb 2007 15:21:51 -0000
@@ -34,7 +34,7 @@
   int exponent, bit;
 
   /* Check for zero, nan and infinity. */
-  if (x != x || x + x == x )
+  if (x != x || x + x == x)
     {
       *exp = 0;
       return x;
@@ -44,7 +44,7 @@
     return -frexpl(-x, exp);
 
   exponent = 0;
-  if (x > 1.0)
+  if (x >= 1.0)
     {
       for (next = exponents, exponents[0] = 2.0L, bit = 1;
           *next <= x + x;
@@ -90,6 +90,7 @@
   x = frexpl(-1.0L / 0.0L, &y); printf ("%.6Lg %d\n", x, y);
   x = frexpl(0.5L, &y); printf ("%.6Lg %d\n", x, y);
   x = frexpl(0.75L, &y); printf ("%.6Lg %d\n", x, y);
+  x = frexpl(1.0L, &y); printf ("%.6Lg %d\n", x, y);
   x = frexpl(3.6L, &y); printf ("%.6Lg %d\n", x, y);
   x = frexpl(17.8L, &y); printf ("%.6Lg %d\n", x, y);
   x = frexpl(8.0L, &y); printf ("%.6Lg %d\n", x, y);