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