www.delorie.com/archives/browse.cgi   search  
Mail Archives: djgpp-workers/2013/07/10/13:40:43

X-Authentication-Warning: delorie.com: mail set sender to djgpp-workers-bounces using -f
X-Recipient: djgpp-workers AT delorie DOT com
Message-ID: <51DD9CCD.6070801@gmx.de>
Date: Wed, 10 Jul 2013 19:41:33 +0200
From: Juan Manuel Guerrero <juan DOT guerrero AT gmx DOT de>
User-Agent: Mozilla/5.0 (X11; Linux x86_64; rv:16.0) Gecko/20121025 Thunderbird/16.0.2
MIME-Version: 1.0
To: djgpp-workers AT delorie DOT com
Subject: Some more strto[d|f|ld] issues.
X-Provags-ID: V03:K0:eh7WGOih6ZhuL7O3cMXLCYZnNZKvI6SxyPEYWKUk7GazZcALwFh
vJ/EwSH3srvKA12Nw2I4MngDblScFW4u2tES7uwxf6gLDZr9wq6O/r2z16AAJH4b+0UpmqI
jjtNrCWDWBYGguEZdbw25whSAShaaasYvbb59l1j54hvWOGG0Kv/X4dxa8Z3QmSyZCBtuUi
D6Leyl31wGEg8MwM+Ix1w==
Reply-To: djgpp-workers AT delorie DOT com

Some years ago I implemented the hex string conversion method
for the strto[d|ld|f] functions.  In those days I decided to
implement rounding towards zero because it was the easiest and
because I did not know it better.  While I was porting Lua,
I noted that that was wrong and I have changed the rounding
to round half up.  This seems to be the same that glibc does
und it lets successfully pass the testsuite of Lua.  That was
not the case with the last port of Lua.

Regards,

Juan M. Guerrero



cvs ci -m"Round half up as glibc does instead of round towards zero." djgpp/src/libc/ansi/stdlib/strtod.c
cvs ci -m"Round half up as glibc does instead of round towards zero." djgpp/src/libc/ansi/stdlib/strtold.c
cvs ci -m"Round half up as glibc does instead of round towards zero." djgpp/src/libc/c99/stdlib/strtof.c


Logging in to :pserver:anonymous AT cvs DOT delorie DOT com:2401/cvs/djgpp
Index: djgpp/src/libc/ansi/stdlib/strtod.c
===================================================================
RCS file: /cvs/djgpp/djgpp/src/libc/ansi/stdlib/strtod.c,v
retrieving revision 1.16
diff -p -U 5 -r1.16 strtod.c
--- djgpp/src/libc/ansi/stdlib/strtod.c    9 Jul 2013 23:00:23 -0000    1.16
+++ djgpp/src/libc/ansi/stdlib/strtod.c    9 Jul 2013 23:25:27 -0000
@@ -14,14 +14,17 @@
  #include <ctype.h>
  #include <string.h>
  #include <libc/unconst.h>
  #include <libc/ieee.h>

+
+#define MANTISSA_SIZE       (52)     /*  Number binary digits in the fractional part of the mantissa.  */
  #define HEX_DIGIT_SIZE      (4)
  #define DOUBLE_BIAS         (0x3FFU)
  #define MAX_BIN_EXPONENT    (1023)   /*  Max. and min. binary exponent (inclusive) as  */
  #define MIN_BIN_EXPONENT    (-1022)  /*  defined in Intel manual (253665.pdf, Table 4.2).  */
+#define IS_ZERO_DIGIT(x)    ((x) == '0')
  #define IS_DEC_DIGIT(x)     (((x) >= '0') && ((x) <= '9'))
  #define IS_HEX_DIGIT(x)     ((((x) >= 'A') && ((x) <= 'F')) || \
                               (((x) >= 'a') && ((x) <= 'f')) || \
                               IS_DEC_DIGIT(x))
  #define IS_DEC_EXPONENT(x)  (((x[0]) == 'E' || (x[0]) == 'e') && \
@@ -120,13 +123,13 @@ strtod(const char *s, char **sret)
    }

    /* Handle 0xH.HHH[p|P][+|-]DDD. */
    if (!strnicmp("0x", s, 2) && (s[2] == '.' || IS_HEX_DIGIT(s[2])))
    {
-    int digits, integer_digits;
-    long int bin_exponent;
-    unsigned long long int mantissa;
+    const int max_digits = 1 + MANTISSA_SIZE / HEX_DIGIT_SIZE + 1; /* Two more digits than fit into mantissa.  */
+    int bin_exponent, digits, integer_digits;
+    unsigned long long int mantissa, msb_mask;
      _double_union_t ieee754;


      /*
       *  Mantissa.
@@ -134,11 +137,11 @@ strtod(const char *s, char **sret)
       */
      bin_exponent = 0;
      integer_digits = 0;
      mantissa = 0x00ULL;
      s += 2;  /*  Skip the hex prefix.  */
-    while (integer_digits < 16 && IS_HEX_DIGIT(*s))
+    while (integer_digits < max_digits && IS_HEX_DIGIT(*s))
      {
        flags = 1;
        mantissa <<= HEX_DIGIT_SIZE;
        mantissa |= IS_DEC_DIGIT(*s) ? *s - '0' :
                    ((*s >= 'A') && (*s <= 'F')) ? *s - 'A' + 10 : *s - 'a' + 10;
@@ -150,27 +153,30 @@ strtod(const char *s, char **sret)
      {
        /*
         *  Compute the binary exponent for a normalized mantissa by
         *  shifting the radix point inside the most significant hex digit.
         */
-      unsigned long long bit = 0x01ULL;

        for (digits = 0; IS_HEX_DIGIT(*s); s++)
          digits++;  /*  Counts hex digits.  */
+
+      msb_mask = 0x01ULL;
        bin_exponent = integer_digits * HEX_DIGIT_SIZE - 1;  /* 2**bin_exponent.  */
-      for (bit <<= bin_exponent; !(mantissa & bit); bin_exponent--)
-        bit >>= 1;
+      for (msb_mask <<= bin_exponent; !(mantissa & msb_mask); msb_mask >>= 1)
+        bin_exponent--;
        bin_exponent += digits * HEX_DIGIT_SIZE;
+      integer_digits += digits;
      }

+    digits = integer_digits;
      if (*s == radix_point)
      {
-      int fraction_zeros = 0;
+      int extra_shifts, fraction_zeros = 0;

        s++;
        digits = integer_digits;
-      while ((digits - fraction_zeros) < 16 && IS_HEX_DIGIT(*s))
+      while ((digits - fraction_zeros) < max_digits && IS_HEX_DIGIT(*s))
        {
          flags = 1;
          digits++;  /*  Counts hex digits.  */
          mantissa <<= HEX_DIGIT_SIZE;
          mantissa |= IS_DEC_DIGIT(*s) ? *s - '0' :
@@ -183,16 +189,18 @@ strtod(const char *s, char **sret)
        {
          /*
           *  Compute the binary exponent for a normalized mantissa by
           *  shifting the radix point inside the most significant hex digit.
           */
-        unsigned long long bit = 0x01ULL;

+        msb_mask = 0x01ULL;
          bin_exponent = -fraction_zeros * HEX_DIGIT_SIZE;  /* 2**bin_exponent.  */
-        for (bit <<= (digits * HEX_DIGIT_SIZE + bin_exponent); !(mantissa & bit); bin_exponent--)
-          bit >>= 1;
+        for (msb_mask <<= (digits * HEX_DIGIT_SIZE + bin_exponent); !(mantissa & msb_mask); msb_mask >>= 1)
+          bin_exponent--;
        }
+      else if ((extra_shifts = digits - integer_digits) > 0)
+        msb_mask <<= extra_shifts * HEX_DIGIT_SIZE;
      }

      if (flags == 0)
      {
        if (sret)
@@ -203,10 +211,28 @@ strtod(const char *s, char **sret)
        }
        errno = EINVAL;  /*  No valid mantissa, no conversion could be performed.  */
        return 0.0L;
      }

+    if (digits >= max_digits)
+    {
+      /*
+       *  Round half towards plus infinity (round half up).
+       */
+      const int lsd = 0x000000000000000FULL & mantissa;  /* Least significant hex digit.  Will be rounded out.  */
+      if (lsd > 0x07)
+      {
+        mantissa += 0x0000000000000010ULL;  /* Smallest float greater than x.  */
+        if (!(mantissa & msb_mask))
+        {
+          /*  Overflow.  */
+          mantissa >>= 1;
+          bin_exponent++;
+        }
+      }
+    }
+
      if (mantissa)
      {
        /*
         *  Normalize mantissa.
         */
Index: djgpp/src/libc/ansi/stdlib/strtold.c
===================================================================
RCS file: /cvs/djgpp/djgpp/src/libc/ansi/stdlib/strtold.c,v
retrieving revision 1.16
diff -p -U 5 -r1.16 strtold.c
--- djgpp/src/libc/ansi/stdlib/strtold.c    9 Jul 2013 23:00:25 -0000    1.16
+++ djgpp/src/libc/ansi/stdlib/strtold.c    9 Jul 2013 23:25:27 -0000
@@ -12,14 +12,17 @@
  #include <math.h>
  #include <string.h>
  #include <libc/unconst.h>
  #include <libc/ieee.h>

+
+#define MANTISSA_SIZE       (64)     /*  Number binary digits in the integer and fractional part of the mantissa.  */
  #define HEX_DIGIT_SIZE      (4)
  #define LONG_DOUBLE_BIAS    (0x3FFFU)
  #define MAX_BIN_EXPONENT    (16383)   /*  Max. and min. binary exponent (inclusive) as  */
  #define MIN_BIN_EXPONENT    (-16382)  /*  defined in Intel manual (253665.pdf, Table 4.2).  */
+#define IS_ZERO_DIGIT(x)    ((x) == '0')
  #define IS_DEC_DIGIT(x)     (((x) >= '0') && ((x) <= '9'))
  #define IS_HEX_DIGIT(x)     ((((x) >= 'A') && ((x) <= 'F')) || \
                               (((x) >= 'a') && ((x) <= 'f')) || \
                               IS_DEC_DIGIT(x))
  #define IS_DEC_EXPONENT(x)  (((x[0]) == 'E' || (x[0]) == 'e') && \
@@ -125,13 +128,13 @@ strtold(const char *s, char **sret)

    /* Handle 0xH.HHH[p|P][+|-]DDD. */
    if (!strnicmp("0x", s, 2) && (s[2] == '.' || IS_HEX_DIGIT(s[2])))
    {
      const char *next_char = NULL;
-    int digits, integer_digits;
-    long long int bin_exponent;
-    unsigned long long int mantissa;
+    const int max_digits = MANTISSA_SIZE / HEX_DIGIT_SIZE;  /* The exact number of digits that fits in mantissa.  */
+    int bin_exponent, digits, integer_digits;
+    unsigned long long int mantissa, msb_mask;
      _longdouble_union_t ieee754;


      /*
       *  Mantissa.
@@ -140,11 +143,11 @@ strtold(const char *s, char **sret)
       */
      bin_exponent = 0;
      integer_digits = 0;
      mantissa = 0x00ULL;
      s += 2;  /*  Skip the hex prefix.  */
-    while (integer_digits < 16 && IS_HEX_DIGIT(*s))
+    while (integer_digits < max_digits && IS_HEX_DIGIT(*s))
      {
        flags = 1;
        mantissa <<= HEX_DIGIT_SIZE;
        mantissa |= IS_DEC_DIGIT(*s) ? *s - '0' :
                    ((*s >= 'A') && (*s <= 'F')) ? *s - 'A' + 10 : *s - 'a' + 10;
@@ -156,28 +159,31 @@ strtold(const char *s, char **sret)
      {
        /*
         *  Compute the binary exponent for a normalized mantissa by
         *  shifting the radix point inside the most significant hex digit.
         */
-      unsigned long long bit = 0x01ULL;

        next_char = s;
        for (digits = 0; IS_HEX_DIGIT(*s); s++)
          digits++;  /*  Counts hex digits.  */
+
+      msb_mask = 0x01ULL;
        bin_exponent = integer_digits * HEX_DIGIT_SIZE - 1;  /* 2**bin_exponent.  */
-      for (bit <<= bin_exponent; !(mantissa & bit); bin_exponent--)
-        bit >>= 1;
+      for (msb_mask <<= bin_exponent; !(mantissa & msb_mask); msb_mask >>= 1)
+        bin_exponent--;
        bin_exponent += digits * HEX_DIGIT_SIZE;
+      integer_digits += digits;
      }

+    digits = integer_digits;
      if (*s == radix_point)
      {
-      int fraction_zeros = 0;
+      int extra_shifts, fraction_zeros = 0;

        s++;
        digits = integer_digits;
-      while ((digits - fraction_zeros) < 16 && IS_HEX_DIGIT(*s))
+      while ((digits - fraction_zeros) < max_digits && IS_HEX_DIGIT(*s))
        {
          flags = 1;
          digits++;  /*  Counts hex digits.  */
          mantissa <<= HEX_DIGIT_SIZE;
          mantissa |= IS_DEC_DIGIT(*s) ? *s - '0' :
@@ -191,16 +197,18 @@ strtold(const char *s, char **sret)
        {
          /*
           *  Compute the binary exponent for a normalized mantissa by
           *  shifting the radix point inside the most significant hex digit.
           */
-        unsigned long long bit = 0x01ULL;

+        msb_mask = 0x01ULL;
          bin_exponent = -fraction_zeros * HEX_DIGIT_SIZE;  /* 2**bin_exponent.  */
-        for (bit <<= (digits * HEX_DIGIT_SIZE + bin_exponent); !(mantissa & bit); bin_exponent--)
-          bit >>= 1;
+        for (msb_mask <<= (digits * HEX_DIGIT_SIZE + bin_exponent); !(mantissa & msb_mask); msb_mask >>= 1)
+          bin_exponent--;
        }
+      else if ((extra_shifts = digits - integer_digits) > 0)
+        msb_mask <<= extra_shifts * HEX_DIGIT_SIZE;
      }

      if (flags == 0)
      {
        if (sret)
@@ -211,38 +219,37 @@ strtold(const char *s, char **sret)
        }
        errno = EINVAL;  /*  No valid mantissa, no conversion could be performed.  */
        return 0.0L;
      }

-    if (mantissa)
+    if (digits >= max_digits)
      {
        /*
-       *  Normalize mantissa.
-       *  If there is still a valid hex character in the string
-       *  its bits will be inserted in the LSB of the mantissa,
-       *  else zeros will be inserted.
+       *  Round half towards plus infinity (round half up).
         */
-      for (digits = 0; !(mantissa & 0x8000000000000000ULL); digits++)
-        mantissa <<= 1;  /*  Shift a binary 1 into the integer part of the mantissa.  */
-      if (IS_HEX_DIGIT(*next_char))
+      const int lsd = IS_DEC_DIGIT(*next_char) ? *next_char - '0' :  /*  Least significant hex digit.  Will be rounded out.  */
+                      ((*next_char >= 'A') && (*next_char <= 'F')) ? *next_char - 'A' + 10 : *next_char - 'a' + 10;
+      if (lsd > 0x07)
        {
-        unsigned long long bits = IS_DEC_DIGIT(*next_char) ? *next_char - '0' :
-                                  ((*next_char >= 'A') && (*next_char <= 'F')) ? *next_char - 'A' + 10 : *next_char - 'a' + 10;
-
-        switch (digits)
+        mantissa += 0x0000000000000001ULL;  /* Smallest float greater than x.  */
+        if (!(mantissa & msb_mask))
          {
-        case 1:
-          mantissa |= (0x01ULL & bits >> 3);
-          break;
-        case 2:
-          mantissa |= (0x03ULL & bits >> 2);
-          break;
-        case 3:
-          mantissa |= (0x07ULL & bits >> 1);
-          break;
+          /*  Overflow.  */
+          mantissa >>= 1;
+          mantissa |= 0x8000000000000000ULL;
+          bin_exponent++;
          }
        }
+    }
+
+    if (mantissa)
+    {
+      /*
+       *  Normalize mantissa.
+       */
+      for (digits = 0; !(mantissa & 0x8000000000000000ULL); digits++)
+        mantissa <<= 1;  /*  Shift a binary 1 into the integer part of the mantissa.  */
        /*  At this point the mantissa is normalized and the exponent has been adjusted accordingly.  */
      }


      /*
Index: djgpp/src/libc/c99/stdlib/strtof.c
===================================================================
RCS file: /cvs/djgpp/djgpp/src/libc/c99/stdlib/strtof.c,v
retrieving revision 1.14
diff -p -U 5 -r1.14 strtof.c
--- djgpp/src/libc/c99/stdlib/strtof.c    9 Jul 2013 23:00:27 -0000    1.14
+++ djgpp/src/libc/c99/stdlib/strtof.c    9 Jul 2013 23:25:27 -0000
@@ -16,14 +16,17 @@
  #include <ctype.h>
  #include <string.h>
  #include <libc/unconst.h>
  #include <libc/ieee.h>

+
+#define MANTISSA_SIZE       (23)     /*  Number binary digits in the fractional part of the mantissa.  */
  #define HEX_DIGIT_SIZE      (4)
  #define FLOAT_BIAS          (0x7FU)
  #define MAX_BIN_EXPONENT    (127)   /*  Max. and min. binary exponent (inclusive) as  */
  #define MIN_BIN_EXPONENT    (-126)  /*  defined in Intel manual (253665.pdf, Table 4.2).  */
+#define IS_ZERO_DIGIT(x)    ((x) == '0')
  #define IS_DEC_DIGIT(x)     (((x) >= '0') && ((x) <= '9'))
  #define IS_HEX_DIGIT(x)     ((((x) >= 'A') && ((x) <= 'F')) || \
                               (((x) >= 'a') && ((x) <= 'f')) || \
                               IS_DEC_DIGIT(x))
  #define IS_DEC_EXPONENT(x)  (((x[0]) == 'E' || (x[0]) == 'e') && \
@@ -120,13 +123,13 @@ strtof(const char *s, char **sret)
    }

    /* Handle 0xH.HHH[p|P][+|-]DDD. */
    if (!strnicmp("0x", s, 2) && (s[2] == '.' || IS_HEX_DIGIT(s[2])))
    {
-    int digits, integer_digits;
-    int bin_exponent;
-    unsigned long int mantissa;
+    const int max_digits = 1 + MANTISSA_SIZE / HEX_DIGIT_SIZE + 1; /* One more digit than fits into mantissa.  */
+    int bin_exponent, digits, integer_digits;
+    unsigned long int mantissa, msb_mask;
      _float_union_t ieee754;


      /*
       *  Mantissa.
@@ -135,11 +138,11 @@ strtof(const char *s, char **sret)
       */
      bin_exponent = 0;
      integer_digits = 0;
      mantissa = 0x00UL;
      s += 2;  /*  Skip the hex prefix.  */
-    while (integer_digits < 8 && IS_HEX_DIGIT(*s))
+    while (integer_digits < max_digits && IS_HEX_DIGIT(*s))
      {
        flags = 1;
        mantissa <<= HEX_DIGIT_SIZE;
        mantissa |= IS_DEC_DIGIT(*s) ? *s - '0' :
                    ((*s >= 'A') && (*s <= 'F')) ? *s - 'A' + 10 : *s - 'a' + 10;
@@ -151,27 +154,29 @@ strtof(const char *s, char **sret)
      {
        /*
         *  Compute the binary exponent for a normalized mantissa by
         *  shifting the radix point inside the most significant hex digit.
         */
-      unsigned long bit = 0x01UL;

        for (digits = 0; IS_HEX_DIGIT(*s); s++)
          digits++;  /*  Counts hex digits.  */
+
+      msb_mask = 0x01UL;
        bin_exponent = integer_digits * HEX_DIGIT_SIZE - 1;  /* 2**bin_exponent.  */
-      for (bit <<= bin_exponent; !(mantissa & bit); bin_exponent--)
-        bit >>= 1;
+      for (msb_mask <<= bin_exponent; !(mantissa & msb_mask); msb_mask >>= 1)
+        bin_exponent--;
        bin_exponent += digits * HEX_DIGIT_SIZE;
+      integer_digits += digits;
      }

+    digits = integer_digits;
      if (*s == radix_point)
      {
-      int fraction_zeros = 0;
+      int extra_shifts, fraction_zeros = 0;

        s++;
-      digits = integer_digits;
-      while ((digits - fraction_zeros) < 8 && IS_HEX_DIGIT(*s))
+      while ((digits - fraction_zeros) < max_digits && IS_HEX_DIGIT(*s))
        {
          flags = 1;
          digits++;  /*  Counts hex digits.  */
          mantissa <<= HEX_DIGIT_SIZE;
          mantissa |= IS_DEC_DIGIT(*s) ? *s - '0' :
@@ -184,15 +189,35 @@ strtof(const char *s, char **sret)
        {
          /*
           *  Compute the binary exponent for a normalized mantissa by
           *  shifting the radix point inside the most significant hex digit.
           */
-        unsigned long bit = 0x01UL;

+        msb_mask = 0x01UL;
          bin_exponent = -fraction_zeros * HEX_DIGIT_SIZE;  /* 2**bin_exponent.  */
-        for (bit <<= (digits * HEX_DIGIT_SIZE + bin_exponent); !(mantissa & bit); bin_exponent--)
-          bit >>= 1;
+        for (msb_mask <<= (digits * HEX_DIGIT_SIZE + bin_exponent); !(mantissa & msb_mask); msb_mask >>= 1)
+          bin_exponent--;
+      }
+      else if ((extra_shifts = digits - integer_digits) > 0)
+        msb_mask <<= extra_shifts * HEX_DIGIT_SIZE;
+    }
+
+    if (digits >= max_digits)
+    {
+      /*
+       *  Round half towards plus infinity (round half up).
+       */
+      const int lsd = 0x0000000FUL & mantissa;  /*  Least significant hex digit.  Will be rounded out.  */
+      if (lsd > 0x07)
+      {
+        mantissa += 0x00000010UL;  /* Smallest float greater than x.  */
+        if (!(mantissa & msb_mask))
+        {
+          /*  Overflow.  */
+          mantissa >>= 1;
+          bin_exponent++;
+        }
        }
      }

      if (flags == 0)
      {

- Raw text -


  webmaster     delorie software   privacy  
  Copyright © 2019   by DJ Delorie     Updated Jul 2019