[BACK]Return to strtodg.c CVS log [TXT][DIR] Up to [cvs.NetBSD.org] / src / lib / libc / gdtoa

Annotation of src/lib/libc/gdtoa/strtodg.c, Revision 1.8

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

CVSweb <webmaster@jp.NetBSD.org>