com.jfinal.template.io.FloatingDecimal Maven / Gradle / Ivy
The newest version!
/*
* Copyright (c) 1996, 2011, Oracle and/or its affiliates. All rights reserved.
* DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER.
*/
package com.jfinal.template.io;
public class FloatingDecimal{
boolean isExceptional;
boolean isNegative;
int decExponent;
char digits[];
int nDigits;
int bigIntExp;
int bigIntNBits;
boolean mustSetRoundDir = false;
boolean fromHex = false;
int roundDir = 0; // set by doubleValue
/*
* Constants of the implementation
* Most are IEEE-754 related.
* (There are more really boring constants at the end.)
*/
static final long signMask = 0x8000000000000000L;
static final long expMask = 0x7ff0000000000000L;
static final long fractMask= ~(signMask|expMask);
static final int expShift = 52;
static final int expBias = 1023;
static final long fractHOB = ( 1L< 0L ) { // i.e. while ((v&highbit) == 0L )
v <<= 1;
}
int n = 0;
while (( v & lowbytes ) != 0L ){
v <<= 8;
n += 8;
}
while ( v != 0L ){
v <<= 1;
n += 1;
}
return n;
}
/*
* Keep big powers of 5 handy for future reference.
*/
private static FDBigInt b5p[];
private static synchronized FDBigInt
big5pow( int p ){
assert p >= 0 : p; // negative power of 5
if ( b5p == null ){
b5p = new FDBigInt[ p+1 ];
}else if (b5p.length <= p ){
FDBigInt t[] = new FDBigInt[ p+1 ];
System.arraycopy( b5p, 0, t, 0, b5p.length );
b5p = t;
}
if ( b5p[p] != null )
return b5p[p];
else if ( p < small5pow.length )
return b5p[p] = new FDBigInt( small5pow[p] );
else if ( p < long5pow.length )
return b5p[p] = new FDBigInt( long5pow[p] );
else {
// construct the value.
// recursively.
int q, r;
// in order to compute 5^p,
// compute its square root, 5^(p/2) and square.
// or, let q = p / 2, r = p -q, then
// 5^p = 5^(q+r) = 5^q * 5^r
q = p >> 1;
r = p - q;
FDBigInt bigq = b5p[q];
if ( bigq == null )
bigq = big5pow ( q );
if ( r < small5pow.length ){
return (b5p[p] = bigq.mult( small5pow[r] ) );
}else{
FDBigInt bigr = b5p[ r ];
if ( bigr == null )
bigr = big5pow( r );
return (b5p[p] = bigq.mult( bigr ) );
}
}
}
//
// a common operation
//
private static FDBigInt
multPow52( FDBigInt v, int p5, int p2 ){
if ( p5 != 0 ){
if ( p5 < small5pow.length ){
v = v.mult( small5pow[p5] );
} else {
v = v.mult( big5pow( p5 ) );
}
}
if ( p2 != 0 ){
v.lshiftMe( p2 );
}
return v;
}
//
// another common operation
//
private static FDBigInt
constructPow52( int p5, int p2 ){
FDBigInt v = new FDBigInt( big5pow( p5 ) );
if ( p2 != 0 ){
v.lshiftMe( p2 );
}
return v;
}
/*
* This is the easy subcase --
* all the significant bits, after scaling, are held in lvalue.
* negSign and decExponent tell us what processing and scaling
* has already been done. Exceptional cases have already been
* stripped out.
* In particular:
* lvalue is a finite number (not Inf, nor NaN)
* lvalue > 0L (not zero, nor negative).
*
* The only reason that we develop the digits here, rather than
* calling on Long.toString() is that we can do it a little faster,
* and besides want to treat trailing 0s specially. If Long.toString
* changes, we should re-evaluate this strategy!
*/
private void
developLongDigits( int decExponent, long lvalue, long insignificant ){
char digits[];
int ndigits;
int digitno;
int c;
//
// Discard non-significant low-order bits, while rounding,
// up to insignificant value.
int i;
for ( i = 0; insignificant >= 10L; i++ )
insignificant /= 10L;
if ( i != 0 ){
long pow10 = long5pow[i] << i; // 10^i == 5^i * 2^i;
long residue = lvalue % pow10;
lvalue /= pow10;
decExponent += i;
if ( residue >= (pow10>>1) ){
// round up based on the low-order bits we're discarding
lvalue++;
}
}
if ( lvalue <= Integer.MAX_VALUE ){
assert lvalue > 0L : lvalue; // lvalue <= 0
// even easier subcase!
// can do int arithmetic rather than long!
int ivalue = (int)lvalue;
ndigits = 10;
digits = (char[])(perThreadBuffer.get());
digitno = ndigits-1;
c = ivalue%10;
ivalue /= 10;
while ( c == 0 ){
decExponent++;
c = ivalue%10;
ivalue /= 10;
}
while ( ivalue != 0){
digits[digitno--] = (char)(c+'0');
decExponent++;
c = ivalue%10;
ivalue /= 10;
}
digits[digitno] = (char)(c+'0');
} else {
// same algorithm as above (same bugs, too )
// but using long arithmetic.
ndigits = 20;
digits = (char[])(perThreadBuffer.get());
digitno = ndigits-1;
c = (int)(lvalue%10L);
lvalue /= 10L;
while ( c == 0 ){
decExponent++;
c = (int)(lvalue%10L);
lvalue /= 10L;
}
while ( lvalue != 0L ){
digits[digitno--] = (char)(c+'0');
decExponent++;
c = (int)(lvalue%10L);
lvalue /= 10;
}
digits[digitno] = (char)(c+'0');
}
char result [];
ndigits -= digitno;
result = new char[ ndigits ];
System.arraycopy( digits, digitno, result, 0, ndigits );
this.digits = result;
this.decExponent = decExponent+1;
this.nDigits = ndigits;
}
//
// add one to the least significant digit.
// in the unlikely event there is a carry out,
// deal with it.
// assert that this will only happen where there
// is only one digit, e.g. (float)1e-44 seems to do it.
//
private void
roundup(){
int i;
int q = digits[ i = (nDigits-1)];
if ( q == '9' ){
while ( q == '9' && i > 0 ){
digits[i] = '0';
q = digits[--i];
}
if ( q == '9' ){
// carryout! High-order 1, rest 0s, larger exp.
decExponent += 1;
digits[0] = '1';
return;
}
// else fall through.
}
digits[i] = (char)(q+1);
}
/*
* FIRST IMPORTANT CONSTRUCTOR: DOUBLE
*/
public FloatingDecimal( double d )
{
long dBits = Double.doubleToLongBits( d );
long fractBits;
int binExp;
int nSignificantBits;
// discover and delete sign
if ( (dBits&signMask) != 0 ){
isNegative = true;
dBits ^= signMask;
} else {
isNegative = false;
}
// Begin to unpack
// Discover obvious special cases of NaN and Infinity.
binExp = (int)( (dBits&expMask) >> expShift );
fractBits = dBits&fractMask;
if ( binExp == (int)(expMask>>expShift) ) {
isExceptional = true;
if ( fractBits == 0L ){
digits = infinity;
} else {
digits = notANumber;
isNegative = false; // NaN has no sign!
}
nDigits = digits.length;
return;
}
isExceptional = false;
// Finish unpacking
// Normalize denormalized numbers.
// Insert assumed high-order bit for normalized numbers.
// Subtract exponent bias.
if ( binExp == 0 ){
if ( fractBits == 0L ){
// not a denorm, just a 0!
decExponent = 0;
digits = zero;
nDigits = 1;
return;
}
while ( (fractBits&fractHOB) == 0L ){
fractBits <<= 1;
binExp -= 1;
}
nSignificantBits = expShift + binExp +1; // recall binExp is - shift count.
binExp += 1;
} else {
fractBits |= fractHOB;
nSignificantBits = expShift+1;
}
binExp -= expBias;
// call the routine that actually does all the hard work.
dtoa( binExp, fractBits, nSignificantBits );
}
/*
* SECOND IMPORTANT CONSTRUCTOR: SINGLE
*/
public FloatingDecimal( float f )
{
int fBits = Float.floatToIntBits( f );
int fractBits;
int binExp;
int nSignificantBits;
// discover and delete sign
if ( (fBits&singleSignMask) != 0 ){
isNegative = true;
fBits ^= singleSignMask;
} else {
isNegative = false;
}
// Begin to unpack
// Discover obvious special cases of NaN and Infinity.
binExp = (int)( (fBits&singleExpMask) >> singleExpShift );
fractBits = fBits&singleFractMask;
if ( binExp == (int)(singleExpMask>>singleExpShift) ) {
isExceptional = true;
if ( fractBits == 0L ){
digits = infinity;
} else {
digits = notANumber;
isNegative = false; // NaN has no sign!
}
nDigits = digits.length;
return;
}
isExceptional = false;
// Finish unpacking
// Normalize denormalized numbers.
// Insert assumed high-order bit for normalized numbers.
// Subtract exponent bias.
if ( binExp == 0 ){
if ( fractBits == 0 ){
// not a denorm, just a 0!
decExponent = 0;
digits = zero;
nDigits = 1;
return;
}
while ( (fractBits&singleFractHOB) == 0 ){
fractBits <<= 1;
binExp -= 1;
}
nSignificantBits = singleExpShift + binExp +1; // recall binExp is - shift count.
binExp += 1;
} else {
fractBits |= singleFractHOB;
nSignificantBits = singleExpShift+1;
}
binExp -= singleExpBias;
// call the routine that actually does all the hard work.
dtoa( binExp, ((long)fractBits)<<(expShift-singleExpShift), nSignificantBits );
}
private void
dtoa( int binExp, long fractBits, int nSignificantBits )
{
int nFractBits; // number of significant bits of fractBits;
int nTinyBits; // number of these to the right of the point.
int decExp;
// Examine number. Determine if it is an easy case,
// which we can do pretty trivially using float/long conversion,
// or whether we must do real work.
nFractBits = countBits( fractBits );
nTinyBits = Math.max( 0, nFractBits - binExp - 1 );
if ( binExp <= maxSmallBinExp && binExp >= minSmallBinExp ){
// Look more closely at the number to decide if,
// with scaling by 10^nTinyBits, the result will fit in
// a long.
if ( (nTinyBits < long5pow.length) && ((nFractBits + n5bits[nTinyBits]) < 64 ) ){
/*
* We can do this:
* take the fraction bits, which are normalized.
* (a) nTinyBits == 0: Shift left or right appropriately
* to align the binary point at the extreme right, i.e.
* where a long int point is expected to be. The integer
* result is easily converted to a string.
* (b) nTinyBits > 0: Shift right by expShift-nFractBits,
* which effectively converts to long and scales by
* 2^nTinyBits. Then multiply by 5^nTinyBits to
* complete the scaling. We know this won't overflow
* because we just counted the number of bits necessary
* in the result. The integer you get from this can
* then be converted to a string pretty easily.
*/
long halfULP;
if ( nTinyBits == 0 ) {
if ( binExp > nSignificantBits ){
halfULP = 1L << ( binExp-nSignificantBits-1);
} else {
halfULP = 0L;
}
if ( binExp >= expShift ){
fractBits <<= (binExp-expShift);
} else {
fractBits >>>= (expShift-binExp) ;
}
developLongDigits( 0, fractBits, halfULP );
return;
}
/*
* The following causes excess digits to be printed
* out in the single-float case. Our manipulation of
* halfULP here is apparently not correct. If we
* better understand how this works, perhaps we can
* use this special case again. But for the time being,
* we do not.
* else {
* fractBits >>>= expShift+1-nFractBits;
* fractBits *= long5pow[ nTinyBits ];
* halfULP = long5pow[ nTinyBits ] >> (1+nSignificantBits-nFractBits);
* developLongDigits( -nTinyBits, fractBits, halfULP );
* return;
* }
*/
}
}
/*
* This is the hard case. We are going to compute large positive
* integers B and S and integer decExp, s.t.
* d = ( B / S ) * 10^decExp
* 1 <= B / S < 10
* Obvious choices are:
* decExp = floor( log10(d) )
* B = d * 2^nTinyBits * 10^max( 0, -decExp )
* S = 10^max( 0, decExp) * 2^nTinyBits
* (noting that nTinyBits has already been forced to non-negative)
* I am also going to compute a large positive integer
* M = (1/2^nSignificantBits) * 2^nTinyBits * 10^max( 0, -decExp )
* i.e. M is (1/2) of the ULP of d, scaled like B.
* When we iterate through dividing B/S and picking off the
* quotient bits, we will know when to stop when the remainder
* is <= M.
*
* We keep track of powers of 2 and powers of 5.
*/
/*
* Estimate decimal exponent. (If it is small-ish,
* we could double-check.)
*
* First, scale the mantissa bits such that 1 <= d2 < 2.
* We are then going to estimate
* log10(d2) ~=~ (d2-1.5)/1.5 + log(1.5)
* and so we can estimate
* log10(d) ~=~ log10(d2) + binExp * log10(2)
* take the floor and call it decExp.
* FIXME -- use more precise constants here. It costs no more.
*/
double d2 = Double.longBitsToDouble(
expOne | ( fractBits &~ fractHOB ) );
decExp = (int)Math.floor(
(d2-1.5D)*0.289529654D + 0.176091259 + (double)binExp * 0.301029995663981 );
int B2, B5; // powers of 2 and powers of 5, respectively, in B
int S2, S5; // powers of 2 and powers of 5, respectively, in S
int M2, M5; // powers of 2 and powers of 5, respectively, in M
int Bbits; // binary digits needed to represent B, approx.
int tenSbits; // binary digits needed to represent 10*S, approx.
FDBigInt Sval, Bval, Mval;
B5 = Math.max( 0, -decExp );
B2 = B5 + nTinyBits + binExp;
S5 = Math.max( 0, decExp );
S2 = S5 + nTinyBits;
M5 = B5;
M2 = B2 - nSignificantBits;
/*
* the long integer fractBits contains the (nFractBits) interesting
* bits from the mantissa of d ( hidden 1 added if necessary) followed
* by (expShift+1-nFractBits) zeros. In the interest of compactness,
* I will shift out those zeros before turning fractBits into a
* FDBigInt. The resulting whole number will be
* d * 2^(nFractBits-1-binExp).
*/
fractBits >>>= (expShift+1-nFractBits);
B2 -= nFractBits-1;
int common2factor = Math.min( B2, S2 );
B2 -= common2factor;
S2 -= common2factor;
M2 -= common2factor;
/*
* HACK!! For exact powers of two, the next smallest number
* is only half as far away as we think (because the meaning of
* ULP changes at power-of-two bounds) for this reason, we
* hack M2. Hope this works.
*/
if ( nFractBits == 1 )
M2 -= 1;
if ( M2 < 0 ){
// oops.
// since we cannot scale M down far enough,
// we must scale the other values up.
B2 -= M2;
S2 -= M2;
M2 = 0;
}
/*
* Construct, Scale, iterate.
* Some day, we'll write a stopping test that takes
* account of the asymmetry of the spacing of floating-point
* numbers below perfect powers of 2
* 26 Sept 96 is not that day.
* So we use a symmetric test.
*/
char digits[] = this.digits = new char[18];
int ndigit = 0;
boolean low, high;
long lowDigitDifference;
int q;
/*
* Detect the special cases where all the numbers we are about
* to compute will fit in int or long integers.
* In these cases, we will avoid doing FDBigInt arithmetic.
* We use the same algorithms, except that we "normalize"
* our FDBigInts before iterating. This is to make division easier,
* as it makes our fist guess (quotient of high-order words)
* more accurate!
*
* Some day, we'll write a stopping test that takes
* account of the asymmetry of the spacing of floating-point
* numbers below perfect powers of 2
* 26 Sept 96 is not that day.
* So we use a symmetric test.
*/
Bbits = nFractBits + B2 + (( B5 < n5bits.length )? n5bits[B5] : ( B5*3 ));
tenSbits = S2+1 + (( (S5+1) < n5bits.length )? n5bits[(S5+1)] : ( (S5+1)*3 ));
if ( Bbits < 64 && tenSbits < 64){
if ( Bbits < 32 && tenSbits < 32){
// wa-hoo! They're all ints!
int b = ((int)fractBits * small5pow[B5] ) << B2;
int s = small5pow[S5] << S2;
int m = small5pow[M5] << M2;
int tens = s * 10;
/*
* Unroll the first iteration. If our decExp estimate
* was too high, our first quotient will be zero. In this
* case, we discard it and decrement decExp.
*/
ndigit = 0;
q = b / s;
b = 10 * ( b % s );
m *= 10;
low = (b < m );
high = (b+m > tens );
assert q < 10 : q; // excessively large digit
if ( (q == 0) && ! high ){
// oops. Usually ignore leading zero.
decExp--;
} else {
digits[ndigit++] = (char)('0' + q);
}
/*
* HACK! Java spec sez that we always have at least
* one digit after the . in either F- or E-form output.
* Thus we will need more than one digit if we're using
* E-form
*/
if ( decExp < -3 || decExp >= 8 ){
high = low = false;
}
while( ! low && ! high ){
q = b / s;
b = 10 * ( b % s );
m *= 10;
assert q < 10 : q; // excessively large digit
if ( m > 0L ){
low = (b < m );
high = (b+m > tens );
} else {
// hack -- m might overflow!
// in this case, it is certainly > b,
// which won't
// and b+m > tens, too, since that has overflowed
// either!
low = true;
high = true;
}
digits[ndigit++] = (char)('0' + q);
}
lowDigitDifference = (b<<1) - tens;
} else {
// still good! they're all longs!
long b = (fractBits * long5pow[B5] ) << B2;
long s = long5pow[S5] << S2;
long m = long5pow[M5] << M2;
long tens = s * 10L;
/*
* Unroll the first iteration. If our decExp estimate
* was too high, our first quotient will be zero. In this
* case, we discard it and decrement decExp.
*/
ndigit = 0;
q = (int) ( b / s );
b = 10L * ( b % s );
m *= 10L;
low = (b < m );
high = (b+m > tens );
assert q < 10 : q; // excessively large digit
if ( (q == 0) && ! high ){
// oops. Usually ignore leading zero.
decExp--;
} else {
digits[ndigit++] = (char)('0' + q);
}
/*
* HACK! Java spec sez that we always have at least
* one digit after the . in either F- or E-form output.
* Thus we will need more than one digit if we're using
* E-form
*/
if ( decExp < -3 || decExp >= 8 ){
high = low = false;
}
while( ! low && ! high ){
q = (int) ( b / s );
b = 10 * ( b % s );
m *= 10;
assert q < 10 : q; // excessively large digit
if ( m > 0L ){
low = (b < m );
high = (b+m > tens );
} else {
// hack -- m might overflow!
// in this case, it is certainly > b,
// which won't
// and b+m > tens, too, since that has overflowed
// either!
low = true;
high = true;
}
digits[ndigit++] = (char)('0' + q);
}
lowDigitDifference = (b<<1) - tens;
}
} else {
FDBigInt tenSval;
int shiftBias;
/*
* We really must do FDBigInt arithmetic.
* Fist, construct our FDBigInt initial values.
*/
Bval = multPow52( new FDBigInt( fractBits ), B5, B2 );
Sval = constructPow52( S5, S2 );
Mval = constructPow52( M5, M2 );
// normalize so that division works better
Bval.lshiftMe( shiftBias = Sval.normalizeMe() );
Mval.lshiftMe( shiftBias );
tenSval = Sval.mult( 10 );
/*
* Unroll the first iteration. If our decExp estimate
* was too high, our first quotient will be zero. In this
* case, we discard it and decrement decExp.
*/
ndigit = 0;
q = Bval.quoRemIteration( Sval );
Mval = Mval.mult( 10 );
low = (Bval.cmp( Mval ) < 0);
high = (Bval.add( Mval ).cmp( tenSval ) > 0 );
assert q < 10 : q; // excessively large digit
if ( (q == 0) && ! high ){
// oops. Usually ignore leading zero.
decExp--;
} else {
digits[ndigit++] = (char)('0' + q);
}
/*
* HACK! Java spec sez that we always have at least
* one digit after the . in either F- or E-form output.
* Thus we will need more than one digit if we're using
* E-form
*/
if ( decExp < -3 || decExp >= 8 ){
high = low = false;
}
while( ! low && ! high ){
q = Bval.quoRemIteration( Sval );
Mval = Mval.mult( 10 );
assert q < 10 : q; // excessively large digit
low = (Bval.cmp( Mval ) < 0);
high = (Bval.add( Mval ).cmp( tenSval ) > 0 );
digits[ndigit++] = (char)('0' + q);
}
if ( high && low ){
Bval.lshiftMe(1);
lowDigitDifference = Bval.cmp(tenSval);
} else
lowDigitDifference = 0L; // this here only for flow analysis!
}
this.decExponent = decExp+1;
this.digits = digits;
this.nDigits = ndigit;
/*
* Last digit gets rounded based on stopping condition.
*/
if ( high ){
if ( low ){
if ( lowDigitDifference == 0L ){
// it's a tie!
// choose based on which digits we like.
if ( (digits[nDigits-1]&1) != 0 ) roundup();
} else if ( lowDigitDifference > 0 ){
roundup();
}
} else {
roundup();
}
}
}
public String
toString(){
// most brain-dead version
StringBuffer result = new StringBuffer( nDigits+8 );
if ( isNegative ){ result.append( '-' ); }
if ( isExceptional ){
result.append( digits, 0, nDigits );
} else {
result.append( "0.");
result.append( digits, 0, nDigits );
result.append('e');
result.append( decExponent );
}
return new String(result);
}
public String toJavaFormatString() {
char result[] = (char[])(perThreadBuffer.get());
int i = getChars(result);
return new String(result, 0, i);
}
public int getChars(char[] result) {
assert nDigits <= 19 : nDigits; // generous bound on size of nDigits
int i = 0;
if (isNegative) { result[0] = '-'; i = 1; }
if (isExceptional) {
System.arraycopy(digits, 0, result, i, nDigits);
i += nDigits;
} else {
if (decExponent > 0 && decExponent < 8) {
// print digits.digits.
int charLength = Math.min(nDigits, decExponent);
System.arraycopy(digits, 0, result, i, charLength);
i += charLength;
if (charLength < decExponent) {
charLength = decExponent-charLength;
System.arraycopy(zero, 0, result, i, charLength);
i += charLength;
result[i++] = '.';
result[i++] = '0';
} else {
result[i++] = '.';
if (charLength < nDigits) {
int t = nDigits - charLength;
System.arraycopy(digits, charLength, result, i, t);
i += t;
} else {
result[i++] = '0';
}
}
} else if (decExponent <=0 && decExponent > -3) {
result[i++] = '0';
result[i++] = '.';
if (decExponent != 0) {
System.arraycopy(zero, 0, result, i, -decExponent);
i -= decExponent;
}
System.arraycopy(digits, 0, result, i, nDigits);
i += nDigits;
} else {
result[i++] = digits[0];
result[i++] = '.';
if (nDigits > 1) {
System.arraycopy(digits, 1, result, i, nDigits-1);
i += nDigits-1;
} else {
result[i++] = '0';
}
result[i++] = 'E';
int e;
if (decExponent <= 0) {
result[i++] = '-';
e = -decExponent+1;
} else {
e = decExponent-1;
}
// decExponent has 1, 2, or 3, digits
if (e <= 9) {
result[i++] = (char)(e+'0');
} else if (e <= 99) {
result[i++] = (char)(e/10 +'0');
result[i++] = (char)(e%10 + '0');
} else {
result[i++] = (char)(e/100+'0');
e %= 100;
result[i++] = (char)(e/10+'0');
result[i++] = (char)(e%10 + '0');
}
}
}
return i;
}
// Per-thread buffer for string/stringbuffer conversion
@SuppressWarnings("rawtypes")
private static ThreadLocal perThreadBuffer = new ThreadLocal() {
protected synchronized Object initialValue() {
return new char[26];
}
};
private static final int small5pow[] = {
1,
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*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*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*5*5*5,
5*5*5*5*5*5*5*5*5*5*5*5*5
};
private static final long long5pow[] = {
1L,
5L,
5L*5,
5L*5*5,
5L*5*5*5,
5L*5*5*5*5,
5L*5*5*5*5*5,
5L*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
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,
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,
};
// approximately ceil( log2( long5pow[i] ) )
private static final int n5bits[] = {
0,
3,
5,
7,
10,
12,
14,
17,
19,
21,
24,
26,
28,
31,
33,
35,
38,
40,
42,
45,
47,
49,
52,
54,
56,
59,
61,
};
private static final char infinity[] = { 'I', 'n', 'f', 'i', 'n', 'i', 't', 'y' };
private static final char notANumber[] = { 'N', 'a', 'N' };
private static final char zero[] = { '0', '0', '0', '0', '0', '0', '0', '0' };
}
/*
* A really, really simple bigint package
* tailored to the needs of floating base conversion.
*/
class FDBigInt {
int nWords; // number of words used
int data[]; // value: data[0] is least significant
public FDBigInt( long v ){
data = new int[2];
data[0] = (int)v;
data[1] = (int)(v>>>32);
nWords = (data[1]==0) ? 1 : 2;
}
public FDBigInt( FDBigInt other ){
data = new int[nWords = other.nWords];
System.arraycopy( other.data, 0, data, 0, nWords );
}
private FDBigInt( int [] d, int n ){
data = d;
nWords = n;
}
/*
* Left shift by c bits.
* Shifts this in place.
*/
public void
lshiftMe( int c )throws IllegalArgumentException {
if ( c <= 0 ){
if ( c == 0 )
return; // silly.
else
throw new IllegalArgumentException("negative shift count");
}
int wordcount = c>>5;
int bitcount = c & 0x1f;
int anticount = 32-bitcount;
int t[] = data;
int s[] = data;
if ( nWords+wordcount+1 > t.length ){
// reallocate.
t = new int[ nWords+wordcount+1 ];
}
int target = nWords+wordcount;
int src = nWords-1;
if ( bitcount == 0 ){
// special hack, since an anticount of 32 won't go!
System.arraycopy( s, 0, t, wordcount, nWords );
target = wordcount-1;
} else {
t[target--] = s[src]>>>anticount;
while ( src >= 1 ){
t[target--] = (s[src]<>>anticount);
}
t[target--] = s[src]<= 0 ){
t[target--] = 0;
}
data = t;
nWords += wordcount + 1;
// may have constructed high-order word of 0.
// if so, trim it
while ( nWords > 1 && data[nWords-1] == 0 )
nWords--;
}
/*
* normalize this number by shifting until
* the MSB of the number is at 0x08000000.
* This is in preparation for quoRemIteration, below.
* The idea is that, to make division easier, we want the
* divisor to be "normalized" -- usually this means shifting
* the MSB into the high words sign bit. But because we know that
* the quotient will be 0 < q < 10, we would like to arrange that
* the dividend not span up into another word of precision.
* (This needs to be explained more clearly!)
*/
public int
normalizeMe() throws IllegalArgumentException {
int src;
int wordcount = 0;
int bitcount = 0;
int v = 0;
for ( src= nWords-1 ; src >= 0 && (v=data[src]) == 0 ; src--){
wordcount += 1;
}
if ( src < 0 ){
// oops. Value is zero. Cannot normalize it!
throw new IllegalArgumentException("zero value");
}
/*
* In most cases, we assume that wordcount is zero. This only
* makes sense, as we try not to maintain any high-order
* words full of zeros. In fact, if there are zeros, we will
* simply SHORTEN our number at this point. Watch closely...
*/
nWords -= wordcount;
/*
* Compute how far left we have to shift v s.t. its highest-
* order bit is in the right place. Then call lshiftMe to
* do the work.
*/
if ( (v & 0xf0000000) != 0 ){
// will have to shift up into the next word.
// too bad.
for( bitcount = 32 ; (v & 0xf0000000) != 0 ; bitcount-- )
v >>>= 1;
} else {
while ( v <= 0x000fffff ){
// hack: byte-at-a-time shifting
v <<= 8;
bitcount += 8;
}
while ( v <= 0x07ffffff ){
v <<= 1;
bitcount += 1;
}
}
if ( bitcount != 0 )
lshiftMe( bitcount );
return bitcount;
}
/*
* Multiply a FDBigInt by an int.
* Result is a new FDBigInt.
*/
public FDBigInt
mult( int iv ) {
long v = iv;
int r[];
long p;
// guess adequate size of r.
r = new int[ ( v * ((long)data[nWords-1]&0xffffffffL) > 0xfffffffL ) ? nWords+1 : nWords ];
p = 0L;
for( int i=0; i < nWords; i++ ) {
p += v * ((long)data[i]&0xffffffffL);
r[i] = (int)p;
p >>>= 32;
}
if ( p == 0L){
return new FDBigInt( r, nWords );
} else {
r[nWords] = (int)p;
return new FDBigInt( r, nWords+1 );
}
}
/*
* Multiply a FDBigInt by another FDBigInt.
* Result is a new FDBigInt.
*/
public FDBigInt
mult( FDBigInt other ){
// crudely guess adequate size for r
int r[] = new int[ nWords + other.nWords ];
int i;
// I think I am promised zeros...
for( i = 0; i < this.nWords; i++ ){
long v = (long)this.data[i] & 0xffffffffL; // UNSIGNED CONVERSION
long p = 0L;
int j;
for( j = 0; j < other.nWords; j++ ){
p += ((long)r[i+j]&0xffffffffL) + v*((long)other.data[j]&0xffffffffL); // UNSIGNED CONVERSIONS ALL 'ROUND.
r[i+j] = (int)p;
p >>>= 32;
}
r[i+j] = (int)p;
}
// compute how much of r we actually needed for all that.
for ( i = r.length-1; i> 0; i--)
if ( r[i] != 0 )
break;
return new FDBigInt( r, i+1 );
}
/*
* Add one FDBigInt to another. Return a FDBigInt
*/
public FDBigInt
add( FDBigInt other ){
int i;
int a[], b[];
int n, m;
long c = 0L;
// arrange such that a.nWords >= b.nWords;
// n = a.nWords, m = b.nWords
if ( this.nWords >= other.nWords ){
a = this.data;
n = this.nWords;
b = other.data;
m = other.nWords;
} else {
a = other.data;
n = other.nWords;
b = this.data;
m = this.nWords;
}
int r[] = new int[ n ];
for ( i = 0; i < n; i++ ){
c += (long)a[i] & 0xffffffffL;
if ( i < m ){
c += (long)b[i] & 0xffffffffL;
}
r[i] = (int) c;
c >>= 32; // signed shift.
}
if ( c != 0L ){
// oops -- carry out -- need longer result.
int s[] = new int[ r.length+1 ];
System.arraycopy( r, 0, s, 0, r.length );
s[i++] = (int)c;
return new FDBigInt( s, i );
}
return new FDBigInt( r, i );
}
/*
* Compare FDBigInt with another FDBigInt. Return an integer
* >0: this > other
* 0: this == other
* <0: this < other
*/
public int
cmp( FDBigInt other ){
int i;
if ( this.nWords > other.nWords ){
// if any of my high-order words is non-zero,
// then the answer is evident
int j = other.nWords-1;
for ( i = this.nWords-1; i > j ; i-- )
if ( this.data[i] != 0 ) return 1;
}else if ( this.nWords < other.nWords ){
// if any of other's high-order words is non-zero,
// then the answer is evident
int j = this.nWords-1;
for ( i = other.nWords-1; i > j ; i-- )
if ( other.data[i] != 0 ) return -1;
} else{
i = this.nWords-1;
}
for ( ; i > 0 ; i-- )
if ( this.data[i] != other.data[i] )
break;
// careful! want unsigned compare!
// use brute force here.
int a = this.data[i];
int b = other.data[i];
if ( a < 0 ){
// a is really big, unsigned
if ( b < 0 ){
return a-b; // both big, negative
} else {
return 1; // b not big, answer is obvious;
}
} else {
// a is not really big
if ( b < 0 ) {
// but b is really big
return -1;
} else {
return a - b;
}
}
}
/*
* Compute
* q = (int)( this / S )
* this = 10 * ( this mod S )
* Return q.
* This is the iteration step of digit development for output.
* We assume that S has been normalized, as above, and that
* "this" has been lshift'ed accordingly.
* Also assume, of course, that the result, q, can be expressed
* as an integer, 0 <= q < 10.
*/
public int
quoRemIteration( FDBigInt S )throws IllegalArgumentException {
// ensure that this and S have the same number of
// digits. If S is properly normalized and q < 10 then
// this must be so.
if ( nWords != S.nWords ){
throw new IllegalArgumentException("disparate values");
}
// estimate q the obvious way. We will usually be
// right. If not, then we're only off by a little and
// will re-add.
int n = nWords-1;
long q = ((long)data[n]&0xffffffffL) / (long)S.data[n];
long diff = 0L;
for ( int i = 0; i <= n ; i++ ){
diff += ((long)data[i]&0xffffffffL) - q*((long)S.data[i]&0xffffffffL);
data[i] = (int)diff;
diff >>= 32; // N.B. SIGNED shift.
}
if ( diff != 0L ) {
// damn, damn, damn. q is too big.
// add S back in until this turns +. This should
// not be very many times!
long sum = 0L;
while ( sum == 0L ){
sum = 0L;
for ( int i = 0; i <= n; i++ ){
sum += ((long)data[i]&0xffffffffL) + ((long)S.data[i]&0xffffffffL);
data[i] = (int) sum;
sum >>= 32; // Signed or unsigned, answer is 0 or 1
}
/*
* Originally the following line read
* "if ( sum !=0 && sum != -1 )"
* but that would be wrong, because of the
* treatment of the two values as entirely unsigned,
* it would be impossible for a carry-out to be interpreted
* as -1 -- it would have to be a single-bit carry-out, or
* +1.
*/
assert sum == 0 || sum == 1 : sum; // carry out of division correction
q -= 1;
}
}
// finally, we can multiply this by 10.
// it cannot overflow, right, as the high-order word has
// at least 4 high-order zeros!
long p = 0L;
for ( int i = 0; i <= n; i++ ){
p += 10*((long)data[i]&0xffffffffL);
data[i] = (int)p;
p >>= 32; // SIGNED shift.
}
assert p == 0L : p; // Carry out of *10
return (int)q;
}
public String
toString() {
StringBuffer r = new StringBuffer(30);
r.append('[');
int i = Math.min( nWords-1, data.length-1) ;
if ( nWords > data.length ){
r.append( "("+data.length+"<"+nWords+"!)" );
}
for( ; i> 0 ; i-- ){
r.append( Integer.toHexString( data[i] ) );
r.append(' ');
}
r.append( Integer.toHexString( data[0] ) );
r.append(']');
return new String( r );
}
}
© 2015 - 2024 Weber Informatics LLC | Privacy Policy