1 |
|
|
/**************************************************************** |
2 |
|
|
|
3 |
|
|
The author of this software is David M. Gay. |
4 |
|
|
|
5 |
|
|
Copyright (C) 1998-2001 by Lucent Technologies |
6 |
|
|
All Rights Reserved |
7 |
|
|
|
8 |
|
|
Permission to use, copy, modify, and distribute this software and |
9 |
|
|
its documentation for any purpose and without fee is hereby |
10 |
|
|
granted, provided that the above copyright notice appear in all |
11 |
|
|
copies and that both that the copyright notice and this |
12 |
|
|
permission notice and warranty disclaimer appear in supporting |
13 |
|
|
documentation, and that the name of Lucent or any of its entities |
14 |
|
|
not be used in advertising or publicity pertaining to |
15 |
|
|
distribution of the software without specific, written prior |
16 |
|
|
permission. |
17 |
|
|
|
18 |
|
|
LUCENT DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE, |
19 |
|
|
INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS. |
20 |
|
|
IN NO EVENT SHALL LUCENT OR ANY OF ITS ENTITIES BE LIABLE FOR ANY |
21 |
|
|
SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES |
22 |
|
|
WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER |
23 |
|
|
IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, |
24 |
|
|
ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF |
25 |
|
|
THIS SOFTWARE. |
26 |
|
|
|
27 |
|
|
****************************************************************/ |
28 |
|
|
|
29 |
|
|
/* Please send bug reports to David M. Gay (dmg at acm dot org, |
30 |
|
|
* with " at " changed at "@" and " dot " changed to "."). */ |
31 |
|
|
|
32 |
|
|
#include "gdtoaimp.h" |
33 |
|
|
#ifndef NO_FENV_H |
34 |
|
|
#include <fenv.h> |
35 |
|
|
#endif |
36 |
|
|
|
37 |
|
|
#ifdef USE_LOCALE |
38 |
|
|
#include "locale.h" |
39 |
|
|
#endif |
40 |
|
|
|
41 |
|
|
#ifdef IEEE_Arith |
42 |
|
|
#ifndef NO_IEEE_Scale |
43 |
|
|
#define Avoid_Underflow |
44 |
|
|
#undef tinytens |
45 |
|
|
/* The factor of 2^106 in tinytens[4] helps us avoid setting the underflow */ |
46 |
|
|
/* flag unnecessarily. It leads to a song and dance at the end of strtod. */ |
47 |
|
|
static CONST double tinytens[] = { 1e-16, 1e-32, 1e-64, 1e-128, |
48 |
|
|
9007199254740992.*9007199254740992.e-256 |
49 |
|
|
}; |
50 |
|
|
#endif |
51 |
|
|
#endif |
52 |
|
|
|
53 |
|
|
#ifdef Honor_FLT_ROUNDS |
54 |
|
|
#undef Check_FLT_ROUNDS |
55 |
|
|
#define Check_FLT_ROUNDS |
56 |
|
|
#else |
57 |
|
|
#define Rounding Flt_Rounds |
58 |
|
|
#endif |
59 |
|
|
|
60 |
|
|
#ifdef Avoid_Underflow /*{*/ |
61 |
|
|
static double |
62 |
|
|
sulp |
63 |
|
|
#ifdef KR_headers |
64 |
|
|
(x, scale) U *x; int scale; |
65 |
|
|
#else |
66 |
|
|
(U *x, int scale) |
67 |
|
|
#endif |
68 |
|
|
{ |
69 |
|
|
U u; |
70 |
|
|
double rv; |
71 |
|
|
int i; |
72 |
|
|
|
73 |
|
|
rv = ulp(x); |
74 |
|
|
if (!scale || (i = 2*P + 1 - ((word0(x) & Exp_mask) >> Exp_shift)) <= 0) |
75 |
|
|
return rv; /* Is there an example where i <= 0 ? */ |
76 |
|
|
word0(&u) = Exp_1 + (i << Exp_shift); |
77 |
|
|
word1(&u) = 0; |
78 |
|
|
return rv * u.d; |
79 |
|
|
} |
80 |
|
|
#endif /*}*/ |
81 |
|
|
|
82 |
|
|
double |
83 |
|
|
strtod |
84 |
|
|
#ifdef KR_headers |
85 |
|
|
(s00, se) CONST char *s00; char **se; |
86 |
|
|
#else |
87 |
|
|
(CONST char *s00, char **se) |
88 |
|
|
#endif |
89 |
|
|
{ |
90 |
|
|
#ifdef Avoid_Underflow |
91 |
|
|
int scale; |
92 |
|
|
#endif |
93 |
|
|
int bb2, bb5, bbe, bd2, bd5, bbbits, bs2, c, decpt, dsign, |
94 |
|
|
e, e1, esign, i, j, k, nd, nd0, nf, nz, nz0, sign; |
95 |
|
|
CONST char *s, *s0, *s1; |
96 |
|
|
double aadj; |
97 |
|
|
Long L; |
98 |
|
|
U adj, aadj1, rv, rv0; |
99 |
|
|
ULong y, z; |
100 |
|
|
Bigint *bb = NULL, *bb1, *bd = NULL, *bd0 = NULL, *bs = NULL, *delta = NULL; |
101 |
|
|
#ifdef Avoid_Underflow |
102 |
|
|
ULong Lsb, Lsb1; |
103 |
|
|
#endif |
104 |
|
|
#ifdef SET_INEXACT |
105 |
|
|
int inexact, oldinexact; |
106 |
|
|
#endif |
107 |
|
|
#ifdef USE_LOCALE /*{{*/ |
108 |
|
|
#ifdef NO_LOCALE_CACHE |
109 |
|
|
char *decimalpoint = localeconv()->decimal_point; |
110 |
|
|
int dplen = strlen(decimalpoint); |
111 |
|
|
#else |
112 |
|
|
char *decimalpoint; |
113 |
|
|
static char *decimalpoint_cache; |
114 |
|
|
static int dplen; |
115 |
|
|
if (!(s0 = decimalpoint_cache)) { |
116 |
|
|
s0 = localeconv()->decimal_point; |
117 |
|
|
decimalpoint_cache = strdup(s0); |
118 |
|
|
dplen = strlen(s0); |
119 |
|
|
} |
120 |
|
|
decimalpoint = (char*)s0; |
121 |
|
|
#endif /*NO_LOCALE_CACHE*/ |
122 |
|
|
#else /*USE_LOCALE}{*/ |
123 |
|
|
#define dplen 1 |
124 |
|
|
#endif /*USE_LOCALE}}*/ |
125 |
|
|
|
126 |
|
|
#ifdef Honor_FLT_ROUNDS /*{*/ |
127 |
|
|
int Rounding; |
128 |
|
|
#ifdef Trust_FLT_ROUNDS /*{{ only define this if FLT_ROUNDS really works! */ |
129 |
|
|
Rounding = Flt_Rounds; |
130 |
|
|
#else /*}{*/ |
131 |
|
|
Rounding = 1; |
132 |
|
|
switch(fegetround()) { |
133 |
|
|
case FE_TOWARDZERO: Rounding = 0; break; |
134 |
|
|
case FE_UPWARD: Rounding = 2; break; |
135 |
|
|
case FE_DOWNWARD: Rounding = 3; |
136 |
|
|
} |
137 |
|
|
#endif /*}}*/ |
138 |
|
|
#endif /*}*/ |
139 |
|
|
|
140 |
|
|
sign = nz0 = nz = decpt = 0; |
141 |
|
|
dval(&rv) = 0.; |
142 |
|
|
for(s = s00;;s++) switch(*s) { |
143 |
|
|
case '-': |
144 |
|
|
sign = 1; |
145 |
|
|
/* no break */ |
146 |
|
|
case '+': |
147 |
|
|
if (*++s) |
148 |
|
|
goto break2; |
149 |
|
|
/* no break */ |
150 |
|
|
case 0: |
151 |
|
|
goto ret0; |
152 |
|
|
case '\t': |
153 |
|
|
case '\n': |
154 |
|
|
case '\v': |
155 |
|
|
case '\f': |
156 |
|
|
case '\r': |
157 |
|
|
case ' ': |
158 |
|
|
continue; |
159 |
|
|
default: |
160 |
|
|
goto break2; |
161 |
|
|
} |
162 |
|
|
break2: |
163 |
|
|
if (*s == '0') { |
164 |
|
|
#ifndef NO_HEX_FP /*{*/ |
165 |
|
|
{ |
166 |
|
|
static FPI fpi = { 53, 1-1023-53+1, 2046-1023-53+1, 1, SI }; |
167 |
|
|
Long exp; |
168 |
|
|
ULong bits[2]; |
169 |
|
|
switch(s[1]) { |
170 |
|
|
case 'x': |
171 |
|
|
case 'X': |
172 |
|
|
{ |
173 |
|
|
#ifdef Honor_FLT_ROUNDS |
174 |
|
|
FPI fpi1 = fpi; |
175 |
|
|
fpi1.rounding = Rounding; |
176 |
|
|
#else |
177 |
|
|
#define fpi1 fpi |
178 |
|
|
#endif |
179 |
|
|
switch((i = gethex(&s, &fpi1, &exp, &bb, sign)) & STRTOG_Retmask) { |
180 |
|
|
case STRTOG_NoMemory: |
181 |
|
|
goto ovfl; |
182 |
|
|
case STRTOG_NoNumber: |
183 |
|
|
s = s00; |
184 |
|
|
sign = 0; |
185 |
|
|
case STRTOG_Zero: |
186 |
|
|
break; |
187 |
|
|
default: |
188 |
|
|
if (bb) { |
189 |
|
|
copybits(bits, fpi.nbits, bb); |
190 |
|
|
Bfree(bb); |
191 |
|
|
} |
192 |
|
|
ULtod(((U*)&rv)->L, bits, exp, i); |
193 |
|
|
}} |
194 |
|
|
goto ret; |
195 |
|
|
} |
196 |
|
|
} |
197 |
|
|
#endif /*}*/ |
198 |
|
|
nz0 = 1; |
199 |
|
|
while(*++s == '0') ; |
200 |
|
|
if (!*s) |
201 |
|
|
goto ret; |
202 |
|
|
} |
203 |
|
|
s0 = s; |
204 |
|
|
y = z = 0; |
205 |
|
|
for(nd = nf = 0; (c = *s) >= '0' && c <= '9'; nd++, s++) |
206 |
|
|
if (nd < 9) |
207 |
|
|
y = 10*y + c - '0'; |
208 |
|
|
else if (nd < 16) |
209 |
|
|
z = 10*z + c - '0'; |
210 |
|
|
nd0 = nd; |
211 |
|
|
#ifdef USE_LOCALE |
212 |
|
|
if (c == *decimalpoint) { |
213 |
|
|
for(i = 1; decimalpoint[i]; ++i) |
214 |
|
|
if (s[i] != decimalpoint[i]) |
215 |
|
|
goto dig_done; |
216 |
|
|
s += i; |
217 |
|
|
c = *s; |
218 |
|
|
#else |
219 |
|
|
if (c == '.') { |
220 |
|
|
c = *++s; |
221 |
|
|
#endif |
222 |
|
|
decpt = 1; |
223 |
|
|
if (!nd) { |
224 |
|
|
for(; c == '0'; c = *++s) |
225 |
|
|
nz++; |
226 |
|
|
if (c > '0' && c <= '9') { |
227 |
|
|
s0 = s; |
228 |
|
|
nf += nz; |
229 |
|
|
nz = 0; |
230 |
|
|
goto have_dig; |
231 |
|
|
} |
232 |
|
|
goto dig_done; |
233 |
|
|
} |
234 |
|
|
for(; c >= '0' && c <= '9'; c = *++s) { |
235 |
|
|
have_dig: |
236 |
|
|
nz++; |
237 |
|
|
if (c -= '0') { |
238 |
|
|
nf += nz; |
239 |
|
|
for(i = 1; i < nz; i++) |
240 |
|
|
if (nd++ < 9) |
241 |
|
|
y *= 10; |
242 |
|
|
else if (nd <= DBL_DIG + 1) |
243 |
|
|
z *= 10; |
244 |
|
|
if (nd++ < 9) |
245 |
|
|
y = 10*y + c; |
246 |
|
|
else if (nd <= DBL_DIG + 1) |
247 |
|
|
z = 10*z + c; |
248 |
|
|
nz = 0; |
249 |
|
|
} |
250 |
|
|
} |
251 |
|
|
}/*}*/ |
252 |
|
|
dig_done: |
253 |
|
|
e = 0; |
254 |
|
|
if (c == 'e' || c == 'E') { |
255 |
|
|
if (!nd && !nz && !nz0) { |
256 |
|
|
goto ret0; |
257 |
|
|
} |
258 |
|
|
s00 = s; |
259 |
|
|
esign = 0; |
260 |
|
|
switch(c = *++s) { |
261 |
|
|
case '-': |
262 |
|
|
esign = 1; |
263 |
|
|
case '+': |
264 |
|
|
c = *++s; |
265 |
|
|
} |
266 |
|
|
if (c >= '0' && c <= '9') { |
267 |
|
|
while(c == '0') |
268 |
|
|
c = *++s; |
269 |
|
|
if (c > '0' && c <= '9') { |
270 |
|
|
L = c - '0'; |
271 |
|
|
s1 = s; |
272 |
|
|
while((c = *++s) >= '0' && c <= '9') |
273 |
|
|
L = 10*L + c - '0'; |
274 |
|
|
if (s - s1 > 8 || L > 19999) |
275 |
|
|
/* Avoid confusion from exponents |
276 |
|
|
* so large that e might overflow. |
277 |
|
|
*/ |
278 |
|
|
e = 19999; /* safe for 16 bit ints */ |
279 |
|
|
else |
280 |
|
|
e = (int)L; |
281 |
|
|
if (esign) |
282 |
|
|
e = -e; |
283 |
|
|
} |
284 |
|
|
else |
285 |
|
|
e = 0; |
286 |
|
|
} |
287 |
|
|
else |
288 |
|
|
s = s00; |
289 |
|
|
} |
290 |
|
|
if (!nd) { |
291 |
|
|
if (!nz && !nz0) { |
292 |
|
|
#ifdef INFNAN_CHECK |
293 |
|
|
/* Check for Nan and Infinity */ |
294 |
|
|
ULong bits[2]; |
295 |
|
|
static FPI fpinan = /* only 52 explicit bits */ |
296 |
|
|
{ 52, 1-1023-53+1, 2046-1023-53+1, 1, SI }; |
297 |
|
|
if (!decpt) |
298 |
|
|
switch(c) { |
299 |
|
|
case 'i': |
300 |
|
|
case 'I': |
301 |
|
|
if (match(&s,"nf")) { |
302 |
|
|
--s; |
303 |
|
|
if (!match(&s,"inity")) |
304 |
|
|
++s; |
305 |
|
|
word0(&rv) = 0x7ff00000; |
306 |
|
|
word1(&rv) = 0; |
307 |
|
|
goto ret; |
308 |
|
|
} |
309 |
|
|
break; |
310 |
|
|
case 'n': |
311 |
|
|
case 'N': |
312 |
|
|
if (match(&s, "an")) { |
313 |
|
|
#ifndef No_Hex_NaN |
314 |
|
|
if (*s == '(' /*)*/ |
315 |
|
|
&& hexnan(&s, &fpinan, bits) |
316 |
|
|
== STRTOG_NaNbits) { |
317 |
|
|
word0(&rv) = 0x7ff00000 | bits[1]; |
318 |
|
|
word1(&rv) = bits[0]; |
319 |
|
|
} |
320 |
|
|
else { |
321 |
|
|
#endif |
322 |
|
|
word0(&rv) = NAN_WORD0; |
323 |
|
|
word1(&rv) = NAN_WORD1; |
324 |
|
|
#ifndef No_Hex_NaN |
325 |
|
|
} |
326 |
|
|
#endif |
327 |
|
|
goto ret; |
328 |
|
|
} |
329 |
|
|
} |
330 |
|
|
#endif /* INFNAN_CHECK */ |
331 |
|
|
ret0: |
332 |
|
|
s = s00; |
333 |
|
|
sign = 0; |
334 |
|
|
} |
335 |
|
|
goto ret; |
336 |
|
|
} |
337 |
|
|
e1 = e -= nf; |
338 |
|
|
|
339 |
|
|
/* Now we have nd0 digits, starting at s0, followed by a |
340 |
|
|
* decimal point, followed by nd-nd0 digits. The number we're |
341 |
|
|
* after is the integer represented by those digits times |
342 |
|
|
* 10**e */ |
343 |
|
|
|
344 |
|
|
if (!nd0) |
345 |
|
|
nd0 = nd; |
346 |
|
|
k = nd < DBL_DIG + 1 ? nd : DBL_DIG + 1; |
347 |
|
|
dval(&rv) = y; |
348 |
|
|
if (k > 9) { |
349 |
|
|
#ifdef SET_INEXACT |
350 |
|
|
if (k > DBL_DIG) |
351 |
|
|
oldinexact = get_inexact(); |
352 |
|
|
#endif |
353 |
|
|
dval(&rv) = tens[k - 9] * dval(&rv) + z; |
354 |
|
|
} |
355 |
|
|
if (nd <= DBL_DIG |
356 |
|
|
#ifndef RND_PRODQUOT |
357 |
|
|
#ifndef Honor_FLT_ROUNDS |
358 |
|
|
&& Flt_Rounds == 1 |
359 |
|
|
#endif |
360 |
|
|
#endif |
361 |
|
|
) { |
362 |
|
|
if (!e) |
363 |
|
|
goto ret; |
364 |
|
|
#ifndef ROUND_BIASED_without_Round_Up |
365 |
|
|
if (e > 0) { |
366 |
|
|
if (e <= Ten_pmax) { |
367 |
|
|
#ifdef VAX |
368 |
|
|
goto vax_ovfl_check; |
369 |
|
|
#else |
370 |
|
|
#ifdef Honor_FLT_ROUNDS |
371 |
|
|
/* round correctly FLT_ROUNDS = 2 or 3 */ |
372 |
|
|
if (sign) { |
373 |
|
|
rv.d = -rv.d; |
374 |
|
|
sign = 0; |
375 |
|
|
} |
376 |
|
|
#endif |
377 |
|
|
/* rv = */ rounded_product(dval(&rv), tens[e]); |
378 |
|
|
goto ret; |
379 |
|
|
#endif |
380 |
|
|
} |
381 |
|
|
i = DBL_DIG - nd; |
382 |
|
|
if (e <= Ten_pmax + i) { |
383 |
|
|
/* A fancier test would sometimes let us do |
384 |
|
|
* this for larger i values. |
385 |
|
|
*/ |
386 |
|
|
#ifdef Honor_FLT_ROUNDS |
387 |
|
|
/* round correctly FLT_ROUNDS = 2 or 3 */ |
388 |
|
|
if (sign) { |
389 |
|
|
rv.d = -rv.d; |
390 |
|
|
sign = 0; |
391 |
|
|
} |
392 |
|
|
#endif |
393 |
|
|
e -= i; |
394 |
|
|
dval(&rv) *= tens[i]; |
395 |
|
|
#ifdef VAX |
396 |
|
|
/* VAX exponent range is so narrow we must |
397 |
|
|
* worry about overflow here... |
398 |
|
|
*/ |
399 |
|
|
vax_ovfl_check: |
400 |
|
|
word0(&rv) -= P*Exp_msk1; |
401 |
|
|
/* rv = */ rounded_product(dval(&rv), tens[e]); |
402 |
|
|
if ((word0(&rv) & Exp_mask) |
403 |
|
|
> Exp_msk1*(DBL_MAX_EXP+Bias-1-P)) |
404 |
|
|
goto ovfl; |
405 |
|
|
word0(&rv) += P*Exp_msk1; |
406 |
|
|
#else |
407 |
|
|
/* rv = */ rounded_product(dval(&rv), tens[e]); |
408 |
|
|
#endif |
409 |
|
|
goto ret; |
410 |
|
|
} |
411 |
|
|
} |
412 |
|
|
#ifndef Inaccurate_Divide |
413 |
|
|
else if (e >= -Ten_pmax) { |
414 |
|
|
#ifdef Honor_FLT_ROUNDS |
415 |
|
|
/* round correctly FLT_ROUNDS = 2 or 3 */ |
416 |
|
|
if (sign) { |
417 |
|
|
rv.d = -rv.d; |
418 |
|
|
sign = 0; |
419 |
|
|
} |
420 |
|
|
#endif |
421 |
|
|
/* rv = */ rounded_quotient(dval(&rv), tens[-e]); |
422 |
|
|
goto ret; |
423 |
|
|
} |
424 |
|
|
#endif |
425 |
|
|
#endif /* ROUND_BIASED_without_Round_Up */ |
426 |
|
|
} |
427 |
|
|
e1 += nd - k; |
428 |
|
|
|
429 |
|
|
#ifdef IEEE_Arith |
430 |
|
|
#ifdef SET_INEXACT |
431 |
|
|
inexact = 1; |
432 |
|
|
if (k <= DBL_DIG) |
433 |
|
|
oldinexact = get_inexact(); |
434 |
|
|
#endif |
435 |
|
|
#ifdef Avoid_Underflow |
436 |
|
|
scale = 0; |
437 |
|
|
#endif |
438 |
|
|
#ifdef Honor_FLT_ROUNDS |
439 |
|
|
if (Rounding >= 2) { |
440 |
|
|
if (sign) |
441 |
|
|
Rounding = Rounding == 2 ? 0 : 2; |
442 |
|
|
else |
443 |
|
|
if (Rounding != 2) |
444 |
|
|
Rounding = 0; |
445 |
|
|
} |
446 |
|
|
#endif |
447 |
|
|
#endif /*IEEE_Arith*/ |
448 |
|
|
|
449 |
|
|
/* Get starting approximation = rv * 10**e1 */ |
450 |
|
|
|
451 |
|
|
if (e1 > 0) { |
452 |
|
|
if ( (i = e1 & 15) !=0) |
453 |
|
|
dval(&rv) *= tens[i]; |
454 |
|
|
if (e1 &= ~15) { |
455 |
|
|
if (e1 > DBL_MAX_10_EXP) { |
456 |
|
|
ovfl: |
457 |
|
|
/* Can't trust HUGE_VAL */ |
458 |
|
|
#ifdef IEEE_Arith |
459 |
|
|
#ifdef Honor_FLT_ROUNDS |
460 |
|
|
switch(Rounding) { |
461 |
|
|
case 0: /* toward 0 */ |
462 |
|
|
case 3: /* toward -infinity */ |
463 |
|
|
word0(&rv) = Big0; |
464 |
|
|
word1(&rv) = Big1; |
465 |
|
|
break; |
466 |
|
|
default: |
467 |
|
|
word0(&rv) = Exp_mask; |
468 |
|
|
word1(&rv) = 0; |
469 |
|
|
} |
470 |
|
|
#else /*Honor_FLT_ROUNDS*/ |
471 |
|
|
word0(&rv) = Exp_mask; |
472 |
|
|
word1(&rv) = 0; |
473 |
|
|
#endif /*Honor_FLT_ROUNDS*/ |
474 |
|
|
#ifdef SET_INEXACT |
475 |
|
|
/* set overflow bit */ |
476 |
|
|
dval(&rv0) = 1e300; |
477 |
|
|
dval(&rv0) *= dval(&rv0); |
478 |
|
|
#endif |
479 |
|
|
#else /*IEEE_Arith*/ |
480 |
|
|
word0(&rv) = Big0; |
481 |
|
|
word1(&rv) = Big1; |
482 |
|
|
#endif /*IEEE_Arith*/ |
483 |
|
|
range_err: |
484 |
|
|
if (bd0) { |
485 |
|
|
Bfree(bb); |
486 |
|
|
Bfree(bd); |
487 |
|
|
Bfree(bs); |
488 |
|
|
Bfree(bd0); |
489 |
|
|
Bfree(delta); |
490 |
|
|
} |
491 |
|
|
#ifndef NO_ERRNO |
492 |
|
|
errno = ERANGE; |
493 |
|
|
#endif |
494 |
|
|
goto ret; |
495 |
|
|
} |
496 |
|
|
e1 >>= 4; |
497 |
|
|
for(j = 0; e1 > 1; j++, e1 >>= 1) |
498 |
|
|
if (e1 & 1) |
499 |
|
|
dval(&rv) *= bigtens[j]; |
500 |
|
|
/* The last multiplication could overflow. */ |
501 |
|
|
word0(&rv) -= P*Exp_msk1; |
502 |
|
|
dval(&rv) *= bigtens[j]; |
503 |
|
|
if ((z = word0(&rv) & Exp_mask) |
504 |
|
|
> Exp_msk1*(DBL_MAX_EXP+Bias-P)) |
505 |
|
|
goto ovfl; |
506 |
|
|
if (z > Exp_msk1*(DBL_MAX_EXP+Bias-1-P)) { |
507 |
|
|
/* set to largest number */ |
508 |
|
|
/* (Can't trust DBL_MAX) */ |
509 |
|
|
word0(&rv) = Big0; |
510 |
|
|
word1(&rv) = Big1; |
511 |
|
|
} |
512 |
|
|
else |
513 |
|
|
word0(&rv) += P*Exp_msk1; |
514 |
|
|
} |
515 |
|
|
} |
516 |
|
|
else if (e1 < 0) { |
517 |
|
|
e1 = -e1; |
518 |
|
|
if ( (i = e1 & 15) !=0) |
519 |
|
|
dval(&rv) /= tens[i]; |
520 |
|
|
if (e1 >>= 4) { |
521 |
|
|
if (e1 >= 1 << n_bigtens) |
522 |
|
|
goto undfl; |
523 |
|
|
#ifdef Avoid_Underflow |
524 |
|
|
if (e1 & Scale_Bit) |
525 |
|
|
scale = 2*P; |
526 |
|
|
for(j = 0; e1 > 0; j++, e1 >>= 1) |
527 |
|
|
if (e1 & 1) |
528 |
|
|
dval(&rv) *= tinytens[j]; |
529 |
|
|
if (scale && (j = 2*P + 1 - ((word0(&rv) & Exp_mask) |
530 |
|
|
>> Exp_shift)) > 0) { |
531 |
|
|
/* scaled rv is denormal; zap j low bits */ |
532 |
|
|
if (j >= 32) { |
533 |
|
|
word1(&rv) = 0; |
534 |
|
|
if (j >= 53) |
535 |
|
|
word0(&rv) = (P+2)*Exp_msk1; |
536 |
|
|
else |
537 |
|
|
word0(&rv) &= 0xffffffff << (j-32); |
538 |
|
|
} |
539 |
|
|
else |
540 |
|
|
word1(&rv) &= 0xffffffff << j; |
541 |
|
|
} |
542 |
|
|
#else |
543 |
|
|
for(j = 0; e1 > 1; j++, e1 >>= 1) |
544 |
|
|
if (e1 & 1) |
545 |
|
|
dval(&rv) *= tinytens[j]; |
546 |
|
|
/* The last multiplication could underflow. */ |
547 |
|
|
dval(&rv0) = dval(&rv); |
548 |
|
|
dval(&rv) *= tinytens[j]; |
549 |
|
|
if (!dval(&rv)) { |
550 |
|
|
dval(&rv) = 2.*dval(&rv0); |
551 |
|
|
dval(&rv) *= tinytens[j]; |
552 |
|
|
#endif |
553 |
|
|
if (!dval(&rv)) { |
554 |
|
|
undfl: |
555 |
|
|
dval(&rv) = 0.; |
556 |
|
|
goto range_err; |
557 |
|
|
} |
558 |
|
|
#ifndef Avoid_Underflow |
559 |
|
|
word0(&rv) = Tiny0; |
560 |
|
|
word1(&rv) = Tiny1; |
561 |
|
|
/* The refinement below will clean |
562 |
|
|
* this approximation up. |
563 |
|
|
*/ |
564 |
|
|
} |
565 |
|
|
#endif |
566 |
|
|
} |
567 |
|
|
} |
568 |
|
|
|
569 |
|
|
/* Now the hard part -- adjusting rv to the correct value.*/ |
570 |
|
|
|
571 |
|
|
/* Put digits into bd: true value = bd * 10^e */ |
572 |
|
|
|
573 |
|
|
bd0 = s2b(s0, nd0, nd, y, dplen); |
574 |
|
|
if (bd0 == NULL) |
575 |
|
|
goto ovfl; |
576 |
|
|
|
577 |
|
|
for(;;) { |
578 |
|
|
bd = Balloc(bd0->k); |
579 |
|
|
if (bd == NULL) |
580 |
|
|
goto ovfl; |
581 |
|
|
Bcopy(bd, bd0); |
582 |
|
|
bb = d2b(dval(&rv), &bbe, &bbbits); /* rv = bb * 2^bbe */ |
583 |
|
|
if (bb == NULL) |
584 |
|
|
goto ovfl; |
585 |
|
|
bs = i2b(1); |
586 |
|
|
if (bs == NULL) |
587 |
|
|
goto ovfl; |
588 |
|
|
|
589 |
|
|
if (e >= 0) { |
590 |
|
|
bb2 = bb5 = 0; |
591 |
|
|
bd2 = bd5 = e; |
592 |
|
|
} |
593 |
|
|
else { |
594 |
|
|
bb2 = bb5 = -e; |
595 |
|
|
bd2 = bd5 = 0; |
596 |
|
|
} |
597 |
|
|
if (bbe >= 0) |
598 |
|
|
bb2 += bbe; |
599 |
|
|
else |
600 |
|
|
bd2 -= bbe; |
601 |
|
|
bs2 = bb2; |
602 |
|
|
#ifdef Honor_FLT_ROUNDS |
603 |
|
|
if (Rounding != 1) |
604 |
|
|
bs2++; |
605 |
|
|
#endif |
606 |
|
|
#ifdef Avoid_Underflow |
607 |
|
|
Lsb = LSB; |
608 |
|
|
Lsb1 = 0; |
609 |
|
|
j = bbe - scale; |
610 |
|
|
i = j + bbbits - 1; /* logb(rv) */ |
611 |
|
|
j = P + 1 - bbbits; |
612 |
|
|
if (i < Emin) { /* denormal */ |
613 |
|
|
i = Emin - i; |
614 |
|
|
j -= i; |
615 |
|
|
if (i < 32) |
616 |
|
|
Lsb <<= i; |
617 |
|
|
else |
618 |
|
|
Lsb1 = Lsb << (i-32); |
619 |
|
|
} |
620 |
|
|
#else /*Avoid_Underflow*/ |
621 |
|
|
#ifdef Sudden_Underflow |
622 |
|
|
#ifdef IBM |
623 |
|
|
j = 1 + 4*P - 3 - bbbits + ((bbe + bbbits - 1) & 3); |
624 |
|
|
#else |
625 |
|
|
j = P + 1 - bbbits; |
626 |
|
|
#endif |
627 |
|
|
#else /*Sudden_Underflow*/ |
628 |
|
|
j = bbe; |
629 |
|
|
i = j + bbbits - 1; /* logb(&rv) */ |
630 |
|
|
if (i < Emin) /* denormal */ |
631 |
|
|
j += P - Emin; |
632 |
|
|
else |
633 |
|
|
j = P + 1 - bbbits; |
634 |
|
|
#endif /*Sudden_Underflow*/ |
635 |
|
|
#endif /*Avoid_Underflow*/ |
636 |
|
|
bb2 += j; |
637 |
|
|
bd2 += j; |
638 |
|
|
#ifdef Avoid_Underflow |
639 |
|
|
bd2 += scale; |
640 |
|
|
#endif |
641 |
|
|
i = bb2 < bd2 ? bb2 : bd2; |
642 |
|
|
if (i > bs2) |
643 |
|
|
i = bs2; |
644 |
|
|
if (i > 0) { |
645 |
|
|
bb2 -= i; |
646 |
|
|
bd2 -= i; |
647 |
|
|
bs2 -= i; |
648 |
|
|
} |
649 |
|
|
if (bb5 > 0) { |
650 |
|
|
bs = pow5mult(bs, bb5); |
651 |
|
|
if (bs == NULL) |
652 |
|
|
goto ovfl; |
653 |
|
|
bb1 = mult(bs, bb); |
654 |
|
|
if (bb1 == NULL) |
655 |
|
|
goto ovfl; |
656 |
|
|
Bfree(bb); |
657 |
|
|
bb = bb1; |
658 |
|
|
} |
659 |
|
|
if (bb2 > 0) { |
660 |
|
|
bb = lshift(bb, bb2); |
661 |
|
|
if (bb == NULL) |
662 |
|
|
goto ovfl; |
663 |
|
|
} |
664 |
|
|
if (bd5 > 0) { |
665 |
|
|
bd = pow5mult(bd, bd5); |
666 |
|
|
if (bd == NULL) |
667 |
|
|
goto ovfl; |
668 |
|
|
} |
669 |
|
|
if (bd2 > 0) { |
670 |
|
|
bd = lshift(bd, bd2); |
671 |
|
|
if (bd == NULL) |
672 |
|
|
goto ovfl; |
673 |
|
|
} |
674 |
|
|
if (bs2 > 0) { |
675 |
|
|
bs = lshift(bs, bs2); |
676 |
|
|
if (bs == NULL) |
677 |
|
|
goto ovfl; |
678 |
|
|
} |
679 |
|
|
delta = diff(bb, bd); |
680 |
|
|
if (delta == NULL) |
681 |
|
|
goto ovfl; |
682 |
|
|
dsign = delta->sign; |
683 |
|
|
delta->sign = 0; |
684 |
|
|
i = cmp(delta, bs); |
685 |
|
|
#ifdef Honor_FLT_ROUNDS |
686 |
|
|
if (Rounding != 1) { |
687 |
|
|
if (i < 0) { |
688 |
|
|
/* Error is less than an ulp */ |
689 |
|
|
if (!delta->x[0] && delta->wds <= 1) { |
690 |
|
|
/* exact */ |
691 |
|
|
#ifdef SET_INEXACT |
692 |
|
|
inexact = 0; |
693 |
|
|
#endif |
694 |
|
|
break; |
695 |
|
|
} |
696 |
|
|
if (Rounding) { |
697 |
|
|
if (dsign) { |
698 |
|
|
dval(&adj) = 1.; |
699 |
|
|
goto apply_adj; |
700 |
|
|
} |
701 |
|
|
} |
702 |
|
|
else if (!dsign) { |
703 |
|
|
dval(&adj) = -1.; |
704 |
|
|
if (!word1(&rv) |
705 |
|
|
&& !(word0(&rv) & Frac_mask)) { |
706 |
|
|
y = word0(&rv) & Exp_mask; |
707 |
|
|
#ifdef Avoid_Underflow |
708 |
|
|
if (!scale || y > 2*P*Exp_msk1) |
709 |
|
|
#else |
710 |
|
|
if (y) |
711 |
|
|
#endif |
712 |
|
|
{ |
713 |
|
|
delta = lshift(delta,Log2P); |
714 |
|
|
if (delta == NULL) |
715 |
|
|
goto ovfl; |
716 |
|
|
if (cmp(delta, bs) <= 0) |
717 |
|
|
dval(&adj) = -0.5; |
718 |
|
|
} |
719 |
|
|
} |
720 |
|
|
apply_adj: |
721 |
|
|
#ifdef Avoid_Underflow |
722 |
|
|
if (scale && (y = word0(&rv) & Exp_mask) |
723 |
|
|
<= 2*P*Exp_msk1) |
724 |
|
|
word0(&adj) += (2*P+1)*Exp_msk1 - y; |
725 |
|
|
#else |
726 |
|
|
#ifdef Sudden_Underflow |
727 |
|
|
if ((word0(&rv) & Exp_mask) <= |
728 |
|
|
P*Exp_msk1) { |
729 |
|
|
word0(&rv) += P*Exp_msk1; |
730 |
|
|
dval(&rv) += adj*ulp(&rv); |
731 |
|
|
word0(&rv) -= P*Exp_msk1; |
732 |
|
|
} |
733 |
|
|
else |
734 |
|
|
#endif /*Sudden_Underflow*/ |
735 |
|
|
#endif /*Avoid_Underflow*/ |
736 |
|
|
dval(&rv) += adj.d*ulp(&rv); |
737 |
|
|
} |
738 |
|
|
break; |
739 |
|
|
} |
740 |
|
|
dval(&adj) = ratio(delta, bs); |
741 |
|
|
if (adj.d < 1.) |
742 |
|
|
dval(&adj) = 1.; |
743 |
|
|
if (adj.d <= 0x7ffffffe) { |
744 |
|
|
/* dval(&adj) = Rounding ? ceil(&adj) : floor(&adj); */ |
745 |
|
|
y = adj.d; |
746 |
|
|
if (y != adj.d) { |
747 |
|
|
if (!((Rounding>>1) ^ dsign)) |
748 |
|
|
y++; |
749 |
|
|
dval(&adj) = y; |
750 |
|
|
} |
751 |
|
|
} |
752 |
|
|
#ifdef Avoid_Underflow |
753 |
|
|
if (scale && (y = word0(&rv) & Exp_mask) <= 2*P*Exp_msk1) |
754 |
|
|
word0(&adj) += (2*P+1)*Exp_msk1 - y; |
755 |
|
|
#else |
756 |
|
|
#ifdef Sudden_Underflow |
757 |
|
|
if ((word0(&rv) & Exp_mask) <= P*Exp_msk1) { |
758 |
|
|
word0(&rv) += P*Exp_msk1; |
759 |
|
|
dval(&adj) *= ulp(&rv); |
760 |
|
|
if (dsign) |
761 |
|
|
dval(&rv) += adj; |
762 |
|
|
else |
763 |
|
|
dval(&rv) -= adj; |
764 |
|
|
word0(&rv) -= P*Exp_msk1; |
765 |
|
|
goto cont; |
766 |
|
|
} |
767 |
|
|
#endif /*Sudden_Underflow*/ |
768 |
|
|
#endif /*Avoid_Underflow*/ |
769 |
|
|
dval(&adj) *= ulp(&rv); |
770 |
|
|
if (dsign) { |
771 |
|
|
if (word0(&rv) == Big0 && word1(&rv) == Big1) |
772 |
|
|
goto ovfl; |
773 |
|
|
dval(&rv) += adj.d; |
774 |
|
|
} |
775 |
|
|
else |
776 |
|
|
dval(&rv) -= adj.d; |
777 |
|
|
goto cont; |
778 |
|
|
} |
779 |
|
|
#endif /*Honor_FLT_ROUNDS*/ |
780 |
|
|
|
781 |
|
|
if (i < 0) { |
782 |
|
|
/* Error is less than half an ulp -- check for |
783 |
|
|
* special case of mantissa a power of two. |
784 |
|
|
*/ |
785 |
|
|
if (dsign || word1(&rv) || word0(&rv) & Bndry_mask |
786 |
|
|
#ifdef IEEE_Arith |
787 |
|
|
#ifdef Avoid_Underflow |
788 |
|
|
|| (word0(&rv) & Exp_mask) <= (2*P+1)*Exp_msk1 |
789 |
|
|
#else |
790 |
|
|
|| (word0(&rv) & Exp_mask) <= Exp_msk1 |
791 |
|
|
#endif |
792 |
|
|
#endif |
793 |
|
|
) { |
794 |
|
|
#ifdef SET_INEXACT |
795 |
|
|
if (!delta->x[0] && delta->wds <= 1) |
796 |
|
|
inexact = 0; |
797 |
|
|
#endif |
798 |
|
|
break; |
799 |
|
|
} |
800 |
|
|
if (!delta->x[0] && delta->wds <= 1) { |
801 |
|
|
/* exact result */ |
802 |
|
|
#ifdef SET_INEXACT |
803 |
|
|
inexact = 0; |
804 |
|
|
#endif |
805 |
|
|
break; |
806 |
|
|
} |
807 |
|
|
delta = lshift(delta,Log2P); |
808 |
|
|
if (delta == NULL) |
809 |
|
|
goto ovfl; |
810 |
|
|
if (cmp(delta, bs) > 0) |
811 |
|
|
goto drop_down; |
812 |
|
|
break; |
813 |
|
|
} |
814 |
|
|
if (i == 0) { |
815 |
|
|
/* exactly half-way between */ |
816 |
|
|
if (dsign) { |
817 |
|
|
if ((word0(&rv) & Bndry_mask1) == Bndry_mask1 |
818 |
|
|
&& word1(&rv) == ( |
819 |
|
|
#ifdef Avoid_Underflow |
820 |
|
|
(scale && (y = word0(&rv) & Exp_mask) <= 2*P*Exp_msk1) |
821 |
|
|
? (0xffffffff & (0xffffffff << (2*P+1-(y>>Exp_shift)))) : |
822 |
|
|
#endif |
823 |
|
|
0xffffffff)) { |
824 |
|
|
/*boundary case -- increment exponent*/ |
825 |
|
|
if (word0(&rv) == Big0 && word1(&rv) == Big1) |
826 |
|
|
goto ovfl; |
827 |
|
|
word0(&rv) = (word0(&rv) & Exp_mask) |
828 |
|
|
+ Exp_msk1 |
829 |
|
|
#ifdef IBM |
830 |
|
|
| Exp_msk1 >> 4 |
831 |
|
|
#endif |
832 |
|
|
; |
833 |
|
|
word1(&rv) = 0; |
834 |
|
|
#ifdef Avoid_Underflow |
835 |
|
|
dsign = 0; |
836 |
|
|
#endif |
837 |
|
|
break; |
838 |
|
|
} |
839 |
|
|
} |
840 |
|
|
else if (!(word0(&rv) & Bndry_mask) && !word1(&rv)) { |
841 |
|
|
drop_down: |
842 |
|
|
/* boundary case -- decrement exponent */ |
843 |
|
|
#ifdef Sudden_Underflow /*{{*/ |
844 |
|
|
L = word0(&rv) & Exp_mask; |
845 |
|
|
#ifdef IBM |
846 |
|
|
if (L < Exp_msk1) |
847 |
|
|
#else |
848 |
|
|
#ifdef Avoid_Underflow |
849 |
|
|
if (L <= (scale ? (2*P+1)*Exp_msk1 : Exp_msk1)) |
850 |
|
|
#else |
851 |
|
|
if (L <= Exp_msk1) |
852 |
|
|
#endif /*Avoid_Underflow*/ |
853 |
|
|
#endif /*IBM*/ |
854 |
|
|
goto undfl; |
855 |
|
|
L -= Exp_msk1; |
856 |
|
|
#else /*Sudden_Underflow}{*/ |
857 |
|
|
#ifdef Avoid_Underflow |
858 |
|
|
if (scale) { |
859 |
|
|
L = word0(&rv) & Exp_mask; |
860 |
|
|
if (L <= (2*P+1)*Exp_msk1) { |
861 |
|
|
if (L > (P+2)*Exp_msk1) |
862 |
|
|
/* round even ==> */ |
863 |
|
|
/* accept rv */ |
864 |
|
|
break; |
865 |
|
|
/* rv = smallest denormal */ |
866 |
|
|
goto undfl; |
867 |
|
|
} |
868 |
|
|
} |
869 |
|
|
#endif /*Avoid_Underflow*/ |
870 |
|
|
L = (word0(&rv) & Exp_mask) - Exp_msk1; |
871 |
|
|
#endif /*Sudden_Underflow}}*/ |
872 |
|
|
word0(&rv) = L | Bndry_mask1; |
873 |
|
|
word1(&rv) = 0xffffffff; |
874 |
|
|
#ifdef IBM |
875 |
|
|
goto cont; |
876 |
|
|
#else |
877 |
|
|
break; |
878 |
|
|
#endif |
879 |
|
|
} |
880 |
|
|
#ifndef ROUND_BIASED |
881 |
|
|
#ifdef Avoid_Underflow |
882 |
|
|
if (Lsb1) { |
883 |
|
|
if (!(word0(&rv) & Lsb1)) |
884 |
|
|
break; |
885 |
|
|
} |
886 |
|
|
else if (!(word1(&rv) & Lsb)) |
887 |
|
|
break; |
888 |
|
|
#else |
889 |
|
|
if (!(word1(&rv) & LSB)) |
890 |
|
|
break; |
891 |
|
|
#endif |
892 |
|
|
#endif |
893 |
|
|
if (dsign) |
894 |
|
|
#ifdef Avoid_Underflow |
895 |
|
|
dval(&rv) += sulp(&rv, scale); |
896 |
|
|
#else |
897 |
|
|
dval(&rv) += ulp(&rv); |
898 |
|
|
#endif |
899 |
|
|
#ifndef ROUND_BIASED |
900 |
|
|
else { |
901 |
|
|
#ifdef Avoid_Underflow |
902 |
|
|
dval(&rv) -= sulp(&rv, scale); |
903 |
|
|
#else |
904 |
|
|
dval(&rv) -= ulp(&rv); |
905 |
|
|
#endif |
906 |
|
|
#ifndef Sudden_Underflow |
907 |
|
|
if (!dval(&rv)) |
908 |
|
|
goto undfl; |
909 |
|
|
#endif |
910 |
|
|
} |
911 |
|
|
#ifdef Avoid_Underflow |
912 |
|
|
dsign = 1 - dsign; |
913 |
|
|
#endif |
914 |
|
|
#endif |
915 |
|
|
break; |
916 |
|
|
} |
917 |
|
|
if ((aadj = ratio(delta, bs)) <= 2.) { |
918 |
|
|
if (dsign) |
919 |
|
|
aadj = dval(&aadj1) = 1.; |
920 |
|
|
else if (word1(&rv) || word0(&rv) & Bndry_mask) { |
921 |
|
|
#ifndef Sudden_Underflow |
922 |
|
|
if (word1(&rv) == Tiny1 && !word0(&rv)) |
923 |
|
|
goto undfl; |
924 |
|
|
#endif |
925 |
|
|
aadj = 1.; |
926 |
|
|
dval(&aadj1) = -1.; |
927 |
|
|
} |
928 |
|
|
else { |
929 |
|
|
/* special case -- power of FLT_RADIX to be */ |
930 |
|
|
/* rounded down... */ |
931 |
|
|
|
932 |
|
|
if (aadj < 2./FLT_RADIX) |
933 |
|
|
aadj = 1./FLT_RADIX; |
934 |
|
|
else |
935 |
|
|
aadj *= 0.5; |
936 |
|
|
dval(&aadj1) = -aadj; |
937 |
|
|
} |
938 |
|
|
} |
939 |
|
|
else { |
940 |
|
|
aadj *= 0.5; |
941 |
|
|
dval(&aadj1) = dsign ? aadj : -aadj; |
942 |
|
|
#ifdef Check_FLT_ROUNDS |
943 |
|
|
switch(Rounding) { |
944 |
|
|
case 2: /* towards +infinity */ |
945 |
|
|
dval(&aadj1) -= 0.5; |
946 |
|
|
break; |
947 |
|
|
case 0: /* towards 0 */ |
948 |
|
|
case 3: /* towards -infinity */ |
949 |
|
|
dval(&aadj1) += 0.5; |
950 |
|
|
} |
951 |
|
|
#else |
952 |
|
|
if (Flt_Rounds == 0) |
953 |
|
|
dval(&aadj1) += 0.5; |
954 |
|
|
#endif /*Check_FLT_ROUNDS*/ |
955 |
|
|
} |
956 |
|
|
y = word0(&rv) & Exp_mask; |
957 |
|
|
|
958 |
|
|
/* Check for overflow */ |
959 |
|
|
|
960 |
|
|
if (y == Exp_msk1*(DBL_MAX_EXP+Bias-1)) { |
961 |
|
|
dval(&rv0) = dval(&rv); |
962 |
|
|
word0(&rv) -= P*Exp_msk1; |
963 |
|
|
dval(&adj) = dval(&aadj1) * ulp(&rv); |
964 |
|
|
dval(&rv) += dval(&adj); |
965 |
|
|
if ((word0(&rv) & Exp_mask) >= |
966 |
|
|
Exp_msk1*(DBL_MAX_EXP+Bias-P)) { |
967 |
|
|
if (word0(&rv0) == Big0 && word1(&rv0) == Big1) |
968 |
|
|
goto ovfl; |
969 |
|
|
word0(&rv) = Big0; |
970 |
|
|
word1(&rv) = Big1; |
971 |
|
|
goto cont; |
972 |
|
|
} |
973 |
|
|
else |
974 |
|
|
word0(&rv) += P*Exp_msk1; |
975 |
|
|
} |
976 |
|
|
else { |
977 |
|
|
#ifdef Avoid_Underflow |
978 |
|
|
if (scale && y <= 2*P*Exp_msk1) { |
979 |
|
|
if (aadj <= 0x7fffffff) { |
980 |
|
|
if ((z = aadj) <= 0) |
981 |
|
|
z = 1; |
982 |
|
|
aadj = z; |
983 |
|
|
dval(&aadj1) = dsign ? aadj : -aadj; |
984 |
|
|
} |
985 |
|
|
word0(&aadj1) += (2*P+1)*Exp_msk1 - y; |
986 |
|
|
} |
987 |
|
|
dval(&adj) = dval(&aadj1) * ulp(&rv); |
988 |
|
|
dval(&rv) += dval(&adj); |
989 |
|
|
#else |
990 |
|
|
#ifdef Sudden_Underflow |
991 |
|
|
if ((word0(&rv) & Exp_mask) <= P*Exp_msk1) { |
992 |
|
|
dval(&rv0) = dval(&rv); |
993 |
|
|
word0(&rv) += P*Exp_msk1; |
994 |
|
|
dval(&adj) = dval(&aadj1) * ulp(&rv); |
995 |
|
|
dval(&rv) += dval(&adj); |
996 |
|
|
#ifdef IBM |
997 |
|
|
if ((word0(&rv) & Exp_mask) < P*Exp_msk1) |
998 |
|
|
#else |
999 |
|
|
if ((word0(&rv) & Exp_mask) <= P*Exp_msk1) |
1000 |
|
|
#endif |
1001 |
|
|
{ |
1002 |
|
|
if (word0(&rv0) == Tiny0 |
1003 |
|
|
&& word1(&rv0) == Tiny1) |
1004 |
|
|
goto undfl; |
1005 |
|
|
word0(&rv) = Tiny0; |
1006 |
|
|
word1(&rv) = Tiny1; |
1007 |
|
|
goto cont; |
1008 |
|
|
} |
1009 |
|
|
else |
1010 |
|
|
word0(&rv) -= P*Exp_msk1; |
1011 |
|
|
} |
1012 |
|
|
else { |
1013 |
|
|
dval(&adj) = dval(&aadj1) * ulp(&rv); |
1014 |
|
|
dval(&rv) += dval(&adj); |
1015 |
|
|
} |
1016 |
|
|
#else /*Sudden_Underflow*/ |
1017 |
|
|
/* Compute dval(&adj) so that the IEEE rounding rules will |
1018 |
|
|
* correctly round rv + dval(&adj) in some half-way cases. |
1019 |
|
|
* If rv * ulp(&rv) is denormalized (i.e., |
1020 |
|
|
* y <= (P-1)*Exp_msk1), we must adjust aadj to avoid |
1021 |
|
|
* trouble from bits lost to denormalization; |
1022 |
|
|
* example: 1.2e-307 . |
1023 |
|
|
*/ |
1024 |
|
|
if (y <= (P-1)*Exp_msk1 && aadj > 1.) { |
1025 |
|
|
dval(&aadj1) = (double)(int)(aadj + 0.5); |
1026 |
|
|
if (!dsign) |
1027 |
|
|
dval(&aadj1) = -dval(&aadj1); |
1028 |
|
|
} |
1029 |
|
|
dval(&adj) = dval(&aadj1) * ulp(&rv); |
1030 |
|
|
dval(&rv) += adj; |
1031 |
|
|
#endif /*Sudden_Underflow*/ |
1032 |
|
|
#endif /*Avoid_Underflow*/ |
1033 |
|
|
} |
1034 |
|
|
z = word0(&rv) & Exp_mask; |
1035 |
|
|
#ifndef SET_INEXACT |
1036 |
|
|
#ifdef Avoid_Underflow |
1037 |
|
|
if (!scale) |
1038 |
|
|
#endif |
1039 |
|
|
if (y == z) { |
1040 |
|
|
/* Can we stop now? */ |
1041 |
|
|
L = (Long)aadj; |
1042 |
|
|
aadj -= L; |
1043 |
|
|
/* The tolerances below are conservative. */ |
1044 |
|
|
if (dsign || word1(&rv) || word0(&rv) & Bndry_mask) { |
1045 |
|
|
if (aadj < .4999999 || aadj > .5000001) |
1046 |
|
|
break; |
1047 |
|
|
} |
1048 |
|
|
else if (aadj < .4999999/FLT_RADIX) |
1049 |
|
|
break; |
1050 |
|
|
} |
1051 |
|
|
#endif |
1052 |
|
|
cont: |
1053 |
|
|
Bfree(bb); |
1054 |
|
|
Bfree(bd); |
1055 |
|
|
Bfree(bs); |
1056 |
|
|
Bfree(delta); |
1057 |
|
|
} |
1058 |
|
|
Bfree(bb); |
1059 |
|
|
Bfree(bd); |
1060 |
|
|
Bfree(bs); |
1061 |
|
|
Bfree(bd0); |
1062 |
|
|
Bfree(delta); |
1063 |
|
|
#ifdef SET_INEXACT |
1064 |
|
|
if (inexact) { |
1065 |
|
|
if (!oldinexact) { |
1066 |
|
|
word0(&rv0) = Exp_1 + (70 << Exp_shift); |
1067 |
|
|
word1(&rv0) = 0; |
1068 |
|
|
dval(&rv0) += 1.; |
1069 |
|
|
} |
1070 |
|
|
} |
1071 |
|
|
else if (!oldinexact) |
1072 |
|
|
clear_inexact(); |
1073 |
|
|
#endif |
1074 |
|
|
#ifdef Avoid_Underflow |
1075 |
|
|
if (scale) { |
1076 |
|
|
word0(&rv0) = Exp_1 - 2*P*Exp_msk1; |
1077 |
|
|
word1(&rv0) = 0; |
1078 |
|
|
dval(&rv) *= dval(&rv0); |
1079 |
|
|
#ifndef NO_ERRNO |
1080 |
|
|
/* try to avoid the bug of testing an 8087 register value */ |
1081 |
|
|
#ifdef IEEE_Arith |
1082 |
|
|
if (!(word0(&rv) & Exp_mask)) |
1083 |
|
|
#else |
1084 |
|
|
if (word0(&rv) == 0 && word1(&rv) == 0) |
1085 |
|
|
#endif |
1086 |
|
|
errno = ERANGE; |
1087 |
|
|
#endif |
1088 |
|
|
} |
1089 |
|
|
#endif /* Avoid_Underflow */ |
1090 |
|
|
#ifdef SET_INEXACT |
1091 |
|
|
if (inexact && !(word0(&rv) & Exp_mask)) { |
1092 |
|
|
/* set underflow bit */ |
1093 |
|
|
dval(&rv0) = 1e-300; |
1094 |
|
|
dval(&rv0) *= dval(&rv0); |
1095 |
|
|
} |
1096 |
|
|
#endif |
1097 |
|
|
ret: |
1098 |
|
|
if (se) |
1099 |
|
|
*se = (char *)s; |
1100 |
|
|
return sign ? -dval(&rv) : dval(&rv); |
1101 |
|
|
} |
1102 |
|
|
DEF_STRONG(strtod); |