Import of the latest gdtoa sources, which include fixes for minorvendor/gdtoa/20050923

problems relating to NaNs and rounding.

1 files changed, 36 insertions, 25 deletions

diff --git a/contrib/gdtoa/strtodg.c b/contrib/gdtoa/strtodg.c
index bf30e33ecf57..cbdf4aadf299 100644
--- a/contrib/gdtoa/strtodg.c
+++ b/contrib/gdtoa/strtodg.c
@@ -26,10 +26,8 @@ THIS SOFTWARE.
 
 ****************************************************************/
 
-/* Please send bug reports to
-	David M. Gay
-	dmg@acm.org
- */
+/* Please send bug reports to David M. Gay (dmg at acm dot org,
+ * with " at " changed at "@" and " dot " changed to ".").	*/
 
 #include "gdtoaimp.h"
 
@@ -239,7 +237,7 @@ rvOK
 			inex = STRTOG_Inexhi;
 			b = increment(b);
 			if ( (j = nb & kmask) !=0)
-				j = 32 - j;
+				j = ULbits - j;
 			if (hi0bits(b->x[b->wds - 1]) != j) {
 				if (!lostbits)
 					lostbits = b->x[0] & 1;
@@ -325,8 +323,8 @@ strtodg
 #endif
 {
 	int abe, abits, asub;
-	int bb0, bb2, bb5, bbe, bd2, bd5, bbbits, bs2;
-	int c, denorm, dsign, e, e1, e2, emin, esign, finished, i, inex, irv;
+	int bb0, bb2, bb5, bbe, bd2, bd5, bbbits, bs2, c, decpt, denorm;
+	int dsign, e, e1, e2, emin, esign, finished, i, inex, irv;
 	int j, k, nbits, nd, nd0, nf, nz, nz0, rd, rvbits, rve, rve1, sign;
 	int sudden_underflow;
 	CONST char *s, *s0, *s1;
@@ -385,7 +383,7 @@ strtodg
 	sudden_underflow = fpi->sudden_underflow;
 	s0 = s;
 	y = z = 0;
-	for(nd = nf = 0; (c = *s) >= '0' && c <= '9'; nd++, s++)
+	for(decpt = nd = nf = 0; (c = *s) >= '0' && c <= '9'; nd++, s++)
 		if (nd < 9)
 			y = 10*y + c - '0';
 		else if (nd < 16)
@@ -397,6 +395,7 @@ strtodg
 	if (c == '.')
 #endif
 		{
+		decpt = 1;
 		c = *++s;
 		if (!nd) {
 			for(; c == '0'; c = *++s)
@@ -471,7 +470,8 @@ strtodg
 		if (!nz && !nz0) {
 #ifdef INFNAN_CHECK
 			/* Check for Nan and Infinity */
-			switch(c) {
+			if (!decpt)
+			 switch(c) {
 			  case 'i':
 			  case 'I':
 				if (match(&s,"nf")) {
@@ -632,7 +632,14 @@ strtodg
 					dval(rv) *= tinytens[j];
 			}
 		}
-
+#ifdef IBM
+	/* e2 is a correction to the (base 2) exponent of the return
+	 * value, reflecting adjustments above to avoid overflow in the
+	 * native arithmetic.  For native IBM (base 16) arithmetic, we
+	 * must multiply e2 by 4 to change from base 16 to 2.
+	 */
+	e2 <<= 2;
+#endif
 	rvb = d2b(dval(rv), &rve, &rvbits);	/* rv = rvb * 2^rve */
 	rve += e2;
 	if ((j = rvbits - nbits) > 0) {
@@ -642,16 +649,8 @@ strtodg
 		}
 	bb0 = 0;	/* trailing zero bits in rvb */
 	e2 = rve + rvbits - nbits;
-	if (e2 > fpi->emax) {
-		rvb->wds = 0;
-		irv = STRTOG_Infinite | STRTOG_Overflow | STRTOG_Inexhi;
-#ifndef NO_ERRNO
-		errno = ERANGE;
-#endif
- infnanexp:
-		*exp = fpi->emax + 1;
-		goto ret;
-		}
+	if (e2 > fpi->emax + 1)
+		goto huge;
 	rve1 = rve + rvbits - nbits;
 	if (e2 < (emin = fpi->emin)) {
 		denorm = 1;
@@ -823,8 +822,8 @@ strtodg
 				break;
 			if (dsign) {
 				rvb = increment(rvb);
-				if ( (j = rvbits >> kshift) !=0)
-					j = 32 - j;
+				if ( (j = rvbits & kmask) !=0)
+					j = ULbits - j;
 				if (hi0bits(rvb->x[(rvb->wds - 1) >> kshift])
 						!= j)
 					rvbits++;
@@ -965,9 +964,11 @@ strtodg
 		Bfree(bs);
 		Bfree(delta);
 		}
-	if (!denorm && rvbits < nbits) {
-		j = nbits - rvbits;
-		rvb = lshift(rvb, j);
+	if (!denorm && (j = nbits - rvbits)) {
+		if (j > 0)
+			rvb = lshift(rvb, j);
+		else
+			rshift(rvb, -j);
 		rve -= j;
 		}
 	*exp = rve;
@@ -976,6 +977,16 @@ strtodg
 	Bfree(bs);
 	Bfree(bd0);
 	Bfree(delta);
+	if (rve > fpi->emax) {
+ huge:
+		rvb->wds = 0;
+		irv = STRTOG_Infinite | STRTOG_Overflow | STRTOG_Inexhi;
+#ifndef NO_ERRNO
+		errno = ERANGE;
+#endif
+ infnanexp:
+		*exp = fpi->emax + 1;
+		}
  ret:
 	if (denorm) {
 		if (sudden_underflow) {