Skip to content

HTTPS clone URL

Subversion checkout URL

You can clone with
or
.
Download ZIP
Browse files

from www.netlib.org

  • Loading branch information...
commit 2d7b777eb320427dd91c754d826789f616c9b966 1 parent d807d49
christos authored
Showing with 19,873 additions and 910 deletions.
  1. +29 −3 lib/libc/gdtoa/README
  2. +0 −2  lib/libc/gdtoa/arithchk.c
  3. +0 −2  lib/libc/gdtoa/dmisc.c
  4. +117 −92 lib/libc/gdtoa/dtoa.c
  5. +10 −7 lib/libc/gdtoa/g_Qfmt.c
  6. +59 −13 lib/libc/gdtoa/g__fmt.c
  7. +29 −15 lib/libc/gdtoa/g_ddfmt.c
  8. +15 −11 lib/libc/gdtoa/g_dfmt.c
  9. +10 −7 lib/libc/gdtoa/g_ffmt.c
  10. +10 −7 lib/libc/gdtoa/g_xLfmt.c
  11. +10 −7 lib/libc/gdtoa/g_xfmt.c
  12. +95 −91 lib/libc/gdtoa/gdtoa.c
  13. +11 −12 lib/libc/gdtoa/gdtoa.h
  14. +18 −0 lib/libc/gdtoa/gdtoa_fltrnds.h
  15. +49 −40 lib/libc/gdtoa/gdtoaimp.h
  16. +124 −23 lib/libc/gdtoa/gethex.c
  17. +0 −2  lib/libc/gdtoa/gmisc.c
  18. +0 −2  lib/libc/gdtoa/hd_init.c
  19. +21 −4 lib/libc/gdtoa/hexnan.c
  20. +27 −18 lib/libc/gdtoa/makefile
  21. +46 −38 lib/libc/gdtoa/misc.c
  22. +8 −0 lib/libc/gdtoa/printf.c
  23. +1,669 −0 lib/libc/gdtoa/printf.c0
  24. +0 −2  lib/libc/gdtoa/qnan.c
  25. +14 −16 lib/libc/gdtoa/smisc.c
  26. +103 −0 lib/libc/gdtoa/stdio1.h
  27. +0 −2  lib/libc/gdtoa/strtoIQ.c
  28. +0 −2  lib/libc/gdtoa/strtoId.c
  29. +0 −2  lib/libc/gdtoa/strtoIdd.c
  30. +0 −2  lib/libc/gdtoa/strtoIf.c
  31. +11 −9 lib/libc/gdtoa/strtoIg.c
  32. +0 −2  lib/libc/gdtoa/strtoIx.c
  33. +0 −2  lib/libc/gdtoa/strtoIxL.c
  34. +309 −219 lib/libc/gdtoa/strtod.c
  35. +42 −48 lib/libc/gdtoa/strtodI.c
  36. +128 −77 lib/libc/gdtoa/strtodg.c
  37. +0 −2  lib/libc/gdtoa/strtodnrp.c
  38. +8 −5 lib/libc/gdtoa/strtof.c
  39. +7 −4 lib/libc/gdtoa/strtopQ.c
  40. +6 −3 lib/libc/gdtoa/strtopd.c
  41. +26 −23 lib/libc/gdtoa/strtopdd.c
  42. +8 −5 lib/libc/gdtoa/strtopf.c
  43. +9 −5 lib/libc/gdtoa/strtopx.c
  44. +9 −5 lib/libc/gdtoa/strtopxL.c
  45. +0 −2  lib/libc/gdtoa/strtorQ.c
  46. +0 −2  lib/libc/gdtoa/strtord.c
  47. +17 −19 lib/libc/gdtoa/strtordd.c
  48. +1 −3 lib/libc/gdtoa/strtorf.c
  49. +2 −3 lib/libc/gdtoa/strtorx.c
  50. +2 −3 lib/libc/gdtoa/strtorxL.c
  51. +0 −2  lib/libc/gdtoa/sum.c
  52. +1,180 −0 lib/libc/gdtoa/test/Q.ou0
  53. +1,468 −0 lib/libc/gdtoa/test/Q.ou1
  54. +178 −0 lib/libc/gdtoa/test/Qtest.c
  55. +81 −0 lib/libc/gdtoa/test/README
  56. +986 −0 lib/libc/gdtoa/test/d.out
  57. +168 −0 lib/libc/gdtoa/test/dI.out
  58. +168 −0 lib/libc/gdtoa/test/dIsi.out
  59. +82 −0 lib/libc/gdtoa/test/dItest.c
  60. +1,356 −0 lib/libc/gdtoa/test/dd.out
  61. +376 −0 lib/libc/gdtoa/test/ddsi.out
  62. +177 −0 lib/libc/gdtoa/test/ddtest.c
  63. +291 −0 lib/libc/gdtoa/test/dt.c
  64. +139 −0 lib/libc/gdtoa/test/dtest.c
  65. +684 −0 lib/libc/gdtoa/test/dtst.out
  66. +988 −0 lib/libc/gdtoa/test/f.out
  67. +140 −0 lib/libc/gdtoa/test/ftest.c
  68. +90 −0 lib/libc/gdtoa/test/getround.c
  69. +172 −0 lib/libc/gdtoa/test/makefile
  70. +6 −0 lib/libc/gdtoa/test/obad/strtodt.out
  71. +1,460 −0 lib/libc/gdtoa/test/obad/xL.out
  72. +13 −0 lib/libc/gdtoa/test/pfLqtestnos
  73. +158 −0 lib/libc/gdtoa/test/pftest.c
  74. +13 −0 lib/libc/gdtoa/test/pftestLq.out
  75. +15 −0 lib/libc/gdtoa/test/pftestQ.out
  76. +15 −0 lib/libc/gdtoa/test/pftestnos
  77. +15 −0 lib/libc/gdtoa/test/pftestx.out
  78. +76 −0 lib/libc/gdtoa/test/rtestnos
  79. +2 −0  lib/libc/gdtoa/test/strtoIdSI.c
  80. +2 −0  lib/libc/gdtoa/test/strtoIddSI.c
  81. +2 −0  lib/libc/gdtoa/test/strtodISI.c
  82. +138 −0 lib/libc/gdtoa/test/strtodt.c
  83. +2 −0  lib/libc/gdtoa/test/strtopddSI.c
  84. +2 −0  lib/libc/gdtoa/test/strtorddSI.c
  85. +28 −0 lib/libc/gdtoa/test/testnos
  86. +20 −0 lib/libc/gdtoa/test/testnos1
  87. +352 −0 lib/libc/gdtoa/test/testnos3
  88. +1,182 −0 lib/libc/gdtoa/test/x.ou0
  89. +1,460 −0 lib/libc/gdtoa/test/x.ou1
  90. +1,182 −0 lib/libc/gdtoa/test/xL.ou0
  91. +1,462 −0 lib/libc/gdtoa/test/xL.ou1
  92. +169 −0 lib/libc/gdtoa/test/xLtest.c
  93. +58 −0 lib/libc/gdtoa/test/xQtest.c
  94. +41 −0 lib/libc/gdtoa/test/xsum0.out
  95. +170 −0 lib/libc/gdtoa/test/xtest.c
  96. +10 −12 lib/libc/gdtoa/ulp.c
  97. +37 −33 lib/libc/gdtoa/xsum0.out
View
32 lib/libc/gdtoa/README
@@ -1,5 +1,3 @@
-$NetBSD: README,v 1.1.1.1 2006/01/25 15:18:40 kleink Exp $
-
This directory contains source for a library of binary -> decimal
and decimal -> binary conversion routines, for single-, double-,
and extended-precision IEEE binary floating-point arithmetic, and
@@ -58,7 +56,9 @@ two letters:
whose sum is the desired value
For decimal -> binary conversions, there are three families of
-helper routines: one for round-nearest:
+helper routines: one for round-nearest (or the current rounding
+mode on IEEE-arithmetic systems that provide the C99 fegetround()
+function, if compiled with -DHonor_FLT_ROUNDS):
strtof
strtod
@@ -193,6 +193,9 @@ in the buffer, if the buffer was long enough, or 0. Other forms of
conversion are easily done with the help of gdtoa(), such as %e or %f
style and conversions with direction of rounding specified (so that, if
desired, the decimal value is either >= or <= the binary value).
+On IEEE-arithmetic systems that provide the C99 fegetround() function,
+if compiled with -DHonor_FLT_ROUNDS, these routines honor the current
+rounding mode.
For an example of more general conversions based on dtoa(), see
netlib's "printf.c from ampl/solvers".
@@ -334,5 +337,28 @@ Compiling g__fmt.c, strtod.c, and strtodg.c with -DUSE_LOCALE causes
the decimal-point character to be taken from the current locale; otherwise
it is '.'.
+Source files dtoa.c and strtod.c in this directory are derived from
+netlib's "dtoa.c from fp" and are meant to function equivalently.
+When compiled with Honor_FLT_ROUNDS #defined (on systems that provide
+FLT_ROUNDS and fegetround() as specified in the C99 standard), they
+honor the current rounding mode. Because FLT_ROUNDS is buggy on some
+(Linux) systems -- not reflecting calls on fesetround(), as the C99
+standard says it should -- when Honor_FLT_ROUNDS is #defined, the
+current rounding mode is obtained from fegetround() rather than from
+FLT_ROUNDS, unless Trust_FLT_ROUNDS is also #defined.
+
+Compile with -DUSE_LOCALE to use the current locale; otherwise
+decimal points are assumed to be '.'. With -DUSE_LOCALE, unless
+you also compile with -DNO_LOCALE_CACHE, the details about the
+current "decimal point" character string are cached and assumed not
+to change during the program's execution.
+
+On machines with a 64-bit long double and perhaps a 113-bit "quad"
+type, you can invoke "make Printf" to add Printf (and variants, such
+as Fprintf) to gdtoa.a. These are analogs, declared in stdio1.h, of
+printf and fprintf, etc. in which %La, %Le, %Lf, and %Lg are for long
+double and (if appropriate) %Lqa, %Lqe, %Lqf, and %Lqg are for quad
+precision printing.
+
Please send comments to David M. Gay (dmg at acm dot org, with " at "
changed at "@" and " dot " changed to ".").
View
2  lib/libc/gdtoa/arithchk.c
@@ -1,5 +1,3 @@
-/* $NetBSD: arithchk.c,v 1.1.1.1 2006/01/25 15:18:40 kleink Exp $ */
-
/****************************************************************
Copyright (C) 1997, 1998 Lucent Technologies
All Rights Reserved
View
2  lib/libc/gdtoa/dmisc.c
@@ -1,5 +1,3 @@
-/* $NetBSD: dmisc.c,v 1.1.1.1 2006/01/25 15:18:40 kleink Exp $ */
-
/****************************************************************
The author of this software is David M. Gay.
View
209 lib/libc/gdtoa/dtoa.c
@@ -1,5 +1,3 @@
-/* $NetBSD: dtoa.c,v 1.1.1.1 2006/01/25 15:18:41 kleink Exp $ */
-
/****************************************************************
The author of this software is David M. Gay.
@@ -68,7 +66,6 @@ THIS SOFTWARE.
*/
#ifdef Honor_FLT_ROUNDS
-#define Rounding rounding
#undef Check_FLT_ROUNDS
#define Check_FLT_ROUNDS
#else
@@ -78,10 +75,10 @@ THIS SOFTWARE.
char *
dtoa
#ifdef KR_headers
- (d, mode, ndigits, decpt, sign, rve)
- double d; int mode, ndigits, *decpt, *sign; char **rve;
+ (d0, mode, ndigits, decpt, sign, rve)
+ double d0; int mode, ndigits, *decpt, *sign; char **rve;
#else
- (double d, int mode, int ndigits, int *decpt, int *sign, char **rve)
+ (double d0, int mode, int ndigits, int *decpt, int *sign, char **rve)
#endif
{
/* Arguments ndigits, decpt, sign are similar to those
@@ -127,14 +124,25 @@ dtoa
ULong x;
#endif
Bigint *b, *b1, *delta, *mlo, *mhi, *S;
- double d2, ds, eps;
+ U d, d2, eps;
+ double ds;
char *s, *s0;
-#ifdef Honor_FLT_ROUNDS
- int rounding;
-#endif
#ifdef SET_INEXACT
int inexact, oldinexact;
#endif
+#ifdef Honor_FLT_ROUNDS /*{*/
+ int Rounding;
+#ifdef Trust_FLT_ROUNDS /*{{ only define this if FLT_ROUNDS really works! */
+ Rounding = Flt_Rounds;
+#else /*}{*/
+ Rounding = 1;
+ switch(fegetround()) {
+ case FE_TOWARDZERO: Rounding = 0; break;
+ case FE_UPWARD: Rounding = 2; break;
+ case FE_DOWNWARD: Rounding = 3;
+ }
+#endif /*}}*/
+#endif /*}*/
#ifndef MULTIPLE_THREADS
if (dtoa_result) {
@@ -142,35 +150,35 @@ dtoa
dtoa_result = 0;
}
#endif
-
- if (word0(d) & Sign_bit) {
+ d.d = d0;
+ if (word0(&d) & Sign_bit) {
/* set sign for everything, including 0's and NaNs */
*sign = 1;
- word0(d) &= ~Sign_bit; /* clear sign bit */
+ word0(&d) &= ~Sign_bit; /* clear sign bit */
}
else
*sign = 0;
#if defined(IEEE_Arith) + defined(VAX)
#ifdef IEEE_Arith
- if ((word0(d) & Exp_mask) == Exp_mask)
+ if ((word0(&d) & Exp_mask) == Exp_mask)
#else
- if (word0(d) == 0x8000)
+ if (word0(&d) == 0x8000)
#endif
{
/* Infinity or NaN */
*decpt = 9999;
#ifdef IEEE_Arith
- if (!word1(d) && !(word0(d) & 0xfffff))
+ if (!word1(&d) && !(word0(&d) & 0xfffff))
return nrv_alloc("Infinity", rve, 8);
#endif
return nrv_alloc("NaN", rve, 3);
}
#endif
#ifdef IBM
- dval(d) += 0; /* normalize */
+ dval(&d) += 0; /* normalize */
#endif
- if (!dval(d)) {
+ if (!dval(&d)) {
*decpt = 1;
return nrv_alloc("0", rve, 1);
}
@@ -180,35 +188,35 @@ dtoa
inexact = 1;
#endif
#ifdef Honor_FLT_ROUNDS
- if ((rounding = Flt_Rounds) >= 2) {
+ if (Rounding >= 2) {
if (*sign)
- rounding = rounding == 2 ? 0 : 2;
+ Rounding = Rounding == 2 ? 0 : 2;
else
- if (rounding != 2)
- rounding = 0;
+ if (Rounding != 2)
+ Rounding = 0;
}
#endif
- b = d2b(dval(d), &be, &bbits);
+ b = d2b(dval(&d), &be, &bbits);
#ifdef Sudden_Underflow
- i = (int)(word0(d) >> Exp_shift1 & (Exp_mask>>Exp_shift1));
+ i = (int)(word0(&d) >> Exp_shift1 & (Exp_mask>>Exp_shift1));
#else
- if (( i = (int)(word0(d) >> Exp_shift1 & (Exp_mask>>Exp_shift1)) )!=0) {
+ if (( i = (int)(word0(&d) >> Exp_shift1 & (Exp_mask>>Exp_shift1)) )!=0) {
#endif
- dval(d2) = dval(d);
- word0(d2) &= Frac_mask1;
- word0(d2) |= Exp_11;
+ dval(&d2) = dval(&d);
+ word0(&d2) &= Frac_mask1;
+ word0(&d2) |= Exp_11;
#ifdef IBM
- if (( j = 11 - hi0bits(word0(d2) & Frac_mask) )!=0)
- dval(d2) /= 1 << j;
+ if (( j = 11 - hi0bits(word0(&d2) & Frac_mask) )!=0)
+ dval(&d2) /= 1 << j;
#endif
/* log(x) ~=~ log(1.5) + (x-1.5)/1.5
* log10(x) = log(x) / log(10)
* ~=~ log(1.5)/log(10) + (x-1.5)/(1.5*log(10))
- * log10(d) = (i-Bias)*log(2)/log(10) + log10(d2)
+ * log10(&d) = (i-Bias)*log(2)/log(10) + log10(&d2)
*
- * This suggests computing an approximation k to log10(d) by
+ * This suggests computing an approximation k to log10(&d) by
*
* k = (i - Bias)*0.301029995663981
* + ( (d2-1.5)*0.289529654602168 + 0.176091259055681 );
@@ -237,21 +245,21 @@ dtoa
/* d is denormalized */
i = bbits + be + (Bias + (P-1) - 1);
- x = i > 32 ? word0(d) << 64 - i | word1(d) >> i - 32
- : word1(d) << 32 - i;
- dval(d2) = x;
- word0(d2) -= 31*Exp_msk1; /* adjust exponent */
+ x = i > 32 ? word0(&d) << (64 - i) | word1(&d) >> (i - 32)
+ : word1(&d) << (32 - i);
+ dval(&d2) = x;
+ word0(&d2) -= 31*Exp_msk1; /* adjust exponent */
i -= (Bias + (P-1) - 1) + 1;
denorm = 1;
}
#endif
- ds = (dval(d2)-1.5)*0.289529654602168 + 0.1760912590558 + i*0.301029995663981;
+ ds = (dval(&d2)-1.5)*0.289529654602168 + 0.1760912590558 + i*0.301029995663981;
k = (int)ds;
if (ds < 0. && ds != k)
k--; /* want k = floor(ds) */
k_check = 1;
if (k >= 0 && k <= Ten_pmax) {
- if (dval(d) < tens[k])
+ if (dval(&d) < tens[k])
k--;
k_check = 0;
}
@@ -290,10 +298,11 @@ dtoa
try_quick = 0;
}
leftright = 1;
+ ilim = ilim1 = -1; /* Values for cases 0 and 1; done here to */
+ /* silence erroneous "gcc -Wall" warning. */
switch(mode) {
case 0:
case 1:
- ilim = ilim1 = -1;
i = 18;
ndigits = 0;
break;
@@ -318,7 +327,7 @@ dtoa
s = s0 = rv_alloc(i);
#ifdef Honor_FLT_ROUNDS
- if (mode > 1 && rounding != 1)
+ if (mode > 1 && Rounding != 1)
leftright = 0;
#endif
@@ -327,7 +336,7 @@ dtoa
/* Try to get by with floating-point arithmetic. */
i = 0;
- dval(d2) = dval(d);
+ dval(&d2) = dval(&d);
k0 = k;
ilim0 = ilim;
ieps = 2; /* conservative */
@@ -337,7 +346,7 @@ dtoa
if (j & Bletch) {
/* prevent overflows */
j &= Bletch - 1;
- dval(d) /= bigtens[n_bigtens-1];
+ dval(&d) /= bigtens[n_bigtens-1];
ieps++;
}
for(; j; j >>= 1, i++)
@@ -345,32 +354,32 @@ dtoa
ieps++;
ds *= bigtens[i];
}
- dval(d) /= ds;
+ dval(&d) /= ds;
}
else if (( j1 = -k )!=0) {
- dval(d) *= tens[j1 & 0xf];
+ dval(&d) *= tens[j1 & 0xf];
for(j = j1 >> 4; j; j >>= 1, i++)
if (j & 1) {
ieps++;
- dval(d) *= bigtens[i];
+ dval(&d) *= bigtens[i];
}
}
- if (k_check && dval(d) < 1. && ilim > 0) {
+ if (k_check && dval(&d) < 1. && ilim > 0) {
if (ilim1 <= 0)
goto fast_failed;
ilim = ilim1;
k--;
- dval(d) *= 10.;
+ dval(&d) *= 10.;
ieps++;
}
- dval(eps) = ieps*dval(d) + 7.;
- word0(eps) -= (P-1)*Exp_msk1;
+ dval(&eps) = ieps*dval(&d) + 7.;
+ word0(&eps) -= (P-1)*Exp_msk1;
if (ilim == 0) {
S = mhi = 0;
- dval(d) -= 5.;
- if (dval(d) > dval(eps))
+ dval(&d) -= 5.;
+ if (dval(&d) > dval(&eps))
goto one_digit;
- if (dval(d) < -dval(eps))
+ if (dval(&d) < -dval(&eps))
goto no_digits;
goto fast_failed;
}
@@ -379,34 +388,34 @@ dtoa
/* Use Steele & White method of only
* generating digits needed.
*/
- dval(eps) = 0.5/tens[ilim-1] - dval(eps);
+ dval(&eps) = 0.5/tens[ilim-1] - dval(&eps);
for(i = 0;;) {
- L = dval(d);
- dval(d) -= L;
+ L = dval(&d);
+ dval(&d) -= L;
*s++ = '0' + (int)L;
- if (dval(d) < dval(eps))
+ if (dval(&d) < dval(&eps))
goto ret1;
- if (1. - dval(d) < dval(eps))
+ if (1. - dval(&d) < dval(&eps))
goto bump_up;
if (++i >= ilim)
break;
- dval(eps) *= 10.;
- dval(d) *= 10.;
+ dval(&eps) *= 10.;
+ dval(&d) *= 10.;
}
}
else {
#endif
/* Generate ilim digits, then fix them up. */
- dval(eps) *= tens[ilim-1];
- for(i = 1;; i++, dval(d) *= 10.) {
- L = (Long)(dval(d));
- if (!(dval(d) -= L))
+ dval(&eps) *= tens[ilim-1];
+ for(i = 1;; i++, dval(&d) *= 10.) {
+ L = (Long)(dval(&d));
+ if (!(dval(&d) -= L))
ilim = i;
*s++ = '0' + (int)L;
if (i == ilim) {
- if (dval(d) > 0.5 + dval(eps))
+ if (dval(&d) > 0.5 + dval(&eps))
goto bump_up;
- else if (dval(d) < 0.5 - dval(eps)) {
+ else if (dval(&d) < 0.5 - dval(&eps)) {
while(*--s == '0');
s++;
goto ret1;
@@ -419,7 +428,7 @@ dtoa
#endif
fast_failed:
s = s0;
- dval(d) = dval(d2);
+ dval(&d) = dval(&d2);
k = k0;
ilim = ilim0;
}
@@ -431,22 +440,22 @@ dtoa
ds = tens[k];
if (ndigits < 0 && ilim <= 0) {
S = mhi = 0;
- if (ilim < 0 || dval(d) <= 5*ds)
+ if (ilim < 0 || dval(&d) <= 5*ds)
goto no_digits;
goto one_digit;
}
- for(i = 1;; i++, dval(d) *= 10.) {
- L = (Long)(dval(d) / ds);
- dval(d) -= L*ds;
+ for(i = 1;; i++, dval(&d) *= 10.) {
+ L = (Long)(dval(&d) / ds);
+ dval(&d) -= L*ds;
#ifdef Check_FLT_ROUNDS
/* If FLT_ROUNDS == 2, L will usually be high by 1 */
- if (dval(d) < 0) {
+ if (dval(&d) < 0) {
L--;
- dval(d) += ds;
+ dval(&d) += ds;
}
#endif
*s++ = '0' + (int)L;
- if (!dval(d)) {
+ if (!dval(&d)) {
#ifdef SET_INEXACT
inexact = 0;
#endif
@@ -455,13 +464,18 @@ dtoa
if (i == ilim) {
#ifdef Honor_FLT_ROUNDS
if (mode > 1)
- switch(rounding) {
+ switch(Rounding) {
case 0: goto ret1;
case 2: goto bump_up;
}
#endif
- dval(d) += dval(d);
- if (dval(d) > ds || dval(d) == ds && L & 1) {
+ dval(&d) += dval(&d);
+#ifdef ROUND_BIASED
+ if (dval(&d) >= ds)
+#else
+ if (dval(&d) > ds || (dval(&d) == ds && L & 1))
+#endif
+ {
bump_up:
while(*--s == '9')
if (s == s0) {
@@ -523,12 +537,12 @@ dtoa
spec_case = 0;
if ((mode < 2 || leftright)
#ifdef Honor_FLT_ROUNDS
- && rounding == 1
+ && Rounding == 1
#endif
) {
- if (!word1(d) && !(word0(d) & Bndry_mask)
+ if (!word1(&d) && !(word0(&d) & Bndry_mask)
#ifndef Sudden_Underflow
- && word0(d) & (Exp_mask & ~Exp_msk1)
+ && word0(&d) & (Exp_mask & ~Exp_msk1)
#endif
) {
/* The special case */
@@ -614,9 +628,9 @@ dtoa
j1 = delta->sign ? 1 : cmp(b, delta);
Bfree(delta);
#ifndef ROUND_BIASED
- if (j1 == 0 && mode != 1 && !(word1(d) & 1)
+ if (j1 == 0 && mode != 1 && !(word1(&d) & 1)
#ifdef Honor_FLT_ROUNDS
- && rounding >= 1
+ && Rounding >= 1
#endif
) {
if (dig == '9')
@@ -631,11 +645,11 @@ dtoa
goto ret;
}
#endif
- if (j < 0 || j == 0 && mode != 1
+ if (j < 0 || (j == 0 && mode != 1
#ifndef ROUND_BIASED
- && !(word1(d) & 1)
+ && !(word1(&d) & 1)
#endif
- ) {
+ )) {
if (!b->x[0] && b->wds <= 1) {
#ifdef SET_INEXACT
inexact = 0;
@@ -644,7 +658,7 @@ dtoa
}
#ifdef Honor_FLT_ROUNDS
if (mode > 1)
- switch(rounding) {
+ switch(Rounding) {
case 0: goto accept_dig;
case 2: goto keep_dig;
}
@@ -652,7 +666,11 @@ dtoa
if (j1 > 0) {
b = lshift(b, 1);
j1 = cmp(b, S);
- if ((j1 > 0 || j1 == 0 && dig & 1)
+#ifdef ROUND_BIASED
+ if (j1 >= 0 /*)*/
+#else
+ if ((j1 > 0 || (j1 == 0 && dig & 1))
+#endif
&& dig++ == '9')
goto round_9_up;
}
@@ -662,7 +680,7 @@ dtoa
}
if (j1 > 0) {
#ifdef Honor_FLT_ROUNDS
- if (!rounding)
+ if (!Rounding)
goto accept_dig;
#endif
if (dig == '9') { /* possible if i == 1 */
@@ -705,14 +723,19 @@ dtoa
/* Round off last digit */
#ifdef Honor_FLT_ROUNDS
- switch(rounding) {
+ switch(Rounding) {
case 0: goto trimzeros;
case 2: goto roundoff;
}
#endif
b = lshift(b, 1);
j = cmp(b, S);
- if (j > 0 || j == 0 && dig & 1) {
+#ifdef ROUND_BIASED
+ if (j >= 0)
+#else
+ if (j > 0 || (j == 0 && dig & 1))
+#endif
+ {
roundoff:
while(*--s == '9')
if (s == s0) {
@@ -723,7 +746,9 @@ dtoa
++*s++;
}
else {
+#ifdef Honor_FLT_ROUNDS
trimzeros:
+#endif
while(*--s == '0');
s++;
}
@@ -738,9 +763,9 @@ dtoa
#ifdef SET_INEXACT
if (inexact) {
if (!oldinexact) {
- word0(d) = Exp_1 + (70 << Exp_shift);
- word1(d) = 0;
- dval(d) += 1.;
+ word0(&d) = Exp_1 + (70 << Exp_shift);
+ word1(&d) = 0;
+ dval(&d) += 1.;
}
}
else if (!oldinexact)
View
17 lib/libc/gdtoa/g_Qfmt.c
@@ -1,5 +1,3 @@
-/* $NetBSD: g_Qfmt.c,v 1.1.1.1 2006/01/25 15:18:41 kleink Exp $ */
-
/****************************************************************
The author of this software is David M. Gay.
@@ -53,15 +51,20 @@ THIS SOFTWARE.
char*
#ifdef KR_headers
-g_Qfmt(buf, V, ndig, bufsize) char *buf; char *V; int ndig; unsigned bufsize;
+g_Qfmt(buf, V, ndig, bufsize) char *buf; char *V; int ndig; size_t bufsize;
#else
-g_Qfmt(char *buf, void *V, int ndig, unsigned bufsize)
+g_Qfmt(char *buf, void *V, int ndig, size_t bufsize)
#endif
{
- static FPI fpi = { 113, 1-16383-113+1, 32766 - 16383 - 113 + 1, 1, 0 };
+ static FPI fpi0 = { 113, 1-16383-113+1, 32766 - 16383 - 113 + 1, 1, 0 };
char *b, *s, *se;
ULong bits[4], *L, sign;
int decpt, ex, i, mode;
+#ifdef Honor_FLT_ROUNDS
+#include "gdtoa_fltrnds.h"
+#else
+#define fpi &fpi0
+#endif
if (ndig < 0)
ndig = 0;
@@ -111,6 +114,6 @@ g_Qfmt(char *buf, void *V, int ndig, unsigned bufsize)
return 0;
mode = 0;
}
- s = gdtoa(&fpi, ex, bits, &i, mode, ndig, &decpt, &se);
- return g__fmt(buf, s, se, decpt, sign);
+ s = gdtoa(fpi, ex, bits, &i, mode, ndig, &decpt, &se);
+ return g__fmt(buf, s, se, decpt, sign, bufsize);
}
View
72 lib/libc/gdtoa/g__fmt.c
@@ -1,5 +1,3 @@
-/* $NetBSD: g__fmt.c,v 1.1.1.1 2006/01/25 15:18:41 kleink Exp $ */
-
/****************************************************************
The author of this software is David M. Gay.
@@ -39,24 +37,51 @@ THIS SOFTWARE.
char *
#ifdef KR_headers
-g__fmt(b, s, se, decpt, sign) char *b; char *s; char *se; int decpt; ULong sign;
+g__fmt(b, s, se, decpt, sign, blen) char *b; char *s; char *se; int decpt; ULong sign; size_t blen;
#else
-g__fmt(char *b, char *s, char *se, int decpt, ULong sign)
+g__fmt(char *b, char *s, char *se, int decpt, ULong sign, size_t blen)
#endif
{
int i, j, k;
- char *s0 = s;
+ char *be, *s0;
+ size_t len;
#ifdef USE_LOCALE
- char decimalpoint = *localeconv()->decimal_point;
+#ifdef NO_LOCALE_CACHE
+ char *decimalpoint = localeconv()->decimal_point;
+ size_t dlen = strlen(decimalpoint);
+#else
+ char *decimalpoint;
+ static char *decimalpoint_cache;
+ static size_t dlen;
+ if (!(s0 = decimalpoint_cache)) {
+ s0 = localeconv()->decimal_point;
+ dlen = strlen(s0);
+ if ((decimalpoint_cache = (char*)MALLOC(strlen(s0) + 1))) {
+ strcpy(decimalpoint_cache, s0);
+ s0 = decimalpoint_cache;
+ }
+ }
+ decimalpoint = s0;
+#endif
#else
-#define decimalpoint '.'
+#define dlen 0
#endif
+ s0 = s;
+ len = (se-s) + dlen + 6; /* 6 = sign + e+dd + trailing null */
+ if (blen < len)
+ goto ret0;
+ be = b + blen - 1;
if (sign)
*b++ = '-';
if (decpt <= -4 || decpt > se - s + 5) {
*b++ = *s++;
if (*s) {
- *b++ = decimalpoint;
+#ifdef USE_LOCALE
+ while((*b = *decimalpoint++))
+ ++b;
+#else
+ *b++ = '.';
+#endif
while((*b = *s++) !=0)
b++;
}
@@ -71,6 +96,8 @@ g__fmt(char *b, char *s, char *se, int decpt, ULong sign)
for(j = 2, k = 10; 10*k <= decpt; j++, k *= 10){}
for(;;) {
i = decpt / k;
+ if (b >= be)
+ goto ret0;
*b++ = i + '0';
if (--j <= 0)
break;
@@ -80,22 +107,41 @@ g__fmt(char *b, char *s, char *se, int decpt, ULong sign)
*b = 0;
}
else if (decpt <= 0) {
- *b++ = decimalpoint;
+#ifdef USE_LOCALE
+ while((*b = *decimalpoint++))
+ ++b;
+#else
+ *b++ = '.';
+#endif
+ if (be < b - decpt + (se - s))
+ goto ret0;
for(; decpt < 0; decpt++)
*b++ = '0';
- while((*b = *s++) !=0)
+ while((*b = *s++) != 0)
b++;
}
else {
- while((*b = *s++) !=0) {
+ while((*b = *s++) != 0) {
b++;
- if (--decpt == 0 && *s)
- *b++ = decimalpoint;
+ if (--decpt == 0 && *s) {
+#ifdef USE_LOCALE
+ while(*b = *decimalpoint++)
+ ++b;
+#else
+ *b++ = '.';
+#endif
+ }
+ }
+ if (b + decpt > be) {
+ ret0:
+ b = 0;
+ goto ret;
}
for(; decpt > 0; decpt--)
*b++ = '0';
*b = 0;
}
+ ret:
freedtoa(s0);
return b;
}
View
44 lib/libc/gdtoa/g_ddfmt.c
@@ -1,5 +1,3 @@
-/* $NetBSD: g_ddfmt.c,v 1.1.1.1 2006/01/25 15:18:44 kleink Exp $ */
-
/****************************************************************
The author of this software is David M. Gay.
@@ -35,9 +33,9 @@ THIS SOFTWARE.
char *
#ifdef KR_headers
-g_ddfmt(buf, dd, ndig, bufsize) char *buf; double *dd; int ndig; unsigned bufsize;
+g_ddfmt(buf, dd0, ndig, bufsize) char *buf; double *dd0; int ndig; size_t bufsize;
#else
-g_ddfmt(char *buf, double *dd, int ndig, unsigned bufsize)
+g_ddfmt(char *buf, double *dd0, int ndig, size_t bufsize)
#endif
{
FPI fpi;
@@ -45,12 +43,28 @@ g_ddfmt(char *buf, double *dd, int ndig, unsigned bufsize)
ULong *L, bits0[4], *bits, *zx;
int bx, by, decpt, ex, ey, i, j, mode;
Bigint *x, *y, *z;
- double ddx[2];
+ U *dd, ddx[2];
+#ifdef Honor_FLT_ROUNDS /*{{*/
+ int Rounding;
+#ifdef Trust_FLT_ROUNDS /*{{ only define this if FLT_ROUNDS really works! */
+ Rounding = Flt_Rounds;
+#else /*}{*/
+ Rounding = 1;
+ switch(fegetround()) {
+ case FE_TOWARDZERO: Rounding = 0; break;
+ case FE_UPWARD: Rounding = 2; break;
+ case FE_DOWNWARD: Rounding = 3;
+ }
+#endif /*}}*/
+#else /*}{*/
+#define Rounding FPI_Round_near
+#endif /*}}*/
if (bufsize < 10 || bufsize < ndig + 8)
return 0;
- L = (ULong*)dd;
+ dd = (U*)dd0;
+ L = dd->L;
if ((L[_0] & 0x7ff00000L) == 0x7ff00000L) {
/* Infinity or NaN */
if (L[_0] & 0xfffff || L[_1]) {
@@ -75,7 +89,7 @@ g_ddfmt(char *buf, double *dd, int ndig, unsigned bufsize)
goto nanret;
goto infret;
}
- if (dd[0] + dd[1] == 0.) {
+ if (dval(&dd[0]) + dval(&dd[1]) == 0.) {
b = buf;
#ifndef IGNORE_ZERO_SIGN
if (L[_0] & L[2+_0] & 0x80000000L)
@@ -86,16 +100,16 @@ g_ddfmt(char *buf, double *dd, int ndig, unsigned bufsize)
return b;
}
if ((L[_0] & 0x7ff00000L) < (L[2+_0] & 0x7ff00000L)) {
- ddx[1] = dd[0];
- ddx[0] = dd[1];
+ dval(&ddx[1]) = dval(&dd[0]);
+ dval(&ddx[0]) = dval(&dd[1]);
dd = ddx;
- L = (ULong*)dd;
+ L = dd->L;
}
- z = d2b(dd[0], &ex, &bx);
- if (dd[1] == 0.)
+ z = d2b(dval(&dd[0]), &ex, &bx);
+ if (dval(&dd[1]) == 0.)
goto no_y;
x = z;
- y = d2b(dd[1], &ey, &by);
+ y = d2b(dval(&dd[1]), &ey, &by);
if ( (i = ex - ey) !=0) {
if (i > 0) {
x = lshift(x, i);
@@ -146,11 +160,11 @@ g_ddfmt(char *buf, double *dd, int ndig, unsigned bufsize)
}
fpi.emin = 1-1023-53+1;
fpi.emax = 2046-1023-106+1;
- fpi.rounding = FPI_Round_near;
+ fpi.rounding = Rounding;
fpi.sudden_underflow = 0;
i = STRTOG_Normal;
s = gdtoa(&fpi, ex, bits, &i, mode, ndig, &decpt, &se);
- b = g__fmt(buf, s, se, decpt, z->sign);
+ b = g__fmt(buf, s, se, decpt, z->sign, bufsize);
Bfree(z);
return b;
}
View
26 lib/libc/gdtoa/g_dfmt.c
@@ -1,5 +1,3 @@
-/* $NetBSD: g_dfmt.c,v 1.1.1.1 2006/01/25 15:18:44 kleink Exp $ */
-
/****************************************************************
The author of this software is David M. Gay.
@@ -35,15 +33,20 @@ THIS SOFTWARE.
char*
#ifdef KR_headers
-g_dfmt(buf, d, ndig, bufsize) char *buf; double *d; int ndig; unsigned bufsize;
+g_dfmt(buf, d, ndig, bufsize) char *buf; double *d; int ndig; size_t bufsize;
#else
-g_dfmt(char *buf, double *d, int ndig, unsigned bufsize)
+g_dfmt(char *buf, double *d, int ndig, size_t bufsize)
#endif
{
- static FPI fpi = { 53, 1-1023-53+1, 2046-1023-53+1, 1, 0 };
+ static FPI fpi0 = { 53, 1-1023-53+1, 2046-1023-53+1, 1, 0 };
char *b, *s, *se;
ULong bits[2], *L, sign;
int decpt, ex, i, mode;
+#ifdef Honor_FLT_ROUNDS
+#include "gdtoa_fltrnds.h"
+#else
+#define fpi &fpi0
+#endif
if (ndig < 0)
ndig = 0;
@@ -54,6 +57,8 @@ g_dfmt(char *buf, double *d, int ndig, unsigned bufsize)
sign = L[_0] & 0x80000000L;
if ((L[_0] & 0x7ff00000) == 0x7ff00000) {
/* Infinity or NaN */
+ if (bufsize < 10)
+ return 0;
if (L[_0] & 0xfffff || L[_1]) {
return strcp(buf, "NaN");
}
@@ -80,12 +85,11 @@ g_dfmt(char *buf, double *d, int ndig, unsigned bufsize)
ex = 1;
ex -= 0x3ff + 52;
mode = 2;
- if (ndig <= 0) {
- if (bufsize < 25)
- return 0;
+ if (ndig <= 0)
mode = 0;
- }
i = STRTOG_Normal;
- s = gdtoa(&fpi, ex, bits, &i, mode, ndig, &decpt, &se);
- return g__fmt(buf, s, se, decpt, sign);
+ if (sign)
+ i = STRTOG_Normal | STRTOG_Neg;
+ s = gdtoa(fpi, ex, bits, &i, mode, ndig, &decpt, &se);
+ return g__fmt(buf, s, se, decpt, sign, bufsize);
}
View
17 lib/libc/gdtoa/g_ffmt.c
@@ -1,5 +1,3 @@
-/* $NetBSD: g_ffmt.c,v 1.1.1.1 2006/01/25 15:18:44 kleink Exp $ */
-
/****************************************************************
The author of this software is David M. Gay.
@@ -35,15 +33,20 @@ THIS SOFTWARE.
char*
#ifdef KR_headers
-g_ffmt(buf, f, ndig, bufsize) char *buf; float *f; int ndig; unsigned bufsize;
+g_ffmt(buf, f, ndig, bufsize) char *buf; float *f; int ndig; size_t bufsize;
#else
-g_ffmt(char *buf, float *f, int ndig, unsigned bufsize)
+g_ffmt(char *buf, float *f, int ndig, size_t bufsize)
#endif
{
- static FPI fpi = { 24, 1-127-24+1, 254-127-24+1, 1, 0 };
+ static FPI fpi0 = { 24, 1-127-24+1, 254-127-24+1, 1, 0 };
char *b, *s, *se;
ULong bits[1], *L, sign;
int decpt, ex, i, mode;
+#ifdef Honor_FLT_ROUNDS
+#include "gdtoa_fltrnds.h"
+#else
+#define fpi &fpi0
+#endif
if (ndig < 0)
ndig = 0;
@@ -85,6 +88,6 @@ g_ffmt(char *buf, float *f, int ndig, unsigned bufsize)
mode = 0;
}
i = STRTOG_Normal;
- s = gdtoa(&fpi, ex, bits, &i, mode, ndig, &decpt, &se);
- return g__fmt(buf, s, se, decpt, sign);
+ s = gdtoa(fpi, ex, bits, &i, mode, ndig, &decpt, &se);
+ return g__fmt(buf, s, se, decpt, sign, bufsize);
}
View
17 lib/libc/gdtoa/g_xLfmt.c
@@ -1,5 +1,3 @@
-/* $NetBSD: g_xLfmt.c,v 1.1.1.1 2006/01/25 15:18:44 kleink Exp $ */
-
/****************************************************************
The author of this software is David M. Gay.
@@ -51,15 +49,20 @@ THIS SOFTWARE.
char*
#ifdef KR_headers
-g_xLfmt(buf, V, ndig, bufsize) char *buf; char *V; int ndig; unsigned bufsize;
+g_xLfmt(buf, V, ndig, bufsize) char *buf; char *V; int ndig; size_t bufsize;
#else
-g_xLfmt(char *buf, void *V, int ndig, unsigned bufsize)
+g_xLfmt(char *buf, void *V, int ndig, size_t bufsize)
#endif
{
- static FPI fpi = { 64, 1-16383-64+1, 32766 - 16383 - 64 + 1, 1, 0 };
+ static FPI fpi0 = { 64, 1-16383-64+1, 32766 - 16383 - 64 + 1, 1, 0 };
char *b, *s, *se;
ULong bits[2], *L, sign;
int decpt, ex, i, mode;
+#ifdef Honor_FLT_ROUNDS
+#include "gdtoa_fltrnds.h"
+#else
+#define fpi &fpi0
+#endif
if (ndig < 0)
ndig = 0;
@@ -105,6 +108,6 @@ g_xLfmt(char *buf, void *V, int ndig, unsigned bufsize)
return 0;
mode = 0;
}
- s = gdtoa(&fpi, ex, bits, &i, mode, ndig, &decpt, &se);
- return g__fmt(buf, s, se, decpt, sign);
+ s = gdtoa(fpi, ex, bits, &i, mode, ndig, &decpt, &se);
+ return g__fmt(buf, s, se, decpt, sign, bufsize);
}
View
17 lib/libc/gdtoa/g_xfmt.c
@@ -1,5 +1,3 @@
-/* $NetBSD: g_xfmt.c,v 1.1.1.1 2006/01/25 15:18:44 kleink Exp $ */
-
/****************************************************************
The author of this software is David M. Gay.
@@ -55,16 +53,21 @@ THIS SOFTWARE.
char*
#ifdef KR_headers
-g_xfmt(buf, V, ndig, bufsize) char *buf; char *V; int ndig; unsigned bufsize;
+g_xfmt(buf, V, ndig, bufsize) char *buf; char *V; int ndig; size_t bufsize;
#else
-g_xfmt(char *buf, void *V, int ndig, unsigned bufsize)
+g_xfmt(char *buf, void *V, int ndig, size_t bufsize)
#endif
{
- static FPI fpi = { 64, 1-16383-64+1, 32766 - 16383 - 64 + 1, 1, 0 };
+ static FPI fpi0 = { 64, 1-16383-64+1, 32766 - 16383 - 64 + 1, 1, 0 };
char *b, *s, *se;
ULong bits[2], sign;
UShort *L;
int decpt, ex, i, mode;
+#ifdef Honor_FLT_ROUNDS
+#include "gdtoa_fltrnds.h"
+#else
+#define fpi &fpi0
+#endif
if (ndig < 0)
ndig = 0;
@@ -111,6 +114,6 @@ g_xfmt(char *buf, void *V, int ndig, unsigned bufsize)
return 0;
mode = 0;
}
- s = gdtoa(&fpi, ex, bits, &i, mode, ndig, &decpt, &se);
- return g__fmt(buf, s, se, decpt, sign);
+ s = gdtoa(fpi, ex, bits, &i, mode, ndig, &decpt, &se);
+ return g__fmt(buf, s, se, decpt, sign, bufsize);
}
View
186 lib/libc/gdtoa/gdtoa.c
@@ -1,5 +1,3 @@
-/* $NetBSD: gdtoa.c,v 1.1.1.1 2006/01/25 15:18:45 kleink Exp $ */
-
/****************************************************************
The author of this software is David M. Gay.
@@ -125,6 +123,7 @@ gdtoa
the returned string. If not null, *rve is set to point
to the end of the return value. If d is +-Infinity or NaN,
then *decpt is set to 9999.
+ be = exponent: value = (integer represented by bits) * (2 to the power of be).
mode:
0 ==> shortest string that yields d when read in
@@ -159,8 +158,9 @@ gdtoa
int rdir, s2, s5, spec_case, try_quick;
Long L;
Bigint *b, *b1, *delta, *mlo, *mhi, *mhi1, *S;
- double d, d2, ds, eps;
+ double d2, ds;
char *s, *s0;
+ U d, eps;
#ifndef MULTIPLE_THREADS
if (dtoa_result) {
@@ -199,21 +199,21 @@ gdtoa
return nrv_alloc("0", rve, 1);
}
- dval(d) = b2d(b, &i);
+ dval(&d) = b2d(b, &i);
i = be + bbits - 1;
- word0(d) &= Frac_mask1;
- word0(d) |= Exp_11;
+ word0(&d) &= Frac_mask1;
+ word0(&d) |= Exp_11;
#ifdef IBM
- if ( (j = 11 - hi0bits(word0(d) & Frac_mask)) !=0)
- dval(d) /= 1 << j;
+ if ( (j = 11 - hi0bits(word0(&d) & Frac_mask)) !=0)
+ dval(&d) /= 1 << j;
#endif
/* log(x) ~=~ log(1.5) + (x-1.5)/1.5
* log10(x) = log(x) / log(10)
* ~=~ log(1.5)/log(10) + (x-1.5)/(1.5*log(10))
- * log10(d) = (i-Bias)*log(2)/log(10) + log10(d2)
+ * log10(&d) = (i-Bias)*log(2)/log(10) + log10(d2)
*
- * This suggests computing an approximation k to log10(d) by
+ * This suggests computing an approximation k to log10(&d) by
*
* k = (i - Bias)*0.301029995663981
* + ( (d2-1.5)*0.289529654602168 + 0.176091259055681 );
@@ -233,7 +233,7 @@ gdtoa
i <<= 2;
i += j;
#endif
- ds = (dval(d)-1.5)*0.289529654602168 + 0.1760912590558 + i*0.301029995663981;
+ ds = (dval(&d)-1.5)*0.289529654602168 + 0.1760912590558 + i*0.301029995663981;
/* correct assumption about exponent range */
if ((j = i) < 0)
@@ -248,13 +248,13 @@ gdtoa
#ifdef IBM
j = be + bbits - 1;
if ( (j1 = j & 3) !=0)
- dval(d) *= 1 << j1;
- word0(d) += j << Exp_shift - 2 & Exp_mask;
+ dval(&d) *= 1 << j1;
+ word0(&d) += j << Exp_shift - 2 & Exp_mask;
#else
- word0(d) += (be + bbits - 1) << Exp_shift;
+ word0(&d) += (be + bbits - 1) << Exp_shift;
#endif
if (k >= 0 && k <= Ten_pmax) {
- if (dval(d) < tens[k])
+ if (dval(&d) < tens[k])
k--;
k_check = 0;
}
@@ -284,11 +284,14 @@ gdtoa
mode -= 4;
try_quick = 0;
}
+ else if (i >= -4 - Emin || i < Emin)
+ try_quick = 0;
leftright = 1;
+ ilim = ilim1 = -1; /* Values for cases 0 and 1; done here to */
+ /* silence erroneous "gcc -Wall" warning. */
switch(mode) {
case 0:
case 1:
- ilim = ilim1 = -1;
i = (int)(nbits * .30103) + 3;
ndigits = 0;
break;
@@ -330,10 +333,10 @@ gdtoa
/* Try to get by with floating-point arithmetic. */
i = 0;
- d2 = dval(d);
+ d2 = dval(&d);
#ifdef IBM
- if ( (j = 11 - hi0bits(word0(d) & Frac_mask)) !=0)
- dval(d) /= 1 << j;
+ if ( (j = 11 - hi0bits(word0(&d) & Frac_mask)) !=0)
+ dval(&d) /= 1 << j;
#endif
k0 = k;
ilim0 = ilim;
@@ -344,7 +347,7 @@ gdtoa
if (j & Bletch) {
/* prevent overflows */
j &= Bletch - 1;
- dval(d) /= bigtens[n_bigtens-1];
+ dval(&d) /= bigtens[n_bigtens-1];
ieps++;
}
for(; j; j >>= 1, i++)
@@ -356,30 +359,30 @@ gdtoa
else {
ds = 1.;
if ( (j1 = -k) !=0) {
- dval(d) *= tens[j1 & 0xf];
+ dval(&d) *= tens[j1 & 0xf];
for(j = j1 >> 4; j; j >>= 1, i++)
if (j & 1) {
ieps++;
- dval(d) *= bigtens[i];
+ dval(&d) *= bigtens[i];
}
}
}
- if (k_check && dval(d) < 1. && ilim > 0) {
+ if (k_check && dval(&d) < 1. && ilim > 0) {
if (ilim1 <= 0)
goto fast_failed;
ilim = ilim1;
k--;
- dval(d) *= 10.;
+ dval(&d) *= 10.;
ieps++;
}
- dval(eps) = ieps*dval(d) + 7.;
- word0(eps) -= (P-1)*Exp_msk1;
+ dval(&eps) = ieps*dval(&d) + 7.;
+ word0(&eps) -= (P-1)*Exp_msk1;
if (ilim == 0) {
S = mhi = 0;
- dval(d) -= 5.;
- if (dval(d) > dval(eps))
+ dval(&d) -= 5.;
+ if (dval(&d) > dval(&eps))
goto one_digit;
- if (dval(d) < -dval(eps))
+ if (dval(&d) < -dval(&eps))
goto no_digits;
goto fast_failed;
}
@@ -388,42 +391,40 @@ gdtoa
/* Use Steele & White method of only
* generating digits needed.
*/
- dval(eps) = ds*0.5/tens[ilim-1] - dval(eps);
+ dval(&eps) = ds*0.5/tens[ilim-1] - dval(&eps);
for(i = 0;;) {
- L = (Long)(dval(d)/ds);
- dval(d) -= L*ds;
+ L = (Long)(dval(&d)/ds);
+ dval(&d) -= L*ds;
*s++ = '0' + (int)L;
- if (dval(d) < dval(eps)) {
- if (dval(d))
+ if (dval(&d) < dval(&eps)) {
+ if (dval(&d))
inex = STRTOG_Inexlo;
goto ret1;
}
- if (ds - dval(d) < dval(eps))
+ if (ds - dval(&d) < dval(&eps))
goto bump_up;
if (++i >= ilim)
break;
- dval(eps) *= 10.;
- dval(d) *= 10.;
+ dval(&eps) *= 10.;
+ dval(&d) *= 10.;
}
}
else {
#endif
/* Generate ilim digits, then fix them up. */
- dval(eps) *= tens[ilim-1];
- for(i = 1;; i++, dval(d) *= 10.) {
- if ( (L = (Long)(dval(d)/ds)) !=0)
- dval(d) -= L*ds;
+ dval(&eps) *= tens[ilim-1];
+ for(i = 1;; i++, dval(&d) *= 10.) {
+ if ( (L = (Long)(dval(&d)/ds)) !=0)
+ dval(&d) -= L*ds;
*s++ = '0' + (int)L;
if (i == ilim) {
ds *= 0.5;
- if (dval(d) > ds + dval(eps))
+ if (dval(&d) > ds + dval(&eps))
goto bump_up;
- else if (dval(d) < ds - dval(eps)) {
- while(*--s == '0'){}
- s++;
- if (dval(d))
+ else if (dval(&d) < ds - dval(&eps)) {
+ if (dval(&d))
inex = STRTOG_Inexlo;
- goto ret1;
+ goto clear_trailing0;
}
break;
}
@@ -433,7 +434,7 @@ gdtoa
#endif
fast_failed:
s = s0;
- dval(d) = d2;
+ dval(&d) = d2;
k = k0;
ilim = ilim0;
}
@@ -445,22 +446,22 @@ gdtoa
ds = tens[k];
if (ndigits < 0 && ilim <= 0) {
S = mhi = 0;
- if (ilim < 0 || dval(d) <= 5*ds)
+ if (ilim < 0 || dval(&d) <= 5*ds)
goto no_digits;
goto one_digit;
}
- for(i = 1;; i++, dval(d) *= 10.) {
- L = dval(d) / ds;
- dval(d) -= L*ds;
+ for(i = 1;; i++, dval(&d) *= 10.) {
+ L = dval(&d) / ds;
+ dval(&d) -= L*ds;
#ifdef Check_FLT_ROUNDS
/* If FLT_ROUNDS == 2, L will usually be high by 1 */
- if (dval(d) < 0) {
+ if (dval(&d) < 0) {
L--;
- dval(d) += ds;
+ dval(&d) += ds;
}
#endif
*s++ = '0' + (int)L;
- if (dval(d) == 0.)
+ if (dval(&d) == 0.)
break;
if (i == ilim) {
if (rdir) {
@@ -469,8 +470,13 @@ gdtoa
inex = STRTOG_Inexlo;
goto ret1;
}
- dval(d) += dval(d);
- if (dval(d) > ds || dval(d) == ds && L & 1) {
+ dval(&d) += dval(&d);
+#ifdef ROUND_BIASED
+ if (dval(&d) >= ds)
+#else
+ if (dval(&d) > ds || (dval(&d) == ds && L & 1))
+#endif
+ {
bump_up:
inex = STRTOG_Inexhi;
while(*--s == '9')
@@ -481,8 +487,12 @@ gdtoa
}
++*s++;
}
- else
+ else {
inex = STRTOG_Inexlo;
+ clear_trailing0:
+ while(*--s == '0'){}
+ ++s;
+ }
break;
}
}
@@ -493,13 +503,15 @@ gdtoa
m5 = b5;
mhi = mlo = 0;
if (leftright) {
- if (mode < 2) {
- i = nbits - bbits;
- if (be - i++ < fpi->emin)
- /* denormal */
- i = be - fpi->emin + 1;
+ i = nbits - bbits;
+ if (be - i++ < fpi->emin && mode != 3 && mode != 5) {
+ /* denormal */
+ i = be - fpi->emin + 1;
+ if (mode >= 2 && ilim > 0 && ilim < i)
+ goto small_ilim;
}
- else {
+ else if (mode >= 2) {
+ small_ilim:
j = ilim - 1;
if (m5 >= j)
m5 -= j;
@@ -560,28 +572,11 @@ gdtoa
* and for all and pass them and a shift to quorem, so it
* can do shifts and ors to compute the numerator for q.
*/
-#ifdef Pack_32
- if ( (i = ((s5 ? 32 - hi0bits(S->x[S->wds-1]) : 1) + s2) & 0x1f) !=0)
- i = 32 - i;
-#else
- if ( (i = ((s5 ? 32 - hi0bits(S->x[S->wds-1]) : 1) + s2) & 0xf) !=0)
- i = 16 - i;
-#endif
- if (i > 4) {
- i -= 4;
- b2 += i;
- m2 += i;
- s2 += i;
- }
- else if (i < 4) {
- i += 28;
- b2 += i;
- m2 += i;
- s2 += i;
- }
- if (b2 > 0)
+ i = ((s5 ? hi0bits(S->x[S->wds-1]) : ULbits - 1) - s2 - 4) & kmask;
+ m2 += i;
+ if ((b2 += i) > 0)
b = lshift(b, b2);
- if (s2 > 0)
+ if ((s2 += i) > 0)
S = lshift(S, s2);
if (k_check) {
if (cmp(b,S) < 0) {
@@ -646,11 +641,11 @@ gdtoa
goto ret;
}
#endif
- if (j < 0 || j == 0 && !mode
+ if (j < 0 || (j == 0 && !mode
#ifndef ROUND_BIASED
&& !(bits[0] & 1)
#endif
- ) {
+ )) {
if (rdir && (b->wds > 1 || b->x[0])) {
if (rdir == 2) {
inex = STRTOG_Inexlo;
@@ -673,7 +668,11 @@ gdtoa
if (j1 > 0) {
b = lshift(b, 1);
j1 = cmp(b, S);
- if ((j1 > 0 || j1 == 0 && dig & 1)
+#ifdef ROUND_BIASED
+ if (j1 >= 0 /*)*/
+#else
+ if ((j1 > 0 || (j1 == 0 && dig & 1))
+#endif
&& dig++ == '9')
goto round_9_up;
inex = STRTOG_Inexhi;
@@ -718,13 +717,18 @@ gdtoa
/* Round off last digit */
if (rdir) {
- if (rdir == 2 || b->wds <= 1 && !b->x[0])
+ if (rdir == 2 || (b->wds <= 1 && !b->x[0]))
goto chopzeros;
goto roundoff;
}
b = lshift(b, 1);
j = cmp(b, S);
- if (j > 0 || j == 0 && dig & 1) {
+#ifdef ROUND_BIASED
+ if (j >= 0)
+#else
+ if (j > 0 || (j == 0 && dig & 1))
+#endif
+ {
roundoff:
inex = STRTOG_Inexhi;
while(*--s == '9')
@@ -740,7 +744,7 @@ gdtoa
if (b->wds > 1 || b->x[0])
inex = STRTOG_Inexlo;
while(*--s == '0'){}
- s++;
+ ++s;
}
ret:
Bfree(S);
View
23 lib/libc/gdtoa/gdtoa.h
@@ -1,5 +1,3 @@
-/* $NetBSD: gdtoa.h,v 1.1.1.1 2006/01/25 15:18:46 kleink Exp $ */
-
/****************************************************************
The author of this software is David M. Gay.
@@ -35,9 +33,10 @@ THIS SOFTWARE.
#define GDTOA_H_INCLUDED
#include "arith.h"
+#include <stddef.h> /* for size_t */
#ifndef Long
-#define Long long
+#define Long int
#endif
#ifndef ULong
typedef unsigned Long ULong;
@@ -76,9 +75,9 @@ typedef unsigned short UShort;
/* The following may be or-ed into one of the above values. */
- STRTOG_Neg = 0x08,
- STRTOG_Inexlo = 0x10,
- STRTOG_Inexhi = 0x20,
+ STRTOG_Neg = 0x08, /* does not affect STRTOG_Inexlo or STRTOG_Inexhi */
+ STRTOG_Inexlo = 0x10, /* returned result rounded toward zero */
+ STRTOG_Inexhi = 0x20, /* returned result rounded away from zero */
STRTOG_Inexact = 0x30,
STRTOG_Underflow= 0x40,
STRTOG_Overflow = 0x80
@@ -113,12 +112,12 @@ extern float strtof ANSI((CONST char *, char **));
extern double strtod ANSI((CONST char *, char **));
extern int strtodg ANSI((CONST char*, char**, FPI*, Long*, ULong*));
-extern char* g_ddfmt ANSI((char*, double*, int, unsigned));
-extern char* g_dfmt ANSI((char*, double*, int, unsigned));
-extern char* g_ffmt ANSI((char*, float*, int, unsigned));
-extern char* g_Qfmt ANSI((char*, void*, int, unsigned));
-extern char* g_xfmt ANSI((char*, void*, int, unsigned));
-extern char* g_xLfmt ANSI((char*, void*, int, unsigned));
+extern char* g_ddfmt ANSI((char*, double*, int, size_t));
+extern char* g_dfmt ANSI((char*, double*, int, size_t));
+extern char* g_ffmt ANSI((char*, float*, int, size_t));
+extern char* g_Qfmt ANSI((char*, void*, int, size_t));
+extern char* g_xfmt ANSI((char*, void*, int, size_t));
+extern char* g_xLfmt ANSI((char*, void*, int, size_t));
extern int strtoId ANSI((CONST char*, char**, double*, double*));
extern int strtoIdd ANSI((CONST char*, char**, double*, double*));
View
18 lib/libc/gdtoa/gdtoa_fltrnds.h
@@ -0,0 +1,18 @@
+ FPI *fpi, fpi1;
+ int Rounding;
+#ifdef Trust_FLT_ROUNDS /*{{ only define this if FLT_ROUNDS really works! */
+ Rounding = Flt_Rounds;
+#else /*}{*/
+ Rounding = 1;
+ switch(fegetround()) {
+ case FE_TOWARDZERO: Rounding = 0; break;
+ case FE_UPWARD: Rounding = 2; break;
+ case FE_DOWNWARD: Rounding = 3;
+ }
+#endif /*}}*/
+ fpi = &fpi0;
+ if (Rounding != 1) {
+ fpi1 = fpi0;
+ fpi = &fpi1;
+ fpi1.rounding = Rounding;
+ }
View
89 lib/libc/gdtoa/gdtoaimp.h
@@ -1,5 +1,3 @@
-/* $NetBSD: gdtoaimp.h,v 1.1.1.1 2006/01/25 15:18:48 kleink Exp $ */
-
/****************************************************************
The author of this software is David M. Gay.
@@ -96,7 +94,12 @@ THIS SOFTWARE.
* #define RND_PRODQUOT to use rnd_prod and rnd_quot (assembly routines
* that use extended-precision instructions to compute rounded
* products and quotients) with IBM.
- * #define ROUND_BIASED for IEEE-format with biased rounding.
+ * #define ROUND_BIASED for IEEE-format with biased rounding and arithmetic
+ * that rounds toward +Infinity.
+ * #define ROUND_BIASED_without_Round_Up for IEEE-format with biased
+ * rounding when the underlying floating-point arithmetic uses
+ * unbiased rounding. This prevent using ordinary floating-point
+ * arithmetic when the result could be computed with one rounding error.
* #define Inaccurate_Divide for IEEE-format with correctly rounded
* products but inaccurate quotients, e.g., for Intel i860.
* #define NO_LONG_LONG on machines that do not have a "long long"
@@ -115,7 +118,12 @@ THIS SOFTWARE.
* #define MALLOC your_malloc, where your_malloc(n) acts like malloc(n)
* if memory is available and otherwise does something you deem
* appropriate. If MALLOC is undefined, malloc will be invoked
- * directly -- and assumed always to succeed.
+ * directly -- and assumed always to succeed. Similarly, if you
+ * want something other than the system's free() to be called to
+ * recycle memory acquired from MALLOC, #define FREE to be the
+ * name of the alternate routine. (FREE or free is only called in
+ * pathological cases, e.g., in a gdtoa call after a gdtoa return in
+ * mode 3 with thousands of digits requested.)
* #define Omit_Private_Memory to omit logic (added Jan. 1998) for making
* memory allocations from a private pool of memory when possible.
* When used, the private pool is PRIVATE_MEM bytes long: 2304 bytes,
@@ -128,17 +136,22 @@ THIS SOFTWARE.
* conversions of IEEE doubles in single-threaded executions with
* 8-byte pointers, PRIVATE_MEM >= 7400 appears to suffice; with
* 4-byte pointers, PRIVATE_MEM >= 7112 appears adequate.
- * #define INFNAN_CHECK on IEEE systems to cause strtod to check for
- * Infinity and NaN (case insensitively).
+ * #define NO_INFNAN_CHECK if you do not wish to have INFNAN_CHECK
+ * #defined automatically on IEEE systems. On such systems,
+ * when INFNAN_CHECK is #defined, strtod checks
+ * for Infinity and NaN (case insensitively).
* When INFNAN_CHECK is #defined and No_Hex_NaN is not #defined,
* strtodg also accepts (case insensitively) strings of the form
- * NaN(x), where x is a string of hexadecimal digits and spaces;
- * if there is only one string of hexadecimal digits, it is taken
- * for the fraction bits of the resulting NaN; if there are two or
- * more strings of hexadecimal digits, each string is assigned
- * to the next available sequence of 32-bit words of fractions
- * bits (starting with the most significant), right-aligned in
- * each sequence.
+ * NaN(x), where x is a string of hexadecimal digits (optionally
+ * preceded by 0x or 0X) and spaces; if there is only one string
+ * of hexadecimal digits, it is taken for the fraction bits of the
+ * resulting NaN; if there are two or more strings of hexadecimal
+ * digits, each string is assigned to the next available sequence
+ * of 32-bit words of fractions bits (starting with the most
+ * significant), right-aligned in each sequence.
+ * Unless GDTOA_NON_PEDANTIC_NANCHECK is #defined, input "NaN(...)"
+ * is consumed even when ... has the wrong form (in which case the
+ * "(...)" is consumed but ignored).
* #define MULTIPLE_THREADS if the system offers preemptively scheduled
* multiple threads. In this case, you must provide (or suitably
* #define) two locks, acquired by ACQUIRE_DTOA_LOCK(n) and freed
@@ -150,7 +163,7 @@ THIS SOFTWARE.
* dtoa. You may do so whether or not MULTIPLE_THREADS is #defined.
* #define IMPRECISE_INEXACT if you do not care about the setting of
* the STRTOG_Inexact bits in the special case of doing IEEE double
- * precision conversions (which could also be done by the strtog in
+ * precision conversions (which could also be done by the strtod in
* dtoa.c).
* #define NO_HEX_FP to disable recognition of C9x's hexadecimal
* floating-point constants.
@@ -159,11 +172,6 @@ THIS SOFTWARE.
* #define NO_STRING_H to use private versions of memcpy.
* On some K&R systems, it may also be necessary to
* #define DECLARE_SIZE_T in this case.
- * #define YES_ALIAS to permit aliasing certain double values with
- * arrays of ULongs. This leads to slightly better code with
- * some compilers and was always used prior to 19990916, but it
- * is not strictly legal and can cause trouble with aggressively
- * optimizing compilers (e.g., gcc 2.95.1 under -O2).
* #define USE_LOCALE to use the current locale's decimal_point value.
*/
@@ -171,6 +179,9 @@ THIS SOFTWARE.
#define GDTOAIMP_H_INCLUDED
#include "gdtoa.h"
#include "gd_qnan.h"
+#ifdef Honor_FLT_ROUNDS
+#include <fenv.h>
+#endif
#ifdef DEBUG
#include "stdio.h"
@@ -264,25 +275,14 @@ Exactly one of IEEE_8087, IEEE_MC68k, VAX, or IBM should be defined.
typedef union { double d; ULong L[2]; } U;
-#ifdef YES_ALIAS
-#define dval(x) x
-#ifdef IEEE_8087
-#define word0(x) ((ULong *)&x)[1]
-#define word1(x) ((ULong *)&x)[0]
-#else
-#define word0(x) ((ULong *)&x)[0]
-#define word1(x) ((ULong *)&x)[1]
-#endif
-#else /* !YES_ALIAS */
#ifdef IEEE_8087
-#define word0(x) ((U*)&x)->L[1]
-#define word1(x) ((U*)&x)->L[0]
+#define word0(x) (x)->L[1]
+#define word1(x) (x)->L[0]
#else
-#define word0(x) ((U*)&x)->L[0]
-#define word1(x) ((U*)&x)->L[1]
+#define word0(x) (x)->L[0]
+#define word1(x) (x)->L[1]
#endif
-#define dval(x) ((U*)&x)->d
-#endif /* YES_ALIAS */
+#define dval(x) (x)->d
/* The following definition of Storeinc is appropriate for MIPS processors.
* An alternative that might be better on some machines is
@@ -396,6 +396,11 @@ typedef union { double d; ULong L[2]; } U;
#ifndef IEEE_Arith
#define ROUND_BIASED
+#else
+#ifdef ROUND_BIASED_without_Round_Up
+#undef ROUND_BIASED
+#define ROUND_BIASED
+#endif
#endif
#ifdef RND_PRODQUOT
@@ -456,7 +461,7 @@ extern double rnd_prod(double, double), rnd_quot(double, double);
#define FREE_DTOA_LOCK(n) /*nothing*/
#endif
-#define Kmax 15
+#define Kmax 9
struct
Bigint {
@@ -540,11 +545,11 @@ extern void memcpy_D2A ANSI((void*, const void*, size_t));
extern int cmp ANSI((Bigint*, Bigint*));
extern void copybits ANSI((ULong*, int, Bigint*));
extern Bigint *d2b ANSI((double, int*, int*));
- extern int decrement ANSI((Bigint*));
+ extern void decrement ANSI((Bigint*));
extern Bigint *diff ANSI((Bigint*, Bigint*));
extern char *dtoa ANSI((double d, int mode, int ndigits,
int *decpt, int *sign, char **rve));
- extern char *g__fmt ANSI((char*, char*, char*, int, ULong));
+ extern char *g__fmt ANSI((char*, char*, char*, int, ULong, size_t));
extern int gethex ANSI((CONST char**, FPI*, Long*, Bigint**, int));
extern void hexdig_init_D2A(Void);
extern int hexnan ANSI((CONST char**, FPI*, ULong*));
@@ -562,14 +567,14 @@ extern void memcpy_D2A ANSI((void*, const void*, size_t));
extern double ratio ANSI((Bigint*, Bigint*));
extern void rshift ANSI((Bigint*, int));
extern char *rv_alloc ANSI((int));
- extern Bigint *s2b ANSI((CONST char*, int, int, ULong));
+ extern Bigint *s2b ANSI((CONST char*, int, int, ULong, int));
extern Bigint *set_ones ANSI((Bigint*, int));
extern char *strcp ANSI((char*, const char*));
extern int strtoIg ANSI((CONST char*, char**, FPI*, Long*, Bigint**, int*));
extern double strtod ANSI((const char *s00, char **se));
extern Bigint *sum ANSI((Bigint*, Bigint*));
extern int trailz ANSI((Bigint*));
- extern double ulp ANSI((double));
+ extern double ulp ANSI((U*));
#ifdef __cplusplus
}
@@ -584,6 +589,10 @@ extern void memcpy_D2A ANSI((void*, const void*, size_t));
* (On HP Series 700/800 machines, -DNAN_WORD0=0x7ff40000 works.)
*/
#ifdef IEEE_Arith
+#ifndef NO_INFNAN_CHECK
+#undef INFNAN_CHECK
+#define INFNAN_CHECK
+#endif
#ifdef IEEE_MC68k
#define _0 0
#define _1 1
View
147 lib/libc/gdtoa/gethex.c
@@ -1,5 +1,3 @@
-/* $NetBSD: gethex.c,v 1.1.1.1 2006/01/25 15:18:48 kleink Exp $ */
-
/****************************************************************
The author of this software is David M. Gay.
@@ -47,17 +45,30 @@ gethex( CONST char **sp, FPI *fpi, Long *exp, Bigint **bp, int sign)
{
Bigint *b;
CONST unsigned char *decpt, *s0, *s, *s1;
- int esign, havedig, irv, k, n, nbits, up, zret;
+ int big, esign, havedig, irv, j, k, n, n0, nbits, up, zret;
ULong L, lostbits, *x;
Long e, e1;
#ifdef USE_LOCALE
- unsigned char decimalpoint = *localeconv()->decimal_point;
+ int i;
+#ifdef NO_LOCALE_CACHE
+ const unsigned char *decimalpoint = (unsigned char*)localeconv()->decimal_point;
#else
-#define decimalpoint '.'
+ const unsigned char *decimalpoint;
+ static unsigned char *decimalpoint_cache;
+ if (!(s0 = decimalpoint_cache)) {
+ s0 = (unsigned char*)localeconv()->decimal_point;
+ if ((decimalpoint_cache = (char*)MALLOC(strlen(s0) + 1))) {
+ strcpy(decimalpoint_cache, s0);
+ s0 = decimalpoint_cache;
+ }
+ }
+ decimalpoint = s0;
+#endif
#endif
if (!hexdig['0'])
hexdig_init_D2A();
+ *bp = 0;
havedig = 0;
s0 = *(CONST unsigned char **)sp + 2;
while(s0[havedig] == '0')
@@ -67,11 +78,21 @@ gethex( CONST char **sp, FPI *fpi, Long *exp, Bigint **bp, int sign)
decpt = 0;
zret = 0;
e = 0;
- if (!hexdig[*s]) {
+ if (hexdig[*s])
+ havedig++;
+ else {
zret = 1;
- if (*s != decimalpoint)
+#ifdef USE_LOCALE
+ for(i = 0; decimalpoint[i]; ++i) {
+ if (s[i] != decimalpoint[i])
+ goto pcheck;
+ }
+ decpt = s += i;
+#else
+ if (*s != '.')
goto pcheck;
decpt = ++s;
+#endif
if (!hexdig[*s])
goto pcheck;
while(*s == '0')
@@ -83,19 +104,28 @@ gethex( CONST char **sp, FPI *fpi, Long *exp, Bigint **bp, int sign)
}
while(hexdig[*s])
s++;
- if (*s == decimalpoint && !decpt) {
+#ifdef USE_LOCALE
+ if (*s == *decimalpoint && !decpt) {
+ for(i = 1; decimalpoint[i]; ++i) {
+ if (s[i] != decimalpoint[i])
+ goto pcheck;
+ }
+ decpt = s += i;
+#else
+ if (*s == '.' && !decpt) {
decpt = ++s;
+#endif
while(hexdig[*s])
s++;
- }
+ }/*}*/
if (decpt)
e = -(((Long)(s-decpt)) << 2);
pcheck:
s1 = s;
+ big = esign = 0;
switch(*s) {
case 'p':
case 'P':
- esign = 0;
switch(*++s) {
case '-':
esign = 1;
@@ -108,26 +138,87 @@ gethex( CONST char **sp, FPI *fpi, Long *exp, Bigint **bp, int sign)
break;
}
e1 = n - 0x10;
- while((n = hexdig[*++s]) !=0 && n <= 0x19)
+ while((n = hexdig[*++s]) !=0 && n <= 0x19) {
+ if (e1 & 0xf8000000)
+ big = 1;
e1 = 10*e1 + n - 0x10;
+ }
if (esign)
e1 = -e1;
e += e1;
}
*sp = (char*)s;
+ if (!havedig)
+ *sp = (char*)s0 - 1;
if (zret)
- return havedig ? STRTOG_Zero : STRTOG_NoNumber;
+ return STRTOG_Zero;
+ if (big) {
+ if (esign) {
+ switch(fpi->rounding) {
+ case FPI_Round_up:
+ if (sign)
+ break;
+ goto ret_tiny;
+ case FPI_Round_down:
+ if (!sign)
+ break;
+ goto ret_tiny;
+ }
+ goto retz;
+ ret_tiny:
+ b = Balloc(0);
+ b->wds = 1;
+ b->x[0] = 1;
+ goto dret;
+ }
+ switch(fpi->rounding) {
+ case FPI_Round_near:
+ goto ovfl1;
+ case FPI_Round_up:
+ if (!sign)
+ goto ovfl1;
+ goto ret_big;
+ case FPI_Round_down:
+ if (sign)
+ goto ovfl1;
+ goto ret_big;
+ }
+ ret_big:
+ nbits = fpi->nbits;
+ n0 = n = nbits >> kshift;
+ if (nbits & kmask)
+ ++n;
+ for(j = n, k = 0; j >>= 1; ++k);
+ *bp = b = Balloc(k);
+ b->wds = n;
+ for(j = 0; j < n0; ++j)
+ b->x[j] = ALL_ON;
+ if (n > n0)
+ b->x[j] = ULbits >> (ULbits - (nbits & kmask));
+ *exp = fpi->emin;
+ return STRTOG_Normal | STRTOG_Inexlo;
+ }
n = s1 - s0 - 1;
- for(k = 0; n > 7; n >>= 1)
+ for(k = 0; n > (1 << (kshift-2)) - 1; n >>= 1)
k++;
b = Balloc(k);
x = b->x;
n = 0;
L = 0;
+#ifdef USE_LOCALE
+ for(i = 0; decimalpoint[i+1]; ++i);
+#endif
while(s1 > s0) {
- if (*--s1 == decimalpoint)
+#ifdef USE_LOCALE
+ if (*--s1 == decimalpoint[i]) {
+ s1 -= i;
continue;
- if (n == 32) {
+ }
+#else
+ if (*--s1 == '.')
+ continue;
+#endif
+ if (n == ULbits) {
*x++ = L;
L = 0;
n = 0;
@@ -137,7 +228,7 @@ gethex( CONST char **sp, FPI *fpi, Long *exp, Bigint **bp, int sign)
}
*x++ = L;
b->wds = n = x - b->x;
- n = 32*n - hi0bits(L);
+ n = ULbits*n - hi0bits(L);
nbits = fpi->nbits;
lostbits = 0;
x = b->x;
@@ -148,7 +239,7 @@ gethex( CONST char **sp, FPI *fpi, Long *exp, Bigint **bp, int sign)
k = n - 1;
if (x[k>>kshift] & 1 << (k & kmask)) {
lostbits = 2;
- if (k > 1 && any_on(b,k-1))
+ if (k > 0 && any_on(b,k))
lostbits = 3;
}
}
@@ -164,7 +255,10 @@ gethex( CONST char **sp, FPI *fpi, Long *exp, Bigint **bp, int sign)
if (e > fpi->emax) {
ovfl:
Bfree(b);
- *bp = 0;
+ ovfl1:
+#ifndef NO_ERRNO
+ errno = ERANGE;
+#endif
return STRTOG_Infinite | STRTOG_Overflow | STRTOG_Inexhi;
}
irv = STRTOG_Normal;
@@ -184,15 +278,22 @@ gethex( CONST char **sp, FPI *fpi, Long *exp, Bigint **bp, int sign)
case FPI_Round_down:
if (sign) {
one_bit:
- *exp = fpi->emin;
x[0] = b->wds = 1;
+ dret:
*bp = b;
+ *exp = fpi->emin;
+#ifndef NO_ERRNO
+ errno = ERANGE;
+#endif
return STRTOG_Denormal | STRTOG_Inexhi
| STRTOG_Underflow;
}
}
Bfree(b);
- *bp = 0;
+ retz:
+#ifndef NO_ERRNO
+ errno = ERANGE;
+#endif
return STRTOG_Zero | STRTOG_Inexlo | STRTOG_Underflow;
}
k = n - 1;
@@ -213,7 +314,7 @@ gethex( CONST char **sp, FPI *fpi, Long *exp, Bigint **bp, int sign)
break;
case FPI_Round_near:
if (lostbits & 2
- && (lostbits & 1) | x[0] & 1)
+ && (lostbits | x[0]) & 1)
up = 1;
break;
case FPI_Round_up:
@@ -232,8 +333,8 @@ gethex( CONST char **sp, FPI *fpi, Long *exp, Bigint **bp, int sign)
irv = STRTOG_Normal;
}
else if (b->wds > k
- || (n = nbits & kmask) !=0
- && hi0bits(x[k-1]) < 32-n) {
+ || ((n = nbits & kmask) !=0
+ && hi0bits(x[k-1]) < 32-n)) {
rshift(b,1);
if (++e > fpi->emax)
goto ovfl;
View
2  lib/libc/gdtoa/gmisc.c
@@ -1,5 +1,3 @@
-/* $NetBSD: gmisc.c,v 1.1.1.1 2006/01/25 15:18:48 kleink Exp $ */
-
/****************************************************************
The author of this software is David M. Gay.
View
2  lib/libc/gdtoa/hd_init.c
@@ -1,5 +1,3 @@
-/* $NetBSD: hd_init.c,v 1.1.1.1 2006/01/25 15:18:48 kleink Exp $ */
-
/****************************************************************
The author of this software is David M. Gay.
View
25 lib/libc/gdtoa/hexnan.c
@@ -1,5 +1,3 @@
-/* $NetBSD: hexnan.c,v 1.1.1.1 2006/01/25 15:18:48 kleink Exp $ */
-
/****************************************************************
The author of this software is David M. Gay.
@@ -73,7 +71,13 @@ hexnan( CONST char **sp, FPI *fpi, ULong *x0)
x1 = xe = x;
havedig = hd0 = i = 0;
s = *sp;
- while(c = *(CONST unsigned char*)++s) {
+ /* allow optional initial 0x or 0X */
+ while((c = *(CONST unsigned char*)(s+1)) && c <= ' ')
+ ++s;
+ if (s[1] == '0' && (s[2] == 'x' || s[2] == 'X')
+ && *(CONST unsigned char*)(s+3) > ' ')
+ s += 2;
+ while((c = *(CONST unsigned char*)++s)) {
if (!(h = hexdig[c])) {
if (c <= ' ') {
if (hd0 < havedig) {
@@ -88,12 +92,25 @@ hexnan( CONST char **sp, FPI *fpi, ULong *x0)
x1 = x;
i = 0;
}
+ while(*(CONST unsigned char*)(s+1) <= ' ')
+ ++s;
+ if (s[1] == '0' && (s[2] == 'x' || s[2] == 'X')
+ && *(CONST unsigned char*)(s+3) > ' ')
+ s += 2;
continue;
}
if (/*(*/ c == ')' && havedig) {
*sp = s + 1;
break;
}
+#ifndef GDTOA_NON_PEDANTIC_NANCHECK
+ do {
+ if (/*(*/ c == ')') {
+ *sp = s + 1;
+ break;
+ }
+ } while((c = *++s));
+#endif
return STRTOG_NaN;
}
havedig++;
@@ -103,7 +120,7 @@ hexnan( CONST char **sp, FPI *fpi, ULong *x0)
i = 1;
*--x = 0;
}
- *x = (*x << 4) | h & 0xf;
+ *x = (*x << 4) | (h & 0xf);
}
if (!havedig)
return STRTOG_NaN;
View
45 lib/libc/gdtoa/makefile
@@ -1,5 +1,3 @@
-/* $NetBSD: makefile,v 1.1.1.1 2006/01/25 15:18:48 kleink Exp $ */
-
# /****************************************************************