1 /* $NetBSD: strtodg.c,v 1.6 2008/03/21 23:13:48 christos Exp $ */ 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 34 #include "gdtoaimp.h" 35 36 #ifdef USE_LOCALE 37 #include "locale.h" 38 #endif 39 40 static CONST int 41 fivesbits[] = { 0, 3, 5, 7, 10, 12, 14, 17, 19, 21, 42 24, 26, 28, 31, 33, 35, 38, 40, 42, 45, 43 47, 49, 52 44 #ifdef VAX 45 , 54, 56 46 #endif 47 }; 48 49 Bigint * 50 #ifdef KR_headers 51 increment(b) Bigint *b; 52 #else 53 increment(Bigint *b) 54 #endif 55 { 56 ULong *x, *xe; 57 Bigint *b1; 58 #ifdef Pack_16 59 ULong carry = 1, y; 60 #endif 61 62 x = b->x; 63 xe = x + b->wds; 64 #ifdef Pack_32 65 do { 66 if (*x < (ULong)0xffffffffL) { 67 ++*x; 68 return b; 69 } 70 *x++ = 0; 71 } while(x < xe); 72 #else 73 do { 74 y = *x + carry; 75 carry = y >> 16; 76 *x++ = y & 0xffff; 77 if (!carry) 78 return b; 79 } while(x < xe); 80 if (carry) 81 #endif 82 { 83 if (b->wds >= b->maxwds) { 84 b1 = Balloc(b->k+1); 85 if (b1 == NULL) 86 return NULL; 87 Bcopy(b1,b); 88 Bfree(b); 89 b = b1; 90 } 91 b->x[b->wds++] = 1; 92 } 93 return b; 94 } 95 96 int 97 #ifdef KR_headers 98 decrement(b) Bigint *b; 99 #else 100 decrement(Bigint *b) 101 #endif 102 { 103 ULong *x, *xe; 104 #ifdef Pack_16 105 ULong borrow = 1, y; 106 #endif 107 108 x = b->x; 109 xe = x + b->wds; 110 #ifdef Pack_32 111 do { 112 if (*x) { 113 --*x; 114 break; 115 } 116 *x++ = 0xffffffffUL; 117 } 118 while(x < xe); 119 #else 120 do { 121 y = *x - borrow; 122 borrow = (y & 0x10000) >> 16; 123 *x++ = y & 0xffff; 124 } while(borrow && x < xe); 125 #endif 126 return STRTOG_Inexlo; 127 } 128 129 static int 130 #ifdef KR_headers 131 all_on(b, n) CONST Bigint *b; int n; 132 #else 133 all_on(CONST Bigint *b, int n) 134 #endif 135 { 136 CONST ULong *x, *xe; 137 138 x = b->x; 139 xe = x + ((unsigned int)n >> kshift); 140 while(x < xe) 141 if ((*x++ & ALL_ON) != ALL_ON) 142 return 0; 143 if (n &= kmask) 144 return ((*x | (ALL_ON << n)) & ALL_ON) == ALL_ON; 145 return 1; 146 } 147 148 Bigint * 149 #ifdef KR_headers 150 set_ones(b, n) Bigint *b; int n; 151 #else 152 set_ones(Bigint *b, int n) 153 #endif 154 { 155 int k; 156 ULong *x, *xe; 157 158 k = (unsigned int)(n + ((1 << kshift) - 1)) >> kshift; 159 if (b->k < k) { 160 Bfree(b); 161 b = Balloc(k); 162 if (b == NULL) 163 return NULL; 164 } 165 k = (unsigned int)n >> kshift; 166 if (n &= kmask) 167 k++; 168 b->wds = k; 169 x = b->x; 170 xe = x + k; 171 while(x < xe) 172 *x++ = ALL_ON; 173 if (n) 174 x[-1] >>= ULbits - n; 175 return b; 176 } 177 178 static int 179 rvOK 180 #ifdef KR_headers 181 (d, fpi, expt, bits, exact, rd, irv) 182 double d; CONST FPI *fpi; Long *expt; ULong *bits; int exact, rd, *irv; 183 #else 184 (double d, CONST FPI *fpi, Long *expt, ULong *bits, int exact, int rd, int *irv) 185 #endif 186 { 187 Bigint *b; 188 ULong carry, inex, lostbits; 189 int bdif, e, j, k, k1, nb, rv; 190 191 carry = rv = 0; 192 b = d2b(d, &e, &bdif); 193 bdif -= nb = fpi->nbits; 194 e += bdif; 195 if (bdif <= 0) { 196 if (exact) 197 goto trunc; 198 goto ret; 199 } 200 if (P == nb) { 201 if ( 202 #ifndef IMPRECISE_INEXACT 203 exact && 204 #endif 205 fpi->rounding == 206 #ifdef RND_PRODQUOT 207 FPI_Round_near 208 #else 209 Flt_Rounds 210 #endif 211 ) goto trunc; 212 goto ret; 213 } 214 switch(rd) { 215 case 1: 216 goto trunc; 217 case 2: 218 break; 219 default: /* round near */ 220 k = bdif - 1; 221 if (!k) { 222 if (!exact) 223 goto ret; 224 if (b->x[0] & 2) 225 break; 226 goto trunc; 227 } 228 if (b->x[(unsigned int)k>>kshift] & ((ULong)1 << (k & kmask))) 229 break; 230 goto trunc; 231 } 232 /* "break" cases: round up 1 bit, then truncate; bdif > 0 */ 233 carry = 1; 234 trunc: 235 inex = lostbits = 0; 236 if (bdif > 0) { 237 if ( (lostbits = any_on(b, bdif)) !=0) 238 inex = STRTOG_Inexlo; 239 rshift(b, bdif); 240 if (carry) { 241 inex = STRTOG_Inexhi; 242 b = increment(b); 243 if ( (j = nb & kmask) !=0) 244 j = ULbits - j; 245 if (hi0bits(b->x[b->wds - 1]) != j) { 246 if (!lostbits) 247 lostbits = b->x[0] & 1; 248 rshift(b, 1); 249 e++; 250 } 251 } 252 } 253 else if (bdif < 0) 254 b = lshift(b, -bdif); 255 if (e < fpi->emin) { 256 k = fpi->emin - e; 257 e = fpi->emin; 258 if (k > nb || fpi->sudden_underflow) { 259 b->wds = inex = 0; 260 *irv = STRTOG_Underflow | STRTOG_Inexlo; 261 } 262 else { 263 k1 = k - 1; 264 if (k1 > 0 && !lostbits) 265 lostbits = any_on(b, k1); 266 if (!lostbits && !exact) 267 goto ret; 268 lostbits |= 269 carry = b->x[(unsigned int)k1>>kshift] & 270 (1 << (k1 & kmask)); 271 rshift(b, k); 272 *irv = STRTOG_Denormal; 273 if (carry) { 274 b = increment(b); 275 inex = STRTOG_Inexhi | STRTOG_Underflow; 276 } 277 else if (lostbits) 278 inex = STRTOG_Inexlo | STRTOG_Underflow; 279 } 280 } 281 else if (e > fpi->emax) { 282 e = fpi->emax + 1; 283 *irv = STRTOG_Infinite | STRTOG_Overflow | STRTOG_Inexhi; 284 #ifndef NO_ERRNO 285 errno = ERANGE; 286 #endif 287 b->wds = inex = 0; 288 } 289 *expt = e; 290 copybits(bits, nb, b); 291 *irv |= inex; 292 rv = 1; 293 ret: 294 Bfree(b); 295 return rv; 296 } 297 298 #ifndef VAX 299 static int 300 #ifdef KR_headers 301 mantbits(d) double d; 302 #else 303 mantbits(double d) 304 #endif 305 { 306 ULong L; 307 #ifdef VAX 308 L = word1(d) << 16 | word1(d) >> 16; 309 if (L) 310 #else 311 if ( (L = word1(d)) !=0) 312 #endif 313 return P - lo0bits(&L); 314 #ifdef VAX 315 L = word0(d) << 16 | word0(d) >> 16 | Exp_msk11; 316 #else 317 L = word0(d) | Exp_msk1; 318 #endif 319 return P - 32 - lo0bits(&L); 320 } 321 #endif /* !VAX */ 322 323 int 324 strtodg 325 #ifdef KR_headers 326 (s00, se, fpi, expt, bits) 327 CONST char *s00; char **se; CONST FPI *fpi; Long *expt; ULong *bits; 328 #else 329 (CONST char *s00, char **se, CONST FPI *fpi, Long *expt, ULong *bits) 330 #endif 331 { 332 int abe, abits, asub; 333 int bb0, bb2, bb5, bbe, bd2, bd5, bbbits, bs2, c, decpt, denorm; 334 int dsign, e, e1, e2, emin, esign, finished, i, inex, irv; 335 int j, k, nbits, nd, nd0, nf, nz, nz0, rd, rvbits, rve, rve1, sign; 336 int sudden_underflow = 0; /* pacify gcc */ 337 CONST char *s, *s0, *s1; 338 double adj, adj0, rv, tol; 339 Long L; 340 ULong y, z; 341 Bigint *ab, *bb, *bb1, *bd, *bd0, *bs, *delta, *rvb, *rvb0; 342 343 e2 = 0; /* XXX gcc */ 344 345 irv = STRTOG_Zero; 346 denorm = sign = nz0 = nz = 0; 347 dval(rv) = 0.; 348 rvb = 0; 349 nbits = fpi->nbits; 350 for(s = s00;;s++) switch(*s) { 351 case '-': 352 sign = 1; 353 /* FALLTHROUGH */ 354 case '+': 355 if (*++s) 356 goto break2; 357 /* FALLTHROUGH */ 358 case 0: 359 sign = 0; 360 irv = STRTOG_NoNumber; 361 s = s00; 362 goto ret; 363 case '\t': 364 case '\n': 365 case '\v': 366 case '\f': 367 case '\r': 368 case ' ': 369 continue; 370 default: 371 goto break2; 372 } 373 break2: 374 if (*s == '0') { 375 #ifndef NO_HEX_FP 376 switch(s[1]) { 377 case 'x': 378 case 'X': 379 irv = gethex(&s, fpi, expt, &rvb, sign); 380 if (irv == STRTOG_NoNumber) { 381 s = s00; 382 sign = 0; 383 } 384 goto ret; 385 } 386 #endif 387 nz0 = 1; 388 while(*++s == '0') ; 389 if (!*s) 390 goto ret; 391 } 392 sudden_underflow = fpi->sudden_underflow; 393 s0 = s; 394 y = z = 0; 395 for(decpt = nd = nf = 0; (c = *s) >= '0' && c <= '9'; nd++, s++) 396 if (nd < 9) 397 y = 10*y + c - '0'; 398 else if (nd < 16) 399 z = 10*z + c - '0'; 400 nd0 = nd; 401 #ifdef USE_LOCALE 402 if (c == *localeconv()->decimal_point) 403 #else 404 if (c == '.') 405 #endif 406 { 407 decpt = 1; 408 c = *++s; 409 if (!nd) { 410 for(; c == '0'; c = *++s) 411 nz++; 412 if (c > '0' && c <= '9') { 413 s0 = s; 414 nf += nz; 415 nz = 0; 416 goto have_dig; 417 } 418 goto dig_done; 419 } 420 for(; c >= '0' && c <= '9'; c = *++s) { 421 have_dig: 422 nz++; 423 if (c -= '0') { 424 nf += nz; 425 for(i = 1; i < nz; i++) 426 if (nd++ < 9) 427 y *= 10; 428 else if (nd <= DBL_DIG + 1) 429 z *= 10; 430 if (nd++ < 9) 431 y = 10*y + c; 432 else if (nd <= DBL_DIG + 1) 433 z = 10*z + c; 434 nz = 0; 435 } 436 } 437 } 438 dig_done: 439 e = 0; 440 if (c == 'e' || c == 'E') { 441 if (!nd && !nz && !nz0) { 442 irv = STRTOG_NoNumber; 443 s = s00; 444 goto ret; 445 } 446 s00 = s; 447 esign = 0; 448 switch(c = *++s) { 449 case '-': 450 esign = 1; 451 /* FALLTHROUGH */ 452 case '+': 453 c = *++s; 454 } 455 if (c >= '0' && c <= '9') { 456 while(c == '0') 457 c = *++s; 458 if (c > '0' && c <= '9') { 459 L = c - '0'; 460 s1 = s; 461 while((c = *++s) >= '0' && c <= '9') 462 L = 10*L + c - '0'; 463 if (s - s1 > 8 || L > 19999) 464 /* Avoid confusion from exponents 465 * so large that e might overflow. 466 */ 467 e = 19999; /* safe for 16 bit ints */ 468 else 469 e = (int)L; 470 if (esign) 471 e = -e; 472 } 473 else 474 e = 0; 475 } 476 else 477 s = s00; 478 } 479 if (!nd) { 480 if (!nz && !nz0) { 481 #ifdef INFNAN_CHECK 482 /* Check for Nan and Infinity */ 483 if (!decpt) 484 switch(c) { 485 case 'i': 486 case 'I': 487 if (match(&s,"nf")) { 488 --s; 489 if (!match(&s,"inity")) 490 ++s; 491 irv = STRTOG_Infinite; 492 goto infnanexp; 493 } 494 break; 495 case 'n': 496 case 'N': 497 if (match(&s, "an")) { 498 irv = STRTOG_NaN; 499 *expt = fpi->emax + 1; 500 #ifndef No_Hex_NaN 501 if (*s == '(') /*)*/ 502 irv = hexnan(&s, fpi, bits); 503 #endif 504 goto infnanexp; 505 } 506 } 507 #endif /* INFNAN_CHECK */ 508 irv = STRTOG_NoNumber; 509 s = s00; 510 } 511 goto ret; 512 } 513 514 irv = STRTOG_Normal; 515 e1 = e -= nf; 516 rd = 0; 517 switch(fpi->rounding & 3) { 518 case FPI_Round_up: 519 rd = 2 - sign; 520 break; 521 case FPI_Round_zero: 522 rd = 1; 523 break; 524 case FPI_Round_down: 525 rd = 1 + sign; 526 } 527 528 /* Now we have nd0 digits, starting at s0, followed by a 529 * decimal point, followed by nd-nd0 digits. The number we're 530 * after is the integer represented by those digits times 531 * 10**e */ 532 533 if (!nd0) 534 nd0 = nd; 535 k = nd < DBL_DIG + 1 ? nd : DBL_DIG + 1; 536 dval(rv) = y; 537 if (k > 9) 538 dval(rv) = tens[k - 9] * dval(rv) + z; 539 bd0 = 0; 540 if (nbits <= P && nd <= DBL_DIG) { 541 if (!e) { 542 if (rvOK(dval(rv), fpi, expt, bits, 1, rd, &irv)) 543 goto ret; 544 } 545 else if (e > 0) { 546 if (e <= Ten_pmax) { 547 #ifdef VAX 548 goto vax_ovfl_check; 549 #else 550 i = fivesbits[e] + mantbits(dval(rv)) <= P; 551 /* rv = */ rounded_product(dval(rv), tens[e]); 552 if (rvOK(dval(rv), fpi, expt, bits, i, rd, &irv)) 553 goto ret; 554 e1 -= e; 555 goto rv_notOK; 556 #endif 557 } 558 i = DBL_DIG - nd; 559 if (e <= Ten_pmax + i) { 560 /* A fancier test would sometimes let us do 561 * this for larger i values. 562 */ 563 e2 = e - i; 564 e1 -= i; 565 dval(rv) *= tens[i]; 566 #ifdef VAX 567 /* VAX exponent range is so narrow we must 568 * worry about overflow here... 569 */ 570 vax_ovfl_check: 571 dval(adj) = dval(rv); 572 word0(adj) -= P*Exp_msk1; 573 /* adj = */ rounded_product(dval(adj), tens[e2]); 574 if ((word0(adj) & Exp_mask) 575 > Exp_msk1*(DBL_MAX_EXP+Bias-1-P)) 576 goto rv_notOK; 577 word0(adj) += P*Exp_msk1; 578 dval(rv) = dval(adj); 579 #else 580 /* rv = */ rounded_product(dval(rv), tens[e2]); 581 #endif 582 if (rvOK(dval(rv), fpi, expt, bits, 0, rd, &irv)) 583 goto ret; 584 e1 -= e2; 585 } 586 } 587 #ifndef Inaccurate_Divide 588 else if (e >= -Ten_pmax) { 589 /* rv = */ rounded_quotient(dval(rv), tens[-e]); 590 if (rvOK(dval(rv), fpi, expt, bits, 0, rd, &irv)) 591 goto ret; 592 e1 -= e; 593 } 594 #endif 595 } 596 rv_notOK: 597 e1 += nd - k; 598 599 /* Get starting approximation = rv * 10**e1 */ 600 601 e2 = 0; 602 if (e1 > 0) { 603 if ( (i = e1 & 15) !=0) 604 dval(rv) *= tens[i]; 605 if (e1 &= ~15) { 606 e1 = (unsigned int)e1 >> 4; 607 while(e1 >= (1 << (n_bigtens-1))) { 608 e2 += ((word0(rv) & Exp_mask) 609 >> Exp_shift1) - Bias; 610 word0(rv) &= ~Exp_mask; 611 word0(rv) |= Bias << Exp_shift1; 612 dval(rv) *= bigtens[n_bigtens-1]; 613 e1 -= 1 << (n_bigtens-1); 614 } 615 e2 += ((word0(rv) & Exp_mask) >> Exp_shift1) - Bias; 616 word0(rv) &= ~Exp_mask; 617 word0(rv) |= Bias << Exp_shift1; 618 for(j = 0; e1 > 0; j++, e1 = (unsigned int)e1 >> 1) 619 if (e1 & 1) 620 dval(rv) *= bigtens[j]; 621 } 622 } 623 else if (e1 < 0) { 624 e1 = -e1; 625 if ( (i = e1 & 15) !=0) 626 dval(rv) /= tens[i]; 627 if (e1 &= ~15) { 628 e1 = (unsigned int)e1 >> 4; 629 while(e1 >= (1 << (n_bigtens-1))) { 630 e2 += ((word0(rv) & Exp_mask) 631 >> Exp_shift1) - Bias; 632 word0(rv) &= ~Exp_mask; 633 word0(rv) |= Bias << Exp_shift1; 634 dval(rv) *= tinytens[n_bigtens-1]; 635 e1 -= 1 << (n_bigtens-1); 636 } 637 e2 += ((word0(rv) & Exp_mask) >> Exp_shift1) - Bias; 638 word0(rv) &= ~Exp_mask; 639 word0(rv) |= Bias << Exp_shift1; 640 for(j = 0; e1 > 0; j++, e1 = (unsigned int)e1 >> 1) 641 if (e1 & 1) 642 dval(rv) *= tinytens[j]; 643 } 644 } 645 #ifdef IBM 646 /* e2 is a correction to the (base 2) exponent of the return 647 * value, reflecting adjustments above to avoid overflow in the 648 * native arithmetic. For native IBM (base 16) arithmetic, we 649 * must multiply e2 by 4 to change from base 16 to 2. 650 */ 651 e2 <<= 2; 652 #endif 653 rvb = d2b(dval(rv), &rve, &rvbits); /* rv = rvb * 2^rve */ 654 if (rvb == NULL) 655 return STRTOG_NoMemory; 656 rve += e2; 657 if ((j = rvbits - nbits) > 0) { 658 rshift(rvb, j); 659 rvbits = nbits; 660 rve += j; 661 } 662 bb0 = 0; /* trailing zero bits in rvb */ 663 e2 = rve + rvbits - nbits; 664 if (e2 > fpi->emax + 1) 665 goto huge; 666 rve1 = rve + rvbits - nbits; 667 if (e2 < (emin = fpi->emin)) { 668 denorm = 1; 669 j = rve - emin; 670 if (j > 0) { 671 rvb = lshift(rvb, j); 672 rvbits += j; 673 } 674 else if (j < 0) { 675 rvbits += j; 676 if (rvbits <= 0) { 677 if (rvbits < -1) { 678 ufl: 679 rvb->wds = 0; 680 rvb->x[0] = 0; 681 *expt = emin; 682 irv = STRTOG_Underflow | STRTOG_Inexlo; 683 goto ret; 684 } 685 rvb->x[0] = rvb->wds = rvbits = 1; 686 } 687 else 688 rshift(rvb, -j); 689 } 690 rve = rve1 = emin; 691 if (sudden_underflow && e2 + 1 < emin) 692 goto ufl; 693 } 694 695 /* Now the hard part -- adjusting rv to the correct value.*/ 696 697 /* Put digits into bd: true value = bd * 10^e */ 698 699 bd0 = s2b(s0, nd0, nd, y); 700 701 for(;;) { 702 bd = Balloc(bd0->k); 703 if (bd == NULL) 704 return STRTOG_NoMemory; 705 Bcopy(bd, bd0); 706 bb = Balloc(rvb->k); 707 if (bb == NULL) 708 return STRTOG_NoMemory; 709 Bcopy(bb, rvb); 710 bbbits = rvbits - bb0; 711 bbe = rve + bb0; 712 bs = i2b(1); 713 if (bs == NULL) 714 return STRTOG_NoMemory; 715 716 if (e >= 0) { 717 bb2 = bb5 = 0; 718 bd2 = bd5 = e; 719 } 720 else { 721 bb2 = bb5 = -e; 722 bd2 = bd5 = 0; 723 } 724 if (bbe >= 0) 725 bb2 += bbe; 726 else 727 bd2 -= bbe; 728 bs2 = bb2; 729 j = nbits + 1 - bbbits; 730 i = bbe + bbbits - nbits; 731 if (i < emin) /* denormal */ 732 j += i - emin; 733 bb2 += j; 734 bd2 += j; 735 i = bb2 < bd2 ? bb2 : bd2; 736 if (i > bs2) 737 i = bs2; 738 if (i > 0) { 739 bb2 -= i; 740 bd2 -= i; 741 bs2 -= i; 742 } 743 if (bb5 > 0) { 744 bs = pow5mult(bs, bb5); 745 if (bs == NULL) 746 return STRTOG_NoMemory; 747 bb1 = mult(bs, bb); 748 if (bb1 == NULL) 749 return STRTOG_NoMemory; 750 Bfree(bb); 751 bb = bb1; 752 } 753 bb2 -= bb0; 754 if (bb2 > 0) { 755 bb = lshift(bb, bb2); 756 if (bb == NULL) 757 return STRTOG_NoMemory; 758 } 759 else if (bb2 < 0) 760 rshift(bb, -bb2); 761 if (bd5 > 0) { 762 bd = pow5mult(bd, bd5); 763 if (bd == NULL) 764 return STRTOG_NoMemory; 765 } 766 if (bd2 > 0) { 767 bd = lshift(bd, bd2); 768 if (bd == NULL) 769 return STRTOG_NoMemory; 770 } 771 if (bs2 > 0) { 772 bs = lshift(bs, bs2); 773 if (bs == NULL) 774 return STRTOG_NoMemory; 775 } 776 asub = 1; 777 inex = STRTOG_Inexhi; 778 delta = diff(bb, bd); 779 if (delta == NULL) 780 return STRTOG_NoMemory; 781 if (delta->wds <= 1 && !delta->x[0]) 782 break; 783 dsign = delta->sign; 784 delta->sign = finished = 0; 785 L = 0; 786 i = cmp(delta, bs); 787 if (rd && i <= 0) { 788 irv = STRTOG_Normal; 789 if ( (finished = dsign ^ (rd&1)) !=0) { 790 if (dsign != 0) { 791 irv |= STRTOG_Inexhi; 792 goto adj1; 793 } 794 irv |= STRTOG_Inexlo; 795 if (rve1 == emin) 796 goto adj1; 797 for(i = 0, j = nbits; j >= ULbits; 798 i++, j -= ULbits) { 799 if (rvb->x[i] & ALL_ON) 800 goto adj1; 801 } 802 if (j > 1 && lo0bits(rvb->x + i) < j - 1) 803 goto adj1; 804 rve = rve1 - 1; 805 rvb = set_ones(rvb, rvbits = nbits); 806 if (rvb == NULL) 807 return STRTOG_NoMemory; 808 break; 809 } 810 irv |= dsign ? STRTOG_Inexlo : STRTOG_Inexhi; 811 break; 812 } 813 if (i < 0) { 814 /* Error is less than half an ulp -- check for 815 * special case of mantissa a power of two. 816 */ 817 irv = dsign 818 ? STRTOG_Normal | STRTOG_Inexlo 819 : STRTOG_Normal | STRTOG_Inexhi; 820 if (dsign || bbbits > 1 || denorm || rve1 == emin) 821 break; 822 delta = lshift(delta,1); 823 if (delta == NULL) 824 return STRTOG_NoMemory; 825 if (cmp(delta, bs) > 0) { 826 irv = STRTOG_Normal | STRTOG_Inexlo; 827 goto drop_down; 828 } 829 break; 830 } 831 if (i == 0) { 832 /* exactly half-way between */ 833 if (dsign) { 834 if (denorm && all_on(rvb, rvbits)) { 835 /*boundary case -- increment exponent*/ 836 rvb->wds = 1; 837 rvb->x[0] = 1; 838 rve = emin + nbits - (rvbits = 1); 839 irv = STRTOG_Normal | STRTOG_Inexhi; 840 denorm = 0; 841 break; 842 } 843 irv = STRTOG_Normal | STRTOG_Inexlo; 844 } 845 else if (bbbits == 1) { 846 irv = STRTOG_Normal; 847 drop_down: 848 /* boundary case -- decrement exponent */ 849 if (rve1 == emin) { 850 irv = STRTOG_Normal | STRTOG_Inexhi; 851 if (rvb->wds == 1 && rvb->x[0] == 1) 852 sudden_underflow = 1; 853 break; 854 } 855 rve -= nbits; 856 rvb = set_ones(rvb, rvbits = nbits); 857 if (rvb == NULL) 858 return STRTOG_NoMemory; 859 break; 860 } 861 else 862 irv = STRTOG_Normal | STRTOG_Inexhi; 863 if ((bbbits < nbits && !denorm) || !(rvb->x[0] & 1)) 864 break; 865 if (dsign) { 866 rvb = increment(rvb); 867 if (rvb == NULL) 868 return STRTOG_NoMemory; 869 if ( (j = rvbits & kmask) !=0) 870 j = ULbits - j; 871 if (hi0bits(rvb->x[(unsigned int)(rvb->wds - 1) 872 >> kshift]) 873 != j) 874 rvbits++; 875 irv = STRTOG_Normal | STRTOG_Inexhi; 876 } 877 else { 878 if (bbbits == 1) 879 goto undfl; 880 decrement(rvb); 881 irv = STRTOG_Normal | STRTOG_Inexlo; 882 } 883 break; 884 } 885 if ((dval(adj) = ratio(delta, bs)) <= 2.) { 886 adj1: 887 inex = STRTOG_Inexlo; 888 if (dsign) { 889 asub = 0; 890 inex = STRTOG_Inexhi; 891 } 892 else if (denorm && bbbits <= 1) { 893 undfl: 894 rvb->wds = 0; 895 rve = emin; 896 irv = STRTOG_Underflow | STRTOG_Inexlo; 897 break; 898 } 899 adj0 = dval(adj) = 1.; 900 } 901 else { 902 adj0 = dval(adj) *= 0.5; 903 if (dsign) { 904 asub = 0; 905 inex = STRTOG_Inexlo; 906 } 907 if (dval(adj) < 2147483647.) { 908 L = adj0; 909 adj0 -= L; 910 switch(rd) { 911 case 0: 912 if (adj0 >= .5) 913 goto inc_L; 914 break; 915 case 1: 916 if (asub && adj0 > 0.) 917 goto inc_L; 918 break; 919 case 2: 920 if (!asub && adj0 > 0.) { 921 inc_L: 922 L++; 923 inex = STRTOG_Inexact - inex; 924 } 925 } 926 dval(adj) = L; 927 } 928 } 929 y = rve + rvbits; 930 931 /* adj *= ulp(dval(rv)); */ 932 /* if (asub) rv -= adj; else rv += adj; */ 933 934 if (!denorm && rvbits < nbits) { 935 rvb = lshift(rvb, j = nbits - rvbits); 936 if (rvb == NULL) 937 return STRTOG_NoMemory; 938 rve -= j; 939 rvbits = nbits; 940 } 941 ab = d2b(dval(adj), &abe, &abits); 942 if (ab == NULL) 943 return STRTOG_NoMemory; 944 if (abe < 0) 945 rshift(ab, -abe); 946 else if (abe > 0) 947 ab = lshift(ab, abe); 948 rvb0 = rvb; 949 if (asub) { 950 /* rv -= adj; */ 951 j = hi0bits(rvb->x[rvb->wds-1]); 952 rvb = diff(rvb, ab); 953 if (rvb == NULL) 954 return STRTOG_NoMemory; 955 k = rvb0->wds - 1; 956 if (denorm) 957 /* do nothing */; 958 else if (rvb->wds <= k 959 || hi0bits( rvb->x[k]) > 960 hi0bits(rvb0->x[k])) { 961 /* unlikely; can only have lost 1 high bit */ 962 if (rve1 == emin) { 963 --rvbits; 964 denorm = 1; 965 } 966 else { 967 rvb = lshift(rvb, 1); 968 if (rvb == NULL) 969 return STRTOG_NoMemory; 970 --rve; 971 --rve1; 972 L = finished = 0; 973 } 974 } 975 } 976 else { 977 rvb = sum(rvb, ab); 978 if (rvb == NULL) 979 return STRTOG_NoMemory; 980 k = rvb->wds - 1; 981 if (k >= rvb0->wds 982 || hi0bits(rvb->x[k]) < hi0bits(rvb0->x[k])) { 983 if (denorm) { 984 if (++rvbits == nbits) 985 denorm = 0; 986 } 987 else { 988 rshift(rvb, 1); 989 rve++; 990 rve1++; 991 L = 0; 992 } 993 } 994 } 995 Bfree(ab); 996 Bfree(rvb0); 997 if (finished) 998 break; 999 1000 z = rve + rvbits; 1001 if (y == z && L) { 1002 /* Can we stop now? */ 1003 tol = dval(adj) * 5e-16; /* > max rel error */ 1004 dval(adj) = adj0 - .5; 1005 if (dval(adj) < -tol) { 1006 if (adj0 > tol) { 1007 irv |= inex; 1008 break; 1009 } 1010 } 1011 else if (dval(adj) > tol && adj0 < 1. - tol) { 1012 irv |= inex; 1013 break; 1014 } 1015 } 1016 bb0 = denorm ? 0 : trailz(rvb); 1017 Bfree(bb); 1018 Bfree(bd); 1019 Bfree(bs); 1020 Bfree(delta); 1021 } 1022 if (!denorm && (j = nbits - rvbits)) { 1023 if (j > 0) 1024 rvb = lshift(rvb, j); 1025 else 1026 rshift(rvb, -j); 1027 rve -= j; 1028 } 1029 *expt = rve; 1030 Bfree(bb); 1031 Bfree(bd); 1032 Bfree(bs); 1033 Bfree(bd0); 1034 Bfree(delta); 1035 if (rve > fpi->emax) { 1036 huge: 1037 rvb->wds = 0; 1038 irv = STRTOG_Infinite | STRTOG_Overflow | STRTOG_Inexhi; 1039 #ifndef NO_ERRNO 1040 errno = ERANGE; 1041 #endif 1042 #ifdef INFNAN_CHECK 1043 infnanexp: 1044 #endif 1045 *expt = fpi->emax + 1; 1046 } 1047 ret: 1048 if (denorm) { 1049 if (sudden_underflow) { 1050 rvb->wds = 0; 1051 irv = STRTOG_Underflow | STRTOG_Inexlo; 1052 } 1053 else { 1054 irv = (irv & ~STRTOG_Retmask) | 1055 (rvb->wds > 0 ? STRTOG_Denormal : STRTOG_Zero); 1056 if (irv & STRTOG_Inexact) 1057 irv |= STRTOG_Underflow; 1058 } 1059 } 1060 if (se) 1061 *se = __UNCONST(s); 1062 if (sign) 1063 irv |= STRTOG_Neg; 1064 if (rvb) { 1065 copybits(bits, nbits, rvb); 1066 Bfree(rvb); 1067 } 1068 return irv; 1069 } 1070