1 /*
   2  * Copyright (c) 2003, 2011, Oracle and/or its affiliates. All rights reserved.
   3  * DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER.
   4  *
   5  * This code is free software; you can redistribute it and/or modify it
   6  * under the terms of the GNU General Public License version 2 only, as
   7  * published by the Free Software Foundation.  Oracle designates this
   8  * particular file as subject to the "Classpath" exception as provided
   9  * by Oracle in the LICENSE file that accompanied this code.
  10  *
  11  * This code is distributed in the hope that it will be useful, but WITHOUT
  12  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
  13  * FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  14  * version 2 for more details (a copy is included in the LICENSE file that
  15  * accompanied this code).
  16  *
  17  * You should have received a copy of the GNU General Public License version
  18  * 2 along with this work; if not, write to the Free Software Foundation,
  19  * Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
  20  *
  21  * Please contact Oracle, 500 Oracle Parkway, Redwood Shores, CA 94065 USA
  22  * or visit www.oracle.com if you need additional information or have any
  23  * questions.
  24  */
  25 
  26 package sun.misc;
  27 
  28 import sun.misc.FpUtils;
  29 import sun.misc.DoubleConsts;
  30 import sun.misc.FloatConsts;
  31 import java.util.regex.*;
  32 
  33 public class FormattedFloatingDecimal{
  34     boolean     isExceptional;
  35     boolean     isNegative;
  36     int         decExponent;  // value set at construction, then immutable
  37     int         decExponentRounded;
  38     char        digits[];
  39     int         nDigits;
  40     int         bigIntExp;
  41     int         bigIntNBits;
  42     boolean     mustSetRoundDir = false;
  43     boolean     fromHex = false;
  44     int         roundDir = 0; // set by doubleValue
  45     int         precision;    // number of digits to the right of decimal
  46 
  47     public enum Form { SCIENTIFIC, COMPATIBLE, DECIMAL_FLOAT, GENERAL };
  48 
  49     private Form form;
  50 
  51     private     FormattedFloatingDecimal( boolean negSign, int decExponent, char []digits, int n, boolean e, int precision, Form form )
  52     {
  53         isNegative = negSign;
  54         isExceptional = e;
  55         this.decExponent = decExponent;
  56         this.digits = digits;
  57         this.nDigits = n;
  58         this.precision = precision;
  59         this.form = form;
  60     }
  61 
  62     /*
  63      * Constants of the implementation
  64      * Most are IEEE-754 related.
  65      * (There are more really boring constants at the end.)
  66      */
  67     static final long   signMask = 0x8000000000000000L;
  68     static final long   expMask  = 0x7ff0000000000000L;
  69     static final long   fractMask= ~(signMask|expMask);
  70     static final int    expShift = 52;
  71     static final int    expBias  = 1023;
  72     static final long   fractHOB = ( 1L<<expShift ); // assumed High-Order bit
  73     static final long   expOne   = ((long)expBias)<<expShift; // exponent of 1.0
  74     static final int    maxSmallBinExp = 62;
  75     static final int    minSmallBinExp = -( 63 / 3 );
  76     static final int    maxDecimalDigits = 15;
  77     static final int    maxDecimalExponent = 308;
  78     static final int    minDecimalExponent = -324;
  79     static final int    bigDecimalExponent = 324; // i.e. abs(minDecimalExponent)
  80 
  81     static final long   highbyte = 0xff00000000000000L;
  82     static final long   highbit  = 0x8000000000000000L;
  83     static final long   lowbytes = ~highbyte;
  84 
  85     static final int    singleSignMask =    0x80000000;
  86     static final int    singleExpMask  =    0x7f800000;
  87     static final int    singleFractMask =   ~(singleSignMask|singleExpMask);
  88     static final int    singleExpShift  =   23;
  89     static final int    singleFractHOB  =   1<<singleExpShift;
  90     static final int    singleExpBias   =   127;
  91     static final int    singleMaxDecimalDigits = 7;
  92     static final int    singleMaxDecimalExponent = 38;
  93     static final int    singleMinDecimalExponent = -45;
  94 
  95     static final int    intDecimalDigits = 9;
  96 
  97 
  98     /*
  99      * count number of bits from high-order 1 bit to low-order 1 bit,
 100      * inclusive.
 101      */
 102     private static int
 103     countBits( long v ){
 104         //
 105         // the strategy is to shift until we get a non-zero sign bit
 106         // then shift until we have no bits left, counting the difference.
 107         // we do byte shifting as a hack. Hope it helps.
 108         //
 109         if ( v == 0L ) return 0;
 110 
 111         while ( ( v & highbyte ) == 0L ){
 112             v <<= 8;
 113         }
 114         while ( v > 0L ) { // i.e. while ((v&highbit) == 0L )
 115             v <<= 1;
 116         }
 117 
 118         int n = 0;
 119         while (( v & lowbytes ) != 0L ){
 120             v <<= 8;
 121             n += 8;
 122         }
 123         while ( v != 0L ){
 124             v <<= 1;
 125             n += 1;
 126         }
 127         return n;
 128     }
 129 
 130     /*
 131      * Keep big powers of 5 handy for future reference.
 132      */
 133     private static FDBigInt b5p[];
 134 
 135     private static synchronized FDBigInt
 136     big5pow( int p ){
 137         assert p >= 0 : p; // negative power of 5
 138         if ( b5p == null ){
 139             b5p = new FDBigInt[ p+1 ];
 140         }else if (b5p.length <= p ){
 141             FDBigInt t[] = new FDBigInt[ p+1 ];
 142             System.arraycopy( b5p, 0, t, 0, b5p.length );
 143             b5p = t;
 144         }
 145         if ( b5p[p] != null )
 146             return b5p[p];
 147         else if ( p < small5pow.length )
 148             return b5p[p] = new FDBigInt( small5pow[p] );
 149         else if ( p < long5pow.length )
 150             return b5p[p] = new FDBigInt( long5pow[p] );
 151         else {
 152             // construct the value.
 153             // recursively.
 154             int q, r;
 155             // in order to compute 5^p,
 156             // compute its square root, 5^(p/2) and square.
 157             // or, let q = p / 2, r = p -q, then
 158             // 5^p = 5^(q+r) = 5^q * 5^r
 159             q = p >> 1;
 160             r = p - q;
 161             FDBigInt bigq =  b5p[q];
 162             if ( bigq == null )
 163                 bigq = big5pow ( q );
 164             if ( r < small5pow.length ){
 165                 return (b5p[p] = bigq.mult( small5pow[r] ) );
 166             }else{
 167                 FDBigInt bigr = b5p[ r ];
 168                 if ( bigr == null )
 169                     bigr = big5pow( r );
 170                 return (b5p[p] = bigq.mult( bigr ) );
 171             }
 172         }
 173     }
 174 
 175     //
 176     // a common operation
 177     //
 178     private static FDBigInt
 179     multPow52( FDBigInt v, int p5, int p2 ){
 180         if ( p5 != 0 ){
 181             if ( p5 < small5pow.length ){
 182                 v = v.mult( small5pow[p5] );
 183             } else {
 184                 v = v.mult( big5pow( p5 ) );
 185             }
 186         }
 187         if ( p2 != 0 ){
 188             v.lshiftMe( p2 );
 189         }
 190         return v;
 191     }
 192 
 193     //
 194     // another common operation
 195     //
 196     private static FDBigInt
 197     constructPow52( int p5, int p2 ){
 198         FDBigInt v = new FDBigInt( big5pow( p5 ) );
 199         if ( p2 != 0 ){
 200             v.lshiftMe( p2 );
 201         }
 202         return v;
 203     }
 204 
 205     /*
 206      * Make a floating double into a FDBigInt.
 207      * This could also be structured as a FDBigInt
 208      * constructor, but we'd have to build a lot of knowledge
 209      * about floating-point representation into it, and we don't want to.
 210      *
 211      * AS A SIDE EFFECT, THIS METHOD WILL SET THE INSTANCE VARIABLES
 212      * bigIntExp and bigIntNBits
 213      *
 214      */
 215     private FDBigInt
 216     doubleToBigInt( double dval ){
 217         long lbits = Double.doubleToLongBits( dval ) & ~signMask;
 218         int binexp = (int)(lbits >>> expShift);
 219         lbits &= fractMask;
 220         if ( binexp > 0 ){
 221             lbits |= fractHOB;
 222         } else {
 223             assert lbits != 0L : lbits; // doubleToBigInt(0.0)
 224             binexp +=1;
 225             while ( (lbits & fractHOB ) == 0L){
 226                 lbits <<= 1;
 227                 binexp -= 1;
 228             }
 229         }
 230         binexp -= expBias;
 231         int nbits = countBits( lbits );
 232         /*
 233          * We now know where the high-order 1 bit is,
 234          * and we know how many there are.
 235          */
 236         int lowOrderZeros = expShift+1-nbits;
 237         lbits >>>= lowOrderZeros;
 238 
 239         bigIntExp = binexp+1-nbits;
 240         bigIntNBits = nbits;
 241         return new FDBigInt( lbits );
 242     }
 243 
 244     /*
 245      * Compute a number that is the ULP of the given value,
 246      * for purposes of addition/subtraction. Generally easy.
 247      * More difficult if subtracting and the argument
 248      * is a normalized a power of 2, as the ULP changes at these points.
 249      */
 250     private static double ulp( double dval, boolean subtracting ){
 251         long lbits = Double.doubleToLongBits( dval ) & ~signMask;
 252         int binexp = (int)(lbits >>> expShift);
 253         double ulpval;
 254         if ( subtracting && ( binexp >= expShift ) && ((lbits&fractMask) == 0L) ){
 255             // for subtraction from normalized, powers of 2,
 256             // use next-smaller exponent
 257             binexp -= 1;
 258         }
 259         if ( binexp > expShift ){
 260             ulpval = Double.longBitsToDouble( ((long)(binexp-expShift))<<expShift );
 261         } else if ( binexp == 0 ){
 262             ulpval = Double.MIN_VALUE;
 263         } else {
 264             ulpval = Double.longBitsToDouble( 1L<<(binexp-1) );
 265         }
 266         if ( subtracting ) ulpval = - ulpval;
 267 
 268         return ulpval;
 269     }
 270 
 271     /*
 272      * Round a double to a float.
 273      * In addition to the fraction bits of the double,
 274      * look at the class instance variable roundDir,
 275      * which should help us avoid double-rounding error.
 276      * roundDir was set in hardValueOf if the estimate was
 277      * close enough, but not exact. It tells us which direction
 278      * of rounding is preferred.
 279      */
 280     float
 281     stickyRound( double dval ){
 282         long lbits = Double.doubleToLongBits( dval );
 283         long binexp = lbits & expMask;
 284         if ( binexp == 0L || binexp == expMask ){
 285             // what we have here is special.
 286             // don't worry, the right thing will happen.
 287             return (float) dval;
 288         }
 289         lbits += (long)roundDir; // hack-o-matic.
 290         return (float)Double.longBitsToDouble( lbits );
 291     }
 292 
 293 
 294     /*
 295      * This is the easy subcase --
 296      * all the significant bits, after scaling, are held in lvalue.
 297      * negSign and decExponent tell us what processing and scaling
 298      * has already been done. Exceptional cases have already been
 299      * stripped out.
 300      * In particular:
 301      * lvalue is a finite number (not Inf, nor NaN)
 302      * lvalue > 0L (not zero, nor negative).
 303      *
 304      * The only reason that we develop the digits here, rather than
 305      * calling on Long.toString() is that we can do it a little faster,
 306      * and besides want to treat trailing 0s specially. If Long.toString
 307      * changes, we should re-evaluate this strategy!
 308      */
 309     private void
 310     developLongDigits( int decExponent, long lvalue, long insignificant ){
 311         char digits[];
 312         int  ndigits;
 313         int  digitno;
 314         int  c;
 315         //
 316         // Discard non-significant low-order bits, while rounding,
 317         // up to insignificant value.
 318         int i;
 319         for ( i = 0; insignificant >= 10L; i++ )
 320             insignificant /= 10L;
 321         if ( i != 0 ){
 322             long pow10 = long5pow[i] << i; // 10^i == 5^i * 2^i;
 323             long residue = lvalue % pow10;
 324             lvalue /= pow10;
 325             decExponent += i;
 326             if ( residue >= (pow10>>1) ){
 327                 // round up based on the low-order bits we're discarding
 328                 lvalue++;
 329             }
 330         }
 331         if ( lvalue <= Integer.MAX_VALUE ){
 332             assert lvalue > 0L : lvalue; // lvalue <= 0
 333             // even easier subcase!
 334             // can do int arithmetic rather than long!
 335             int  ivalue = (int)lvalue;
 336             ndigits = 10;
 337             digits = (char[])(perThreadBuffer.get());
 338             digitno = ndigits-1;
 339             c = ivalue%10;
 340             ivalue /= 10;
 341             while ( c == 0 ){
 342                 decExponent++;
 343                 c = ivalue%10;
 344                 ivalue /= 10;
 345             }
 346             while ( ivalue != 0){
 347                 digits[digitno--] = (char)(c+'0');
 348                 decExponent++;
 349                 c = ivalue%10;
 350                 ivalue /= 10;
 351             }
 352             digits[digitno] = (char)(c+'0');
 353         } else {
 354             // same algorithm as above (same bugs, too )
 355             // but using long arithmetic.
 356             ndigits = 20;
 357             digits = (char[])(perThreadBuffer.get());
 358             digitno = ndigits-1;
 359             c = (int)(lvalue%10L);
 360             lvalue /= 10L;
 361             while ( c == 0 ){
 362                 decExponent++;
 363                 c = (int)(lvalue%10L);
 364                 lvalue /= 10L;
 365             }
 366             while ( lvalue != 0L ){
 367                 digits[digitno--] = (char)(c+'0');
 368                 decExponent++;
 369                 c = (int)(lvalue%10L);
 370                 lvalue /= 10;
 371             }
 372             digits[digitno] = (char)(c+'0');
 373         }
 374         char result [];
 375         ndigits -= digitno;
 376         result = new char[ ndigits ];
 377         System.arraycopy( digits, digitno, result, 0, ndigits );
 378         this.digits = result;
 379         this.decExponent = decExponent+1;
 380         this.nDigits = ndigits;
 381     }
 382 
 383     //
 384     // add one to the least significant digit.
 385     // in the unlikely event there is a carry out,
 386     // deal with it.
 387     // assert that this will only happen where there
 388     // is only one digit, e.g. (float)1e-44 seems to do it.
 389     //
 390     private void
 391     roundup(){
 392         int i;
 393         int q = digits[ i = (nDigits-1)];
 394         if ( q == '9' ){
 395             while ( q == '9' && i > 0 ){
 396                 digits[i] = '0';
 397                 q = digits[--i];
 398             }
 399             if ( q == '9' ){
 400                 // carryout! High-order 1, rest 0s, larger exp.
 401                 decExponent += 1;
 402                 digits[0] = '1';
 403                 return;
 404             }
 405             // else fall through.
 406         }
 407         digits[i] = (char)(q+1);
 408     }
 409 
 410     // Given the desired number of digits predict the result's exponent.
 411     private int checkExponent(int length) {
 412         if (length >= nDigits || length < 0)
 413             return decExponent;
 414 
 415         for (int i = 0; i < length; i++)
 416             if (digits[i] != '9')
 417                 // a '9' anywhere in digits will absorb the round
 418                 return decExponent;
 419         return decExponent + (digits[length] >= '5' ? 1 : 0);
 420     }
 421 
 422     // Unlike roundup(), this method does not modify digits.  It also
 423     // rounds at a particular precision.
 424     private char [] applyPrecision(int length) {
 425         char [] result = new char[nDigits];
 426         for (int i = 0; i < result.length; i++) result[i] = '0';
 427 
 428         if (length >= nDigits || length < 0) {
 429             // no rounding necessary
 430             System.arraycopy(digits, 0, result, 0, nDigits);
 431             return result;
 432         }
 433         if (length == 0) {
 434             // only one digit (0 or 1) is returned because the precision
 435             // excludes all significant digits
 436             if (digits[0] >= '5') {
 437                 result[0] = '1';
 438             }
 439             return result;
 440         }
 441 
 442         int i = length;
 443         int q = digits[i];
 444         if (q >= '5' && i > 0) {
 445             q = digits[--i];
 446             if ( q == '9' ) {
 447                 while ( q == '9' && i > 0 ){
 448                     q = digits[--i];
 449                 }
 450                 if ( q == '9' ){
 451                     // carryout! High-order 1, rest 0s, larger exp.
 452                     result[0] = '1';
 453                     return result;
 454                 }
 455             }
 456             result[i] = (char)(q + 1);
 457         }
 458         while (--i >= 0) {
 459             result[i] = digits[i];
 460         }
 461         return result;
 462     }
 463 
 464     /*
 465      * FIRST IMPORTANT CONSTRUCTOR: DOUBLE
 466      */
 467     public FormattedFloatingDecimal( double d )
 468     {
 469         this(d, Integer.MAX_VALUE, Form.COMPATIBLE);
 470     }
 471 
 472     public FormattedFloatingDecimal( double d, int precision, Form form )
 473     {
 474         long    dBits = Double.doubleToLongBits( d );
 475         long    fractBits;
 476         int     binExp;
 477         int     nSignificantBits;
 478 
 479         this.precision = precision;
 480         this.form      = form;
 481 
 482         // discover and delete sign
 483         if ( (dBits&signMask) != 0 ){
 484             isNegative = true;
 485             dBits ^= signMask;
 486         } else {
 487             isNegative = false;
 488         }
 489         // Begin to unpack
 490         // Discover obvious special cases of NaN and Infinity.
 491         binExp = (int)( (dBits&expMask) >> expShift );
 492         fractBits = dBits&fractMask;
 493         if ( binExp == (int)(expMask>>expShift) ) {
 494             isExceptional = true;
 495             if ( fractBits == 0L ){
 496                 digits =  infinity;
 497             } else {
 498                 digits = notANumber;
 499                 isNegative = false; // NaN has no sign!
 500             }
 501             nDigits = digits.length;
 502             return;
 503         }
 504         isExceptional = false;
 505         // Finish unpacking
 506         // Normalize denormalized numbers.
 507         // Insert assumed high-order bit for normalized numbers.
 508         // Subtract exponent bias.
 509         if ( binExp == 0 ){
 510             if ( fractBits == 0L ){
 511                 // not a denorm, just a 0!
 512                 decExponent = 0;
 513                 digits = zero;
 514                 nDigits = 1;
 515                 return;
 516             }
 517             while ( (fractBits&fractHOB) == 0L ){
 518                 fractBits <<= 1;
 519                 binExp -= 1;
 520             }
 521             nSignificantBits = expShift + binExp +1; // recall binExp is  - shift count.
 522             binExp += 1;
 523         } else {
 524             fractBits |= fractHOB;
 525             nSignificantBits = expShift+1;
 526         }
 527         binExp -= expBias;
 528         // call the routine that actually does all the hard work.
 529         dtoa( binExp, fractBits, nSignificantBits );
 530     }
 531 
 532     /*
 533      * SECOND IMPORTANT CONSTRUCTOR: SINGLE
 534      */
 535     public FormattedFloatingDecimal( float f )
 536     {
 537         this(f, Integer.MAX_VALUE, Form.COMPATIBLE);
 538     }
 539     public FormattedFloatingDecimal( float f, int precision, Form form )
 540     {
 541         int     fBits = Float.floatToIntBits( f );
 542         int     fractBits;
 543         int     binExp;
 544         int     nSignificantBits;
 545 
 546         this.precision = precision;
 547         this.form      = form;
 548 
 549         // discover and delete sign
 550         if ( (fBits&singleSignMask) != 0 ){
 551             isNegative = true;
 552             fBits ^= singleSignMask;
 553         } else {
 554             isNegative = false;
 555         }
 556         // Begin to unpack
 557         // Discover obvious special cases of NaN and Infinity.
 558         binExp = (int)( (fBits&singleExpMask) >> singleExpShift );
 559         fractBits = fBits&singleFractMask;
 560         if ( binExp == (int)(singleExpMask>>singleExpShift) ) {
 561             isExceptional = true;
 562             if ( fractBits == 0L ){
 563                 digits =  infinity;
 564             } else {
 565                 digits = notANumber;
 566                 isNegative = false; // NaN has no sign!
 567             }
 568             nDigits = digits.length;
 569             return;
 570         }
 571         isExceptional = false;
 572         // Finish unpacking
 573         // Normalize denormalized numbers.
 574         // Insert assumed high-order bit for normalized numbers.
 575         // Subtract exponent bias.
 576         if ( binExp == 0 ){
 577             if ( fractBits == 0 ){
 578                 // not a denorm, just a 0!
 579                 decExponent = 0;
 580                 digits = zero;
 581                 nDigits = 1;
 582                 return;
 583             }
 584             while ( (fractBits&singleFractHOB) == 0 ){
 585                 fractBits <<= 1;
 586                 binExp -= 1;
 587             }
 588             nSignificantBits = singleExpShift + binExp +1; // recall binExp is  - shift count.
 589             binExp += 1;
 590         } else {
 591             fractBits |= singleFractHOB;
 592             nSignificantBits = singleExpShift+1;
 593         }
 594         binExp -= singleExpBias;
 595         // call the routine that actually does all the hard work.
 596         dtoa( binExp, ((long)fractBits)<<(expShift-singleExpShift), nSignificantBits );
 597     }
 598 
 599     private void
 600     dtoa( int binExp, long fractBits, int nSignificantBits )
 601     {
 602         int     nFractBits; // number of significant bits of fractBits;
 603         int     nTinyBits;  // number of these to the right of the point.
 604         int     decExp;
 605 
 606         // Examine number. Determine if it is an easy case,
 607         // which we can do pretty trivially using float/long conversion,
 608         // or whether we must do real work.
 609         nFractBits = countBits( fractBits );
 610         nTinyBits = Math.max( 0, nFractBits - binExp - 1 );
 611         if ( binExp <= maxSmallBinExp && binExp >= minSmallBinExp ){
 612             // Look more closely at the number to decide if,
 613             // with scaling by 10^nTinyBits, the result will fit in
 614             // a long.
 615             if ( (nTinyBits < long5pow.length) && ((nFractBits + n5bits[nTinyBits]) < 64 ) ){
 616                 /*
 617                  * We can do this:
 618                  * take the fraction bits, which are normalized.
 619                  * (a) nTinyBits == 0: Shift left or right appropriately
 620                  *     to align the binary point at the extreme right, i.e.
 621                  *     where a long int point is expected to be. The integer
 622                  *     result is easily converted to a string.
 623                  * (b) nTinyBits > 0: Shift right by expShift-nFractBits,
 624                  *     which effectively converts to long and scales by
 625                  *     2^nTinyBits. Then multiply by 5^nTinyBits to
 626                  *     complete the scaling. We know this won't overflow
 627                  *     because we just counted the number of bits necessary
 628                  *     in the result. The integer you get from this can
 629                  *     then be converted to a string pretty easily.
 630                  */
 631                 long halfULP;
 632                 if ( nTinyBits == 0 ) {
 633                     if ( binExp > nSignificantBits ){
 634                         halfULP = 1L << ( binExp-nSignificantBits-1);
 635                     } else {
 636                         halfULP = 0L;
 637                     }
 638                     if ( binExp >= expShift ){
 639                         fractBits <<= (binExp-expShift);
 640                     } else {
 641                         fractBits >>>= (expShift-binExp) ;
 642                     }
 643                     developLongDigits( 0, fractBits, halfULP );
 644                     return;
 645                 }
 646                 /*
 647                  * The following causes excess digits to be printed
 648                  * out in the single-float case. Our manipulation of
 649                  * halfULP here is apparently not correct. If we
 650                  * better understand how this works, perhaps we can
 651                  * use this special case again. But for the time being,
 652                  * we do not.
 653                  * else {
 654                  *     fractBits >>>= expShift+1-nFractBits;
 655                  *     fractBits *= long5pow[ nTinyBits ];
 656                  *     halfULP = long5pow[ nTinyBits ] >> (1+nSignificantBits-nFractBits);
 657                  *     developLongDigits( -nTinyBits, fractBits, halfULP );
 658                  *     return;
 659                  * }
 660                  */
 661             }
 662         }
 663         /*
 664          * This is the hard case. We are going to compute large positive
 665          * integers B and S and integer decExp, s.t.
 666          *      d = ( B / S ) * 10^decExp
 667          *      1 <= B / S < 10
 668          * Obvious choices are:
 669          *      decExp = floor( log10(d) )
 670          *      B      = d * 2^nTinyBits * 10^max( 0, -decExp )
 671          *      S      = 10^max( 0, decExp) * 2^nTinyBits
 672          * (noting that nTinyBits has already been forced to non-negative)
 673          * I am also going to compute a large positive integer
 674          *      M      = (1/2^nSignificantBits) * 2^nTinyBits * 10^max( 0, -decExp )
 675          * i.e. M is (1/2) of the ULP of d, scaled like B.
 676          * When we iterate through dividing B/S and picking off the
 677          * quotient bits, we will know when to stop when the remainder
 678          * is <= M.
 679          *
 680          * We keep track of powers of 2 and powers of 5.
 681          */
 682 
 683         /*
 684          * Estimate decimal exponent. (If it is small-ish,
 685          * we could double-check.)
 686          *
 687          * First, scale the mantissa bits such that 1 <= d2 < 2.
 688          * We are then going to estimate
 689          *          log10(d2) ~=~  (d2-1.5)/1.5 + log(1.5)
 690          * and so we can estimate
 691          *      log10(d) ~=~ log10(d2) + binExp * log10(2)
 692          * take the floor and call it decExp.
 693          * FIXME -- use more precise constants here. It costs no more.
 694          */
 695         double d2 = Double.longBitsToDouble(
 696             expOne | ( fractBits &~ fractHOB ) );
 697         decExp = (int)Math.floor(
 698             (d2-1.5D)*0.289529654D + 0.176091259 + (double)binExp * 0.301029995663981 );
 699         int B2, B5; // powers of 2 and powers of 5, respectively, in B
 700         int S2, S5; // powers of 2 and powers of 5, respectively, in S
 701         int M2, M5; // powers of 2 and powers of 5, respectively, in M
 702         int Bbits; // binary digits needed to represent B, approx.
 703         int tenSbits; // binary digits needed to represent 10*S, approx.
 704         FDBigInt Sval, Bval, Mval;
 705 
 706         B5 = Math.max( 0, -decExp );
 707         B2 = B5 + nTinyBits + binExp;
 708 
 709         S5 = Math.max( 0, decExp );
 710         S2 = S5 + nTinyBits;
 711 
 712         M5 = B5;
 713         M2 = B2 - nSignificantBits;
 714 
 715         /*
 716          * the long integer fractBits contains the (nFractBits) interesting
 717          * bits from the mantissa of d ( hidden 1 added if necessary) followed
 718          * by (expShift+1-nFractBits) zeros. In the interest of compactness,
 719          * I will shift out those zeros before turning fractBits into a
 720          * FDBigInt. The resulting whole number will be
 721          *      d * 2^(nFractBits-1-binExp).
 722          */
 723         fractBits >>>= (expShift+1-nFractBits);
 724         B2 -= nFractBits-1;
 725         int common2factor = Math.min( B2, S2 );
 726         B2 -= common2factor;
 727         S2 -= common2factor;
 728         M2 -= common2factor;
 729 
 730         /*
 731          * HACK!! For exact powers of two, the next smallest number
 732          * is only half as far away as we think (because the meaning of
 733          * ULP changes at power-of-two bounds) for this reason, we
 734          * hack M2. Hope this works.
 735          */
 736         if ( nFractBits == 1 )
 737             M2 -= 1;
 738 
 739         if ( M2 < 0 ){
 740             // oops.
 741             // since we cannot scale M down far enough,
 742             // we must scale the other values up.
 743             B2 -= M2;
 744             S2 -= M2;
 745             M2 =  0;
 746         }
 747         /*
 748          * Construct, Scale, iterate.
 749          * Some day, we'll write a stopping test that takes
 750          * account of the assymetry of the spacing of floating-point
 751          * numbers below perfect powers of 2
 752          * 26 Sept 96 is not that day.
 753          * So we use a symmetric test.
 754          */
 755         char digits[] = this.digits = new char[18];
 756         int  ndigit = 0;
 757         boolean low, high;
 758         long lowDigitDifference;
 759         int  q;
 760 
 761         /*
 762          * Detect the special cases where all the numbers we are about
 763          * to compute will fit in int or long integers.
 764          * In these cases, we will avoid doing FDBigInt arithmetic.
 765          * We use the same algorithms, except that we "normalize"
 766          * our FDBigInts before iterating. This is to make division easier,
 767          * as it makes our fist guess (quotient of high-order words)
 768          * more accurate!
 769          *
 770          * Some day, we'll write a stopping test that takes
 771          * account of the assymetry of the spacing of floating-point
 772          * numbers below perfect powers of 2
 773          * 26 Sept 96 is not that day.
 774          * So we use a symmetric test.
 775          */
 776         Bbits = nFractBits + B2 + (( B5 < n5bits.length )? n5bits[B5] : ( B5*3 ));
 777         tenSbits = S2+1 + (( (S5+1) < n5bits.length )? n5bits[(S5+1)] : ( (S5+1)*3 ));
 778         if ( Bbits < 64 && tenSbits < 64){
 779             if ( Bbits < 32 && tenSbits < 32){
 780                 // wa-hoo! They're all ints!
 781                 int b = ((int)fractBits * small5pow[B5] ) << B2;
 782                 int s = small5pow[S5] << S2;
 783                 int m = small5pow[M5] << M2;
 784                 int tens = s * 10;
 785                 /*
 786                  * Unroll the first iteration. If our decExp estimate
 787                  * was too high, our first quotient will be zero. In this
 788                  * case, we discard it and decrement decExp.
 789                  */
 790                 ndigit = 0;
 791                 q = b / s;
 792                 b = 10 * ( b % s );
 793                 m *= 10;
 794                 low  = (b <  m );
 795                 high = (b+m > tens );
 796                 assert q < 10 : q; // excessively large digit
 797                 if ( (q == 0) && ! high ){
 798                     // oops. Usually ignore leading zero.
 799                     decExp--;
 800                 } else {
 801                     digits[ndigit++] = (char)('0' + q);
 802                 }
 803                 /*
 804                  * HACK! Java spec sez that we always have at least
 805                  * one digit after the . in either F- or E-form output.
 806                  * Thus we will need more than one digit if we're using
 807                  * E-form
 808                  */
 809                 if (! (form == Form.COMPATIBLE && -3 < decExp && decExp < 8)) {
 810                     high = low = false;
 811                 }
 812                 while( ! low && ! high ){
 813                     q = b / s;
 814                     b = 10 * ( b % s );
 815                     m *= 10;
 816                     assert q < 10 : q; // excessively large digit
 817                     if ( m > 0L ){
 818                         low  = (b <  m );
 819                         high = (b+m > tens );
 820                     } else {
 821                         // hack -- m might overflow!
 822                         // in this case, it is certainly > b,
 823                         // which won't
 824                         // and b+m > tens, too, since that has overflowed
 825                         // either!
 826                         low = true;
 827                         high = true;
 828                     }
 829                     digits[ndigit++] = (char)('0' + q);
 830                 }
 831                 lowDigitDifference = (b<<1) - tens;
 832             } else {
 833                 // still good! they're all longs!
 834                 long b = (fractBits * long5pow[B5] ) << B2;
 835                 long s = long5pow[S5] << S2;
 836                 long m = long5pow[M5] << M2;
 837                 long tens = s * 10L;
 838                 /*
 839                  * Unroll the first iteration. If our decExp estimate
 840                  * was too high, our first quotient will be zero. In this
 841                  * case, we discard it and decrement decExp.
 842                  */
 843                 ndigit = 0;
 844                 q = (int) ( b / s );
 845                 b = 10L * ( b % s );
 846                 m *= 10L;
 847                 low  = (b <  m );
 848                 high = (b+m > tens );
 849                 assert q < 10 : q; // excessively large digit
 850                 if ( (q == 0) && ! high ){
 851                     // oops. Usually ignore leading zero.
 852                     decExp--;
 853                 } else {
 854                     digits[ndigit++] = (char)('0' + q);
 855                 }
 856                 /*
 857                  * HACK! Java spec sez that we always have at least
 858                  * one digit after the . in either F- or E-form output.
 859                  * Thus we will need more than one digit if we're using
 860                  * E-form
 861                  */
 862                 if (! (form == Form.COMPATIBLE && -3 < decExp && decExp < 8)) {
 863                     high = low = false;
 864                 }
 865                 while( ! low && ! high ){
 866                     q = (int) ( b / s );
 867                     b = 10 * ( b % s );
 868                     m *= 10;
 869                     assert q < 10 : q;  // excessively large digit
 870                     if ( m > 0L ){
 871                         low  = (b <  m );
 872                         high = (b+m > tens );
 873                     } else {
 874                         // hack -- m might overflow!
 875                         // in this case, it is certainly > b,
 876                         // which won't
 877                         // and b+m > tens, too, since that has overflowed
 878                         // either!
 879                         low = true;
 880                         high = true;
 881                     }
 882                     digits[ndigit++] = (char)('0' + q);
 883                 }
 884                 lowDigitDifference = (b<<1) - tens;
 885             }
 886         } else {
 887             FDBigInt tenSval;
 888             int  shiftBias;
 889 
 890             /*
 891              * We really must do FDBigInt arithmetic.
 892              * Fist, construct our FDBigInt initial values.
 893              */
 894             Bval = multPow52( new FDBigInt( fractBits  ), B5, B2 );
 895             Sval = constructPow52( S5, S2 );
 896             Mval = constructPow52( M5, M2 );
 897 
 898 
 899             // normalize so that division works better
 900             Bval.lshiftMe( shiftBias = Sval.normalizeMe() );
 901             Mval.lshiftMe( shiftBias );
 902             tenSval = Sval.mult( 10 );
 903             /*
 904              * Unroll the first iteration. If our decExp estimate
 905              * was too high, our first quotient will be zero. In this
 906              * case, we discard it and decrement decExp.
 907              */
 908             ndigit = 0;
 909             q = Bval.quoRemIteration( Sval );
 910             Mval = Mval.mult( 10 );
 911             low  = (Bval.cmp( Mval ) < 0);
 912             high = (Bval.add( Mval ).cmp( tenSval ) > 0 );
 913             assert q < 10 : q; // excessively large digit
 914             if ( (q == 0) && ! high ){
 915                 // oops. Usually ignore leading zero.
 916                 decExp--;
 917             } else {
 918                 digits[ndigit++] = (char)('0' + q);
 919             }
 920             /*
 921              * HACK! Java spec sez that we always have at least
 922              * one digit after the . in either F- or E-form output.
 923              * Thus we will need more than one digit if we're using
 924              * E-form
 925              */
 926             if (! (form == Form.COMPATIBLE && -3 < decExp && decExp < 8)) {
 927                 high = low = false;
 928             }
 929             while( ! low && ! high ){
 930                 q = Bval.quoRemIteration( Sval );
 931                 Mval = Mval.mult( 10 );
 932                 assert q < 10 : q;  // excessively large digit
 933                 low  = (Bval.cmp( Mval ) < 0);
 934                 high = (Bval.add( Mval ).cmp( tenSval ) > 0 );
 935                 digits[ndigit++] = (char)('0' + q);
 936             }
 937             if ( high && low ){
 938                 Bval.lshiftMe(1);
 939                 lowDigitDifference = Bval.cmp(tenSval);
 940             } else
 941                 lowDigitDifference = 0L; // this here only for flow analysis!
 942         }
 943         this.decExponent = decExp+1;
 944         this.digits = digits;
 945         this.nDigits = ndigit;
 946         /*
 947          * Last digit gets rounded based on stopping condition.
 948          */
 949         if ( high ){
 950             if ( low ){
 951                 if ( lowDigitDifference == 0L ){
 952                     // it's a tie!
 953                     // choose based on which digits we like.
 954                     if ( (digits[nDigits-1]&1) != 0 ) roundup();
 955                 } else if ( lowDigitDifference > 0 ){
 956                     roundup();
 957                 }
 958             } else {
 959                 roundup();
 960             }
 961         }
 962     }
 963 
 964     public String
 965     toString(){
 966         // most brain-dead version
 967         StringBuffer result = new StringBuffer( nDigits+8 );
 968         if ( isNegative ){ result.append( '-' ); }
 969         if ( isExceptional ){
 970             result.append( digits, 0, nDigits );
 971         } else {
 972             result.append( "0.");
 973             result.append( digits, 0, nDigits );
 974             result.append('e');
 975             result.append( decExponent );
 976         }
 977         return new String(result);
 978     }
 979 
 980     // returns the exponent before rounding
 981     public int getExponent() {
 982         return decExponent - 1;
 983     }
 984 
 985     // returns the exponent after rounding has been done by applyPrecision
 986     public int getExponentRounded() {
 987         return decExponentRounded - 1;
 988     }
 989 
 990     public int getChars(char[] result) {
 991         assert nDigits <= 19 : nDigits; // generous bound on size of nDigits
 992         int i = 0;
 993         if (isNegative) { result[0] = '-'; i = 1; }
 994         if (isExceptional) {
 995             System.arraycopy(digits, 0, result, i, nDigits);
 996             i += nDigits;
 997         } else {
 998             char digits [] = this.digits;
 999             int exp = decExponent;
1000             switch (form) {
1001             case COMPATIBLE:
1002                 break;
1003             case DECIMAL_FLOAT:
1004                 exp = checkExponent(decExponent + precision);
1005                 digits = applyPrecision(decExponent + precision);
1006                 break;
1007             case SCIENTIFIC:
1008                 exp = checkExponent(precision + 1);
1009                 digits = applyPrecision(precision + 1);
1010                 break;
1011             case GENERAL:
1012                 exp = checkExponent(precision);
1013                 digits = applyPrecision(precision);
1014                 // adjust precision to be the number of digits to right of decimal
1015                 // the real exponent to be output is actually exp - 1, not exp
1016                 if (exp - 1 < -4 || exp - 1 >= precision) {
1017                     form = Form.SCIENTIFIC;
1018                     precision--;
1019                 } else {
1020                     form = Form.DECIMAL_FLOAT;
1021                     precision = precision - exp;
1022                 }
1023                 break;
1024             default:
1025                 assert false;
1026             }
1027             decExponentRounded = exp;
1028 
1029             if (exp > 0
1030                 && ((form == Form.COMPATIBLE && (exp < 8))
1031                     || (form == Form.DECIMAL_FLOAT)))
1032             {
1033                 // print digits.digits.
1034                 int charLength = Math.min(nDigits, exp);
1035                 System.arraycopy(digits, 0, result, i, charLength);
1036                 i += charLength;
1037                 if (charLength < exp) {
1038                     charLength = exp-charLength;
1039                     for (int nz = 0; nz < charLength; nz++)
1040                         result[i++] = '0';
1041                     // Do not append ".0" for formatted floats since the user
1042                     // may request that it be omitted. It is added as necessary
1043                     // by the Formatter.
1044                     if (form == Form.COMPATIBLE) {
1045                         result[i++] = '.';
1046                         result[i++] = '0';
1047                     }
1048                 } else {
1049                     // Do not append ".0" for formatted floats since the user
1050                     // may request that it be omitted. It is added as necessary
1051                     // by the Formatter.
1052                     if (form == Form.COMPATIBLE) {
1053                         result[i++] = '.';
1054                         if (charLength < nDigits) {
1055                             int t = Math.min(nDigits - charLength, precision);
1056                             System.arraycopy(digits, charLength, result, i, t);
1057                             i += t;
1058                         } else {
1059                             result[i++] = '0';
1060                         }
1061                     } else {
1062                         int t = Math.min(nDigits - charLength, precision);
1063                         if (t > 0) {
1064                             result[i++] = '.';
1065                             System.arraycopy(digits, charLength, result, i, t);
1066                             i += t;
1067                         }
1068                     }
1069                 }
1070             } else if (exp <= 0
1071                        && ((form == Form.COMPATIBLE && exp > -3)
1072                        || (form == Form.DECIMAL_FLOAT)))
1073             {
1074                 // print 0.0* digits
1075                 result[i++] = '0';
1076                 if (exp != 0) {
1077                     // write '0' s before the significant digits
1078                     int t = Math.min(-exp, precision);
1079                     if (t > 0) {
1080                         result[i++] = '.';
1081                         for (int nz = 0; nz < t; nz++)
1082                             result[i++] = '0';
1083                     }
1084                 }
1085                 int t = Math.min(digits.length, precision + exp);
1086                 if (t > 0) {
1087                     if (i == 1)
1088                         result[i++] = '.';
1089                     // copy only when significant digits are within the precision
1090                     System.arraycopy(digits, 0, result, i, t);
1091                     i += t;
1092                 }
1093             } else {
1094                 result[i++] = digits[0];
1095                 if (form == Form.COMPATIBLE) {
1096                     result[i++] = '.';
1097                     if (nDigits > 1) {
1098                         System.arraycopy(digits, 1, result, i, nDigits-1);
1099                         i += nDigits-1;
1100                     } else {
1101                         result[i++] = '0';
1102                     }
1103                     result[i++] = 'E';
1104                 } else {
1105                     if (nDigits > 1) {
1106                         int t = Math.min(nDigits -1, precision);
1107                         if (t > 0) {
1108                             result[i++] = '.';
1109                             System.arraycopy(digits, 1, result, i, t);
1110                             i += t;
1111                         }
1112                     }
1113                     result[i++] = 'e';
1114                 }
1115                 int e;
1116                 if (exp <= 0) {
1117                     result[i++] = '-';
1118                     e = -exp+1;
1119                 } else {
1120                     if (form != Form.COMPATIBLE)
1121                         result[i++] = '+';
1122                     e = exp-1;
1123                 }
1124                 // decExponent has 1, 2, or 3, digits
1125                 if (e <= 9) {
1126                     if (form != Form.COMPATIBLE)
1127                         result[i++] = '0';
1128                     result[i++] = (char)(e+'0');
1129                 } else if (e <= 99) {
1130                     result[i++] = (char)(e/10 +'0');
1131                     result[i++] = (char)(e%10 + '0');
1132                 } else {
1133                     result[i++] = (char)(e/100+'0');
1134                     e %= 100;
1135                     result[i++] = (char)(e/10+'0');
1136                     result[i++] = (char)(e%10 + '0');
1137                 }
1138             }
1139         }
1140         return i;
1141     }
1142 
1143     // Per-thread buffer for string/stringbuffer conversion
1144     private static ThreadLocal perThreadBuffer = new ThreadLocal() {
1145             protected synchronized Object initialValue() {
1146                 return new char[26];
1147             }
1148         };
1149 
1150     /*
1151      * Take a FormattedFloatingDecimal, which we presumably just scanned in,
1152      * and find out what its value is, as a double.
1153      *
1154      * AS A SIDE EFFECT, SET roundDir TO INDICATE PREFERRED
1155      * ROUNDING DIRECTION in case the result is really destined
1156      * for a single-precision float.
1157      */
1158 
1159     public strictfp double doubleValue(){
1160         int     kDigits = Math.min( nDigits, maxDecimalDigits+1 );
1161         long    lValue;
1162         double  dValue;
1163         double  rValue, tValue;
1164 
1165         // First, check for NaN and Infinity values
1166         if(digits == infinity || digits == notANumber) {
1167             if(digits == notANumber)
1168                 return Double.NaN;
1169             else
1170                 return (isNegative?Double.NEGATIVE_INFINITY:Double.POSITIVE_INFINITY);
1171         }
1172         else {
1173             if (mustSetRoundDir) {
1174                 roundDir = 0;
1175             }
1176             /*
1177              * convert the lead kDigits to a long integer.
1178              */
1179             // (special performance hack: start to do it using int)
1180             int iValue = (int)digits[0]-(int)'0';
1181             int iDigits = Math.min( kDigits, intDecimalDigits );
1182             for ( int i=1; i < iDigits; i++ ){
1183                 iValue = iValue*10 + (int)digits[i]-(int)'0';
1184             }
1185             lValue = (long)iValue;
1186             for ( int i=iDigits; i < kDigits; i++ ){
1187                 lValue = lValue*10L + (long)((int)digits[i]-(int)'0');
1188             }
1189             dValue = (double)lValue;
1190             int exp = decExponent-kDigits;
1191             /*
1192              * lValue now contains a long integer with the value of
1193              * the first kDigits digits of the number.
1194              * dValue contains the (double) of the same.
1195              */
1196 
1197             if ( nDigits <= maxDecimalDigits ){
1198                 /*
1199                  * possibly an easy case.
1200                  * We know that the digits can be represented
1201                  * exactly. And if the exponent isn't too outrageous,
1202                  * the whole thing can be done with one operation,
1203                  * thus one rounding error.
1204                  * Note that all our constructors trim all leading and
1205                  * trailing zeros, so simple values (including zero)
1206                  * will always end up here
1207                  */
1208                 if (exp == 0 || dValue == 0.0)
1209                     return (isNegative)? -dValue : dValue; // small floating integer
1210                 else if ( exp >= 0 ){
1211                     if ( exp <= maxSmallTen ){
1212                         /*
1213                          * Can get the answer with one operation,
1214                          * thus one roundoff.
1215                          */
1216                         rValue = dValue * small10pow[exp];
1217                         if ( mustSetRoundDir ){
1218                             tValue = rValue / small10pow[exp];
1219                             roundDir = ( tValue ==  dValue ) ? 0
1220                                 :( tValue < dValue ) ? 1
1221                                 : -1;
1222                         }
1223                         return (isNegative)? -rValue : rValue;
1224                     }
1225                     int slop = maxDecimalDigits - kDigits;
1226                     if ( exp <= maxSmallTen+slop ){
1227                         /*
1228                          * We can multiply dValue by 10^(slop)
1229                          * and it is still "small" and exact.
1230                          * Then we can multiply by 10^(exp-slop)
1231                          * with one rounding.
1232                          */
1233                         dValue *= small10pow[slop];
1234                         rValue = dValue * small10pow[exp-slop];
1235 
1236                         if ( mustSetRoundDir ){
1237                             tValue = rValue / small10pow[exp-slop];
1238                             roundDir = ( tValue ==  dValue ) ? 0
1239                                 :( tValue < dValue ) ? 1
1240                                 : -1;
1241                         }
1242                         return (isNegative)? -rValue : rValue;
1243                     }
1244                     /*
1245                      * Else we have a hard case with a positive exp.
1246                      */
1247                 } else {
1248                     if ( exp >= -maxSmallTen ){
1249                         /*
1250                          * Can get the answer in one division.
1251                          */
1252                         rValue = dValue / small10pow[-exp];
1253                         tValue = rValue * small10pow[-exp];
1254                         if ( mustSetRoundDir ){
1255                             roundDir = ( tValue ==  dValue ) ? 0
1256                                 :( tValue < dValue ) ? 1
1257                                 : -1;
1258                         }
1259                         return (isNegative)? -rValue : rValue;
1260                     }
1261                     /*
1262                      * Else we have a hard case with a negative exp.
1263                      */
1264                 }
1265             }
1266 
1267             /*
1268              * Harder cases:
1269              * The sum of digits plus exponent is greater than
1270              * what we think we can do with one error.
1271              *
1272              * Start by approximating the right answer by,
1273              * naively, scaling by powers of 10.
1274              */
1275             if ( exp > 0 ){
1276                 if ( decExponent > maxDecimalExponent+1 ){
1277                     /*
1278                      * Lets face it. This is going to be
1279                      * Infinity. Cut to the chase.
1280                      */
1281                     return (isNegative)? Double.NEGATIVE_INFINITY : Double.POSITIVE_INFINITY;
1282                 }
1283                 if ( (exp&15) != 0 ){
1284                     dValue *= small10pow[exp&15];
1285                 }
1286                 if ( (exp>>=4) != 0 ){
1287                     int j;
1288                     for( j = 0; exp > 1; j++, exp>>=1 ){
1289                         if ( (exp&1)!=0)
1290                             dValue *= big10pow[j];
1291                     }
1292                     /*
1293                      * The reason for the weird exp > 1 condition
1294                      * in the above loop was so that the last multiply
1295                      * would get unrolled. We handle it here.
1296                      * It could overflow.
1297                      */
1298                     double t = dValue * big10pow[j];
1299                     if ( Double.isInfinite( t ) ){
1300                         /*
1301                          * It did overflow.
1302                          * Look more closely at the result.
1303                          * If the exponent is just one too large,
1304                          * then use the maximum finite as our estimate
1305                          * value. Else call the result infinity
1306                          * and punt it.
1307                          * ( I presume this could happen because
1308                          * rounding forces the result here to be
1309                          * an ULP or two larger than
1310                          * Double.MAX_VALUE ).
1311                          */
1312                         t = dValue / 2.0;
1313                         t *= big10pow[j];
1314                         if ( Double.isInfinite( t ) ){
1315                             return (isNegative)? Double.NEGATIVE_INFINITY : Double.POSITIVE_INFINITY;
1316                         }
1317                         t = Double.MAX_VALUE;
1318                     }
1319                     dValue = t;
1320                 }
1321             } else if ( exp < 0 ){
1322                 exp = -exp;
1323                 if ( decExponent < minDecimalExponent-1 ){
1324                     /*
1325                      * Lets face it. This is going to be
1326                      * zero. Cut to the chase.
1327                      */
1328                     return (isNegative)? -0.0 : 0.0;
1329                 }
1330                 if ( (exp&15) != 0 ){
1331                     dValue /= small10pow[exp&15];
1332                 }
1333                 if ( (exp>>=4) != 0 ){
1334                     int j;
1335                     for( j = 0; exp > 1; j++, exp>>=1 ){
1336                         if ( (exp&1)!=0)
1337                             dValue *= tiny10pow[j];
1338                     }
1339                     /*
1340                      * The reason for the weird exp > 1 condition
1341                      * in the above loop was so that the last multiply
1342                      * would get unrolled. We handle it here.
1343                      * It could underflow.
1344                      */
1345                     double t = dValue * tiny10pow[j];
1346                     if ( t == 0.0 ){
1347                         /*
1348                          * It did underflow.
1349                          * Look more closely at the result.
1350                          * If the exponent is just one too small,
1351                          * then use the minimum finite as our estimate
1352                          * value. Else call the result 0.0
1353                          * and punt it.
1354                          * ( I presume this could happen because
1355                          * rounding forces the result here to be
1356                          * an ULP or two less than
1357                          * Double.MIN_VALUE ).
1358                          */
1359                         t = dValue * 2.0;
1360                         t *= tiny10pow[j];
1361                         if ( t == 0.0 ){
1362                             return (isNegative)? -0.0 : 0.0;
1363                         }
1364                         t = Double.MIN_VALUE;
1365                     }
1366                     dValue = t;
1367                 }
1368             }
1369 
1370             /*
1371              * dValue is now approximately the result.
1372              * The hard part is adjusting it, by comparison
1373              * with FDBigInt arithmetic.
1374              * Formulate the EXACT big-number result as
1375              * bigD0 * 10^exp
1376              */
1377             FDBigInt bigD0 = new FDBigInt( lValue, digits, kDigits, nDigits );
1378             exp   = decExponent - nDigits;
1379 
1380             correctionLoop:
1381             while(true){
1382                 /* AS A SIDE EFFECT, THIS METHOD WILL SET THE INSTANCE VARIABLES
1383                  * bigIntExp and bigIntNBits
1384                  */
1385                 FDBigInt bigB = doubleToBigInt( dValue );
1386 
1387                 /*
1388                  * Scale bigD, bigB appropriately for
1389                  * big-integer operations.
1390                  * Naively, we multipy by powers of ten
1391                  * and powers of two. What we actually do
1392                  * is keep track of the powers of 5 and
1393                  * powers of 2 we would use, then factor out
1394                  * common divisors before doing the work.
1395                  */
1396                 int B2, B5; // powers of 2, 5 in bigB
1397                 int     D2, D5; // powers of 2, 5 in bigD
1398                 int Ulp2;   // powers of 2 in halfUlp.
1399                 if ( exp >= 0 ){
1400                     B2 = B5 = 0;
1401                     D2 = D5 = exp;
1402                 } else {
1403                     B2 = B5 = -exp;
1404                     D2 = D5 = 0;
1405                 }
1406                 if ( bigIntExp >= 0 ){
1407                     B2 += bigIntExp;
1408                 } else {
1409                     D2 -= bigIntExp;
1410                 }
1411                 Ulp2 = B2;
1412                 // shift bigB and bigD left by a number s. t.
1413                 // halfUlp is still an integer.
1414                 int hulpbias;
1415                 if ( bigIntExp+bigIntNBits <= -expBias+1 ){
1416                     // This is going to be a denormalized number
1417                     // (if not actually zero).
1418                     // half an ULP is at 2^-(expBias+expShift+1)
1419                     hulpbias = bigIntExp+ expBias + expShift;
1420                 } else {
1421                     hulpbias = expShift + 2 - bigIntNBits;
1422                 }
1423                 B2 += hulpbias;
1424                 D2 += hulpbias;
1425                 // if there are common factors of 2, we might just as well
1426                 // factor them out, as they add nothing useful.
1427                 int common2 = Math.min( B2, Math.min( D2, Ulp2 ) );
1428                 B2 -= common2;
1429                 D2 -= common2;
1430                 Ulp2 -= common2;
1431                 // do multiplications by powers of 5 and 2
1432                 bigB = multPow52( bigB, B5, B2 );
1433                 FDBigInt bigD = multPow52( new FDBigInt( bigD0 ), D5, D2 );
1434                 //
1435                 // to recap:
1436                 // bigB is the scaled-big-int version of our floating-point
1437                 // candidate.
1438                 // bigD is the scaled-big-int version of the exact value
1439                 // as we understand it.
1440                 // halfUlp is 1/2 an ulp of bigB, except for special cases
1441                 // of exact powers of 2
1442                 //
1443                 // the plan is to compare bigB with bigD, and if the difference
1444                 // is less than halfUlp, then we're satisfied. Otherwise,
1445                 // use the ratio of difference to halfUlp to calculate a fudge
1446                 // factor to add to the floating value, then go 'round again.
1447                 //
1448                 FDBigInt diff;
1449                 int cmpResult;
1450                 boolean overvalue;
1451                 if ( (cmpResult = bigB.cmp( bigD ) ) > 0 ){
1452                     overvalue = true; // our candidate is too big.
1453                     diff = bigB.sub( bigD );
1454                     if ( (bigIntNBits == 1) && (bigIntExp > -expBias) ){
1455                         // candidate is a normalized exact power of 2 and
1456                         // is too big. We will be subtracting.
1457                         // For our purposes, ulp is the ulp of the
1458                         // next smaller range.
1459                         Ulp2 -= 1;
1460                         if ( Ulp2 < 0 ){
1461                             // rats. Cannot de-scale ulp this far.
1462                             // must scale diff in other direction.
1463                             Ulp2 = 0;
1464                             diff.lshiftMe( 1 );
1465                         }
1466                     }
1467                 } else if ( cmpResult < 0 ){
1468                     overvalue = false; // our candidate is too small.
1469                     diff = bigD.sub( bigB );
1470                 } else {
1471                     // the candidate is exactly right!
1472                     // this happens with surprising fequency
1473                     break correctionLoop;
1474                 }
1475                 FDBigInt halfUlp = constructPow52( B5, Ulp2 );
1476                 if ( (cmpResult = diff.cmp( halfUlp ) ) < 0 ){
1477                     // difference is small.
1478                     // this is close enough
1479                     if (mustSetRoundDir) {
1480                         roundDir = overvalue ? -1 : 1;
1481                     }
1482                     break correctionLoop;
1483                 } else if ( cmpResult == 0 ){
1484                     // difference is exactly half an ULP
1485                     // round to some other value maybe, then finish
1486                     dValue += 0.5*ulp( dValue, overvalue );
1487                     // should check for bigIntNBits == 1 here??
1488                     if (mustSetRoundDir) {
1489                         roundDir = overvalue ? -1 : 1;
1490                     }
1491                     break correctionLoop;
1492                 } else {
1493                     // difference is non-trivial.
1494                     // could scale addend by ratio of difference to
1495                     // halfUlp here, if we bothered to compute that difference.
1496                     // Most of the time ( I hope ) it is about 1 anyway.
1497                     dValue += ulp( dValue, overvalue );
1498                     if ( dValue == 0.0 || dValue == Double.POSITIVE_INFINITY )
1499                         break correctionLoop; // oops. Fell off end of range.
1500                     continue; // try again.
1501                 }
1502 
1503             }
1504             return (isNegative)? -dValue : dValue;
1505         }
1506     }
1507 
1508     /*
1509      * Take a FormattedFloatingDecimal, which we presumably just scanned in,
1510      * and find out what its value is, as a float.
1511      * This is distinct from doubleValue() to avoid the extremely
1512      * unlikely case of a double rounding error, wherein the converstion
1513      * to double has one rounding error, and the conversion of that double
1514      * to a float has another rounding error, IN THE WRONG DIRECTION,
1515      * ( because of the preference to a zero low-order bit ).
1516      */
1517 
1518     public strictfp float floatValue(){
1519         int     kDigits = Math.min( nDigits, singleMaxDecimalDigits+1 );
1520         int     iValue;
1521         float   fValue;
1522 
1523         // First, check for NaN and Infinity values
1524         if(digits == infinity || digits == notANumber) {
1525             if(digits == notANumber)
1526                 return Float.NaN;
1527             else
1528                 return (isNegative?Float.NEGATIVE_INFINITY:Float.POSITIVE_INFINITY);
1529         }
1530         else {
1531             /*
1532              * convert the lead kDigits to an integer.
1533              */
1534             iValue = (int)digits[0]-(int)'0';
1535             for ( int i=1; i < kDigits; i++ ){
1536                 iValue = iValue*10 + (int)digits[i]-(int)'0';
1537             }
1538             fValue = (float)iValue;
1539             int exp = decExponent-kDigits;
1540             /*
1541              * iValue now contains an integer with the value of
1542              * the first kDigits digits of the number.
1543              * fValue contains the (float) of the same.
1544              */
1545 
1546             if ( nDigits <= singleMaxDecimalDigits ){
1547                 /*
1548                  * possibly an easy case.
1549                  * We know that the digits can be represented
1550                  * exactly. And if the exponent isn't too outrageous,
1551                  * the whole thing can be done with one operation,
1552                  * thus one rounding error.
1553                  * Note that all our constructors trim all leading and
1554                  * trailing zeros, so simple values (including zero)
1555                  * will always end up here.
1556                  */
1557                 if (exp == 0 || fValue == 0.0f)
1558                     return (isNegative)? -fValue : fValue; // small floating integer
1559                 else if ( exp >= 0 ){
1560                     if ( exp <= singleMaxSmallTen ){
1561                         /*
1562                          * Can get the answer with one operation,
1563                          * thus one roundoff.
1564                          */
1565                         fValue *= singleSmall10pow[exp];
1566                         return (isNegative)? -fValue : fValue;
1567                     }
1568                     int slop = singleMaxDecimalDigits - kDigits;
1569                     if ( exp <= singleMaxSmallTen+slop ){
1570                         /*
1571                          * We can multiply dValue by 10^(slop)
1572                          * and it is still "small" and exact.
1573                          * Then we can multiply by 10^(exp-slop)
1574                          * with one rounding.
1575                          */
1576                         fValue *= singleSmall10pow[slop];
1577                         fValue *= singleSmall10pow[exp-slop];
1578                         return (isNegative)? -fValue : fValue;
1579                     }
1580                     /*
1581                      * Else we have a hard case with a positive exp.
1582                      */
1583                 } else {
1584                     if ( exp >= -singleMaxSmallTen ){
1585                         /*
1586                          * Can get the answer in one division.
1587                          */
1588                         fValue /= singleSmall10pow[-exp];
1589                         return (isNegative)? -fValue : fValue;
1590                     }
1591                     /*
1592                      * Else we have a hard case with a negative exp.
1593                      */
1594                 }
1595             } else if ( (decExponent >= nDigits) && (nDigits+decExponent <= maxDecimalDigits) ){
1596                 /*
1597                  * In double-precision, this is an exact floating integer.
1598                  * So we can compute to double, then shorten to float
1599                  * with one round, and get the right answer.
1600                  *
1601                  * First, finish accumulating digits.
1602                  * Then convert that integer to a double, multiply
1603                  * by the appropriate power of ten, and convert to float.
1604                  */
1605                 long lValue = (long)iValue;
1606                 for ( int i=kDigits; i < nDigits; i++ ){
1607                     lValue = lValue*10L + (long)((int)digits[i]-(int)'0');
1608                 }
1609                 double dValue = (double)lValue;
1610                 exp = decExponent-nDigits;
1611                 dValue *= small10pow[exp];
1612                 fValue = (float)dValue;
1613                 return (isNegative)? -fValue : fValue;
1614 
1615             }
1616             /*
1617              * Harder cases:
1618              * The sum of digits plus exponent is greater than
1619              * what we think we can do with one error.
1620              *
1621              * Start by weeding out obviously out-of-range
1622              * results, then convert to double and go to
1623              * common hard-case code.
1624              */
1625             if ( decExponent > singleMaxDecimalExponent+1 ){
1626                 /*
1627                  * Lets face it. This is going to be
1628                  * Infinity. Cut to the chase.
1629                  */
1630                 return (isNegative)? Float.NEGATIVE_INFINITY : Float.POSITIVE_INFINITY;
1631             } else if ( decExponent < singleMinDecimalExponent-1 ){
1632                 /*
1633                  * Lets face it. This is going to be
1634                  * zero. Cut to the chase.
1635                  */
1636                 return (isNegative)? -0.0f : 0.0f;
1637             }
1638 
1639             /*
1640              * Here, we do 'way too much work, but throwing away
1641              * our partial results, and going and doing the whole
1642              * thing as double, then throwing away half the bits that computes
1643              * when we convert back to float.
1644              *
1645              * The alternative is to reproduce the whole multiple-precision
1646              * algorythm for float precision, or to try to parameterize it
1647              * for common usage. The former will take about 400 lines of code,
1648              * and the latter I tried without success. Thus the semi-hack
1649              * answer here.
1650              */
1651             mustSetRoundDir = !fromHex;
1652             double dValue = doubleValue();
1653             return stickyRound( dValue );
1654         }
1655     }
1656 
1657 
1658     /*
1659      * All the positive powers of 10 that can be
1660      * represented exactly in double/float.
1661      */
1662     private static final double small10pow[] = {
1663         1.0e0,
1664         1.0e1, 1.0e2, 1.0e3, 1.0e4, 1.0e5,
1665         1.0e6, 1.0e7, 1.0e8, 1.0e9, 1.0e10,
1666         1.0e11, 1.0e12, 1.0e13, 1.0e14, 1.0e15,
1667         1.0e16, 1.0e17, 1.0e18, 1.0e19, 1.0e20,
1668         1.0e21, 1.0e22
1669     };
1670 
1671     private static final float singleSmall10pow[] = {
1672         1.0e0f,
1673         1.0e1f, 1.0e2f, 1.0e3f, 1.0e4f, 1.0e5f,
1674         1.0e6f, 1.0e7f, 1.0e8f, 1.0e9f, 1.0e10f
1675     };
1676 
1677     private static final double big10pow[] = {
1678         1e16, 1e32, 1e64, 1e128, 1e256 };
1679     private static final double tiny10pow[] = {
1680         1e-16, 1e-32, 1e-64, 1e-128, 1e-256 };
1681 
1682     private static final int maxSmallTen = small10pow.length-1;
1683     private static final int singleMaxSmallTen = singleSmall10pow.length-1;
1684 
1685     private static final int small5pow[] = {
1686         1,
1687         5,
1688         5*5,
1689         5*5*5,
1690         5*5*5*5,
1691         5*5*5*5*5,
1692         5*5*5*5*5*5,
1693         5*5*5*5*5*5*5,
1694         5*5*5*5*5*5*5*5,
1695         5*5*5*5*5*5*5*5*5,
1696         5*5*5*5*5*5*5*5*5*5,
1697         5*5*5*5*5*5*5*5*5*5*5,
1698         5*5*5*5*5*5*5*5*5*5*5*5,
1699         5*5*5*5*5*5*5*5*5*5*5*5*5
1700     };
1701 
1702 
1703     private static final long long5pow[] = {
1704         1L,
1705         5L,
1706         5L*5,
1707         5L*5*5,
1708         5L*5*5*5,
1709         5L*5*5*5*5,
1710         5L*5*5*5*5*5,
1711         5L*5*5*5*5*5*5,
1712         5L*5*5*5*5*5*5*5,
1713         5L*5*5*5*5*5*5*5*5,
1714         5L*5*5*5*5*5*5*5*5*5,
1715         5L*5*5*5*5*5*5*5*5*5*5,
1716         5L*5*5*5*5*5*5*5*5*5*5*5,
1717         5L*5*5*5*5*5*5*5*5*5*5*5*5,
1718         5L*5*5*5*5*5*5*5*5*5*5*5*5*5,
1719         5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
1720         5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
1721         5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
1722         5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
1723         5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
1724         5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
1725         5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
1726         5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
1727         5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
1728         5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
1729         5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
1730         5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
1731     };
1732 
1733     // approximately ceil( log2( long5pow[i] ) )
1734     private static final int n5bits[] = {
1735         0,
1736         3,
1737         5,
1738         7,
1739         10,
1740         12,
1741         14,
1742         17,
1743         19,
1744         21,
1745         24,
1746         26,
1747         28,
1748         31,
1749         33,
1750         35,
1751         38,
1752         40,
1753         42,
1754         45,
1755         47,
1756         49,
1757         52,
1758         54,
1759         56,
1760         59,
1761         61,
1762     };
1763 
1764     private static final char infinity[] = { 'I', 'n', 'f', 'i', 'n', 'i', 't', 'y' };
1765     private static final char notANumber[] = { 'N', 'a', 'N' };
1766     private static final char zero[] = { '0', '0', '0', '0', '0', '0', '0', '0' };
1767 }
--- EOF ---