1 /**************************************************************** 2 3 The author of this software (_dtoa, strtod) is David M. Gay. 4 Please send bug reports to 5 David M. Gay 6 Bell Laboratories, Room 2C-463 7 600 Mountain Avenue 8 Murray Hill, NJ 07974-2070 9 U.S.A. 10 dmg@research.bell-labs.com 11 */ 12 #include <stdlib.h> 13 #include <string.h> 14 #define _RESEARCH_SOURCE 15 #include <errno.h> 16 #include <float.h> 17 #include <math.h> 18 #define CONST const 19 20 /* 21 * #define IEEE_8087 for IEEE-arithmetic machines where the least 22 * significant byte has the lowest address. 23 * #define IEEE_MC68k for IEEE-arithmetic machines where the most 24 * significant byte has the lowest address. 25 * #define Sudden_Underflow for IEEE-format machines without gradual 26 * underflow (i.e., that flush to zero on underflow). 27 * #define IBM for IBM mainframe-style floating-point arithmetic. 28 * #define VAX for VAX-style floating-point arithmetic. 29 * #define Unsigned_Shifts if >> does treats its left operand as unsigned. 30 * #define No_leftright to omit left-right logic in fast floating-point 31 * computation of dtoa. 32 * #define Check_FLT_ROUNDS if FLT_ROUNDS can assume the values 2 or 3. 33 * #define RND_PRODQUOT to use rnd_prod and rnd_quot (assembly routines 34 * that use extended-precision instructions to compute rounded 35 * products and quotients) with IBM. 36 * #define ROUND_BIASED for IEEE-format with biased rounding. 37 * #define Inaccurate_Divide for IEEE-format with correctly rounded 38 * products but inaccurate quotients, e.g., for Intel i860. 39 * #define Just_16 to store 16 bits per 32-bit long when doing high-precision 40 * integer arithmetic. Whether this speeds things up or slows things 41 * down depends on the machine and the number being converted. 42 */ 43 44 #if defined(IEEE_8087) + defined(IEEE_MC68k) + defined(VAX) + defined(IBM) != 1 45 Exactly one of IEEE_8087, IEEE_MC68k, VAX, or IBM should be defined. 46 #endif 47 48 typedef union { 49 double d; 50 unsigned long ul[2]; 51 } Dul; 52 53 #ifdef IEEE_8087 54 #define word0(x) ((x).ul[1]) 55 #define word1(x) ((x).ul[0]) 56 #else 57 #define word0(x) ((x).ul[0]) 58 #define word1(x) ((x).ul[1]) 59 #endif 60 61 #ifdef Unsigned_Shifts 62 #define Sign_Extend(a,b) if (b < 0) a |= 0xffff0000; 63 #else 64 #define Sign_Extend(a,b) /*no-op*/ 65 #endif 66 67 /* The following definition of Storeinc is appropriate for MIPS processors. 68 * An alternative that might be better on some machines is 69 * #define Storeinc(a,b,c) (*a++ = b << 16 | c & 0xffff) 70 */ 71 #if defined(IEEE_8087) + defined(VAX) 72 #define Storeinc(a,b,c) (((unsigned short *)a)[1] = (unsigned short)b, \ 73 ((unsigned short *)a)[0] = (unsigned short)c, a++) 74 #else 75 #define Storeinc(a,b,c) (((unsigned short *)a)[0] = (unsigned short)b, \ 76 ((unsigned short *)a)[1] = (unsigned short)c, a++) 77 #endif 78 79 /* #define P DBL_MANT_DIG */ 80 /* Ten_pmax = floor(P*log(2)/log(5)) */ 81 /* Bletch = (highest power of 2 < DBL_MAX_10_EXP) / 16 */ 82 /* Quick_max = floor((P-1)*log(FLT_RADIX)/log(10) - 1) */ 83 /* Int_max = floor(P*log(FLT_RADIX)/log(10) - 1) */ 84 85 #if defined(IEEE_8087) + defined(IEEE_MC68k) 86 #define Exp_shift 20 87 #define Exp_shift1 20 88 #define Exp_msk1 0x100000 89 #define Exp_msk11 0x100000 90 #define Exp_mask 0x7ff00000 91 #define P 53 92 #define Bias 1023 93 #define IEEE_Arith 94 #define Emin (-1022) 95 #define Exp_1 0x3ff00000 96 #define Exp_11 0x3ff00000 97 #define Ebits 11 98 #define Frac_mask 0xfffff 99 #define Frac_mask1 0xfffff 100 #define Ten_pmax 22 101 #define Bletch 0x10 102 #define Bndry_mask 0xfffff 103 #define Bndry_mask1 0xfffff 104 #define LSB 1 105 #define Sign_bit 0x80000000 106 #define Log2P 1 107 #define Tiny0 0 108 #define Tiny1 1 109 #define Quick_max 14 110 #define Int_max 14 111 #define Infinite(x) (word0(x) == 0x7ff00000) /* sufficient test for here */ 112 #else 113 #undef Sudden_Underflow 114 #define Sudden_Underflow 115 #ifdef IBM 116 #define Exp_shift 24 117 #define Exp_shift1 24 118 #define Exp_msk1 0x1000000 119 #define Exp_msk11 0x1000000 120 #define Exp_mask 0x7f000000 121 #define P 14 122 #define Bias 65 123 #define Exp_1 0x41000000 124 #define Exp_11 0x41000000 125 #define Ebits 8 /* exponent has 7 bits, but 8 is the right value in b2d */ 126 #define Frac_mask 0xffffff 127 #define Frac_mask1 0xffffff 128 #define Bletch 4 129 #define Ten_pmax 22 130 #define Bndry_mask 0xefffff 131 #define Bndry_mask1 0xffffff 132 #define LSB 1 133 #define Sign_bit 0x80000000 134 #define Log2P 4 135 #define Tiny0 0x100000 136 #define Tiny1 0 137 #define Quick_max 14 138 #define Int_max 15 139 #else /* VAX */ 140 #define Exp_shift 23 141 #define Exp_shift1 7 142 #define Exp_msk1 0x80 143 #define Exp_msk11 0x800000 144 #define Exp_mask 0x7f80 145 #define P 56 146 #define Bias 129 147 #define Exp_1 0x40800000 148 #define Exp_11 0x4080 149 #define Ebits 8 150 #define Frac_mask 0x7fffff 151 #define Frac_mask1 0xffff007f 152 #define Ten_pmax 24 153 #define Bletch 2 154 #define Bndry_mask 0xffff007f 155 #define Bndry_mask1 0xffff007f 156 #define LSB 0x10000 157 #define Sign_bit 0x8000 158 #define Log2P 1 159 #define Tiny0 0x80 160 #define Tiny1 0 161 #define Quick_max 15 162 #define Int_max 15 163 #endif 164 #endif 165 166 #ifndef IEEE_Arith 167 #define ROUND_BIASED 168 #endif 169 170 #define Big0 (Frac_mask1 | Exp_msk1*(DBL_MAX_EXP+Bias-1)) 171 #define Big1 0xffffffff 172 173 #ifndef Just_16 174 /* When Pack_32 is not defined, we store 16 bits per 32-bit long. 175 * This makes some inner loops simpler and sometimes saves work 176 * during multiplications, but it often seems to make things slightly 177 * slower. Hence the default is now to store 32 bits per long. 178 */ 179 #ifndef Pack_32 180 #define Pack_32 181 #endif 182 #endif 183 184 #define Kmax 15 185 186 struct 187 Bigint { 188 struct Bigint *next; 189 int k, maxwds, sign, wds; 190 unsigned long x[1]; 191 }; 192 193 typedef struct Bigint Bigint; 194 195 /* This routines shouldn't be visible externally */ 196 extern Bigint *_Balloc(int); 197 extern void _Bfree(Bigint *); 198 extern Bigint *_multadd(Bigint *, int, int); 199 extern int _hi0bits(unsigned long); 200 extern Bigint *_mult(Bigint *, Bigint *); 201 extern Bigint *_pow5mult(Bigint *, int); 202 extern Bigint *_lshift(Bigint *, int); 203 extern int _cmp(Bigint *, Bigint *); 204 extern Bigint *_diff(Bigint *, Bigint *); 205 extern Bigint *_d2b(double, int *, int *); 206 extern Bigint *_i2b(int); 207 208 extern double _tens[], _bigtens[], _tinytens[]; 209 210 #ifdef IEEE_Arith 211 #define n_bigtens 5 212 #else 213 #ifdef IBM 214 #define n_bigtens 3 215 #else 216 #define n_bigtens 2 217 #endif 218 #endif 219 220 #define Balloc(x) _Balloc(x) 221 #define Bfree(x) _Bfree(x) 222 #define Bcopy(x,y) memcpy((char *)&x->sign, (char *)&y->sign, \ 223 y->wds*sizeof(long) + 2*sizeof(int)) 224 #define multadd(x,y,z) _multadd(x,y,z) 225 #define hi0bits(x) _hi0bits(x) 226 #define i2b(x) _i2b(x) 227 #define mult(x,y) _mult(x,y) 228 #define pow5mult(x,y) _pow5mult(x,y) 229 #define lshift(x,y) _lshift(x,y) 230 #define cmp(x,y) _cmp(x,y) 231 #define diff(x,y) _diff(x,y) 232 #define d2b(x,y,z) _d2b(x,y,z) 233 234 #define tens _tens 235 #define bigtens _bigtens 236 #define tinytens _tinytens 237