/*
* DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER.
*
* under the terms of the GNU General Public License version 2 only, as
* published by the Free Software Foundation. Oracle designates this
* particular file as subject to the "Classpath" exception as provided
* by Oracle in the LICENSE file that accompanied this code.
*
* This code is distributed in the hope that it will be useful, but WITHOUT
* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
* version 2 for more details (a copy is included in the LICENSE file that
* accompanied this code).
*
* You should have received a copy of the GNU General Public License version
* 2 along with this work; if not, write to the Free Software Foundation,
* Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
*
* Please contact Oracle, 500 Oracle Parkway, Redwood Shores, CA 94065 USA
* or visit www.oracle.com if you need additional information or have any
* questions.
*/
public class FloatingDecimal{
boolean isExceptional;
boolean isNegative;
int decExponent;
char digits[];
int nDigits;
int bigIntExp;
int bigIntNBits;
boolean mustSetRoundDir = false;
boolean fromHex = false;
{
isExceptional = e;
this.decExponent = decExponent;
this.nDigits = n;
}
/*
* Constants of the implementation
* Most are IEEE-754 related.
* (There are more really boring constants at the end.)
*/
/*
* count number of bits from high-order 1 bit to low-order 1 bit,
* inclusive.
*/
private static int
countBits( long v ){
//
// the strategy is to shift until we get a non-zero sign bit
// then shift until we have no bits left, counting the difference.
// we do byte shifting as a hack. Hope it helps.
//
if ( v == 0L ) return 0;
while ( ( v & highbyte ) == 0L ){
v <<= 8;
}
while ( v > 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 synchronized FDBigInt
big5pow( int p ){
assert p >= 0 : p; // negative power of 5
b5p = t;
}
return b5p[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;
}else{
}
}
}
//
// a common operation
//
private static FDBigInt
if ( p5 != 0 ){
} else {
}
}
if ( p2 != 0 ){
}
return v;
}
//
// another common operation
//
private static FDBigInt
if ( p2 != 0 ){
}
return v;
}
/*
* Make a floating double into a FDBigInt.
* This could also be structured as a FDBigInt
* constructor, but we'd have to build a lot of knowledge
* about floating-point representation into it, and we don't want to.
*
* AS A SIDE EFFECT, THIS METHOD WILL SET THE INSTANCE VARIABLES
* bigIntExp and bigIntNBits
*
*/
private FDBigInt
if ( binexp > 0 ){
} else {
binexp +=1;
lbits <<= 1;
binexp -= 1;
}
}
/*
* We now know where the high-order 1 bit is,
* and we know how many there are.
*/
lbits >>>= lowOrderZeros;
bigIntNBits = nbits;
}
/*
* Compute a number that is the ULP of the given value,
* for purposes of addition/subtraction. Generally easy.
* More difficult if subtracting and the argument
* is a normalized a power of 2, as the ULP changes at these points.
*/
double ulpval;
// for subtraction from normalized, powers of 2,
// use next-smaller exponent
binexp -= 1;
}
} else if ( binexp == 0 ){
} else {
}
return ulpval;
}
/*
* Round a double to a float.
* In addition to the fraction bits of the double,
* look at the class instance variable roundDir,
* which should help us avoid double-rounding error.
* roundDir was set in hardValueOf if the estimate was
* close enough, but not exact. It tells us which direction
* of rounding is preferred.
*/
float
// what we have here is special.
// don't worry, the right thing will happen.
return (float) dval;
}
}
/*
* 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
char digits[];
int ndigits;
int digitno;
int c;
//
// Discard non-significant low-order bits, while rounding,
// up to insignificant value.
int i;
insignificant /= 10L;
if ( i != 0 ){
decExponent += i;
// round up based on the low-order bits we're discarding
lvalue++;
}
}
// even easier subcase!
// can do int arithmetic rather than long!
ndigits = 10;
c = ivalue%10;
ivalue /= 10;
while ( c == 0 ){
decExponent++;
c = ivalue%10;
ivalue /= 10;
}
while ( ivalue != 0){
decExponent++;
c = ivalue%10;
ivalue /= 10;
}
} else {
// same algorithm as above (same bugs, too )
// but using long arithmetic.
ndigits = 20;
c = (int)(lvalue%10L);
lvalue /= 10L;
while ( c == 0 ){
decExponent++;
c = (int)(lvalue%10L);
lvalue /= 10L;
}
while ( lvalue != 0L ){
decExponent++;
c = (int)(lvalue%10L);
lvalue /= 10;
}
}
char result [];
}
//
// 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;
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;
return;
}
// else fall through.
}
digits[i] = (char)(q+1);
}
/*
* FIRST IMPORTANT CONSTRUCTOR: DOUBLE
*/
public FloatingDecimal( double d )
{
long fractBits;
int binExp;
int nSignificantBits;
// discover and delete sign
isNegative = true;
} else {
isNegative = false;
}
// Begin to unpack
// Discover obvious special cases of NaN and Infinity.
isExceptional = true;
if ( fractBits == 0L ){
} else {
digits = notANumber;
isNegative = false; // NaN has no sign!
}
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;
nDigits = 1;
return;
}
fractBits <<= 1;
binExp -= 1;
}
binExp += 1;
} else {
}
// call the routine that actually does all the hard work.
}
/*
* SECOND IMPORTANT CONSTRUCTOR: SINGLE
*/
public FloatingDecimal( float f )
{
int fractBits;
int binExp;
int nSignificantBits;
// discover and delete sign
isNegative = true;
fBits ^= singleSignMask;
} else {
isNegative = false;
}
// Begin to unpack
// Discover obvious special cases of NaN and Infinity.
isExceptional = true;
if ( fractBits == 0L ){
} else {
digits = notANumber;
isNegative = false; // NaN has no sign!
}
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;
nDigits = 1;
return;
}
fractBits <<= 1;
binExp -= 1;
}
binExp += 1;
} else {
}
binExp -= singleExpBias;
// call the routine that actually does all the hard work.
}
private void
{
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,
// or whether we must do real work.
// Look more closely at the number to decide if,
// with scaling by 10^nTinyBits, the result will fit in
// a long.
/*
* 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 ){
} else {
halfULP = 0L;
}
} else {
}
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.
*/
int Bbits; // binary digits needed to represent B, approx.
int tenSbits; // binary digits needed to represent 10*S, approx.
/*
* 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).
*/
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.
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.
*/
int ndigit = 0;
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.
*/
// wa-hoo! They're all ints!
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 );
assert q < 10 : q; // excessively large digit
if ( (q == 0) && ! high ){
// oops. Usually ignore leading zero.
decExp--;
} else {
}
/*
* 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
*/
}
q = b / s;
b = 10 * ( b % s );
m *= 10;
assert q < 10 : q; // excessively large digit
if ( m > 0L ){
low = (b < m );
} 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;
}
}
} else {
// still good! they're all longs!
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 );
assert q < 10 : q; // excessively large digit
if ( (q == 0) && ! high ){
// oops. Usually ignore leading zero.
decExp--;
} else {
}
/*
* 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
*/
}
q = (int) ( b / s );
b = 10 * ( b % s );
m *= 10;
assert q < 10 : q; // excessively large digit
if ( m > 0L ){
low = (b < m );
} 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;
}
}
}
} else {
int shiftBias;
/*
* We really must do FDBigInt arithmetic.
* Fist, construct our FDBigInt initial values.
*/
// normalize so that division works better
/*
* 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;
assert q < 10 : q; // excessively large digit
if ( (q == 0) && ! high ){
// oops. Usually ignore leading zero.
decExp--;
} else {
}
/*
* 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
*/
}
assert q < 10 : q; // excessively large digit
}
} else
}
/*
* 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.
} else if ( lowDigitDifference > 0 ){
roundup();
}
} else {
roundup();
}
}
}
public String
toString(){
// most brain-dead version
if ( isExceptional ){
} else {
}
}
}
int i = 0;
if (isExceptional) {
i += nDigits;
} else {
// print digits.digits.
i += charLength;
if (charLength < decExponent) {
i += charLength;
result[i++] = '.';
result[i++] = '0';
} else {
result[i++] = '.';
if (charLength < nDigits) {
int t = nDigits - charLength;
i += t;
} else {
result[i++] = '0';
}
}
result[i++] = '0';
result[i++] = '.';
if (decExponent != 0) {
i -= decExponent;
}
i += nDigits;
} else {
result[i++] = '.';
if (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) {
} else {
e %= 100;
}
}
}
return i;
}
// Per-thread buffer for string/stringbuffer conversion
protected synchronized Object initialValue() {
return new char[26];
}
};
if (buf instanceof StringBuilder)
else if (buf instanceof StringBuffer)
else
assert false;
}
public static FloatingDecimal
boolean isNegative = false;
boolean signSeen = false;
int decExp;
char c;
try{
// throws NullPointerException if null
int i = 0;
case '-':
isNegative = true;
//FALLTHROUGH
case '+':
i++;
signSeen = true;
}
// Check for NaN and Infinity strings
if(c == 'N' || c == 'I') { // possible NaN or infinity
boolean potentialNaN = false;
if(c == 'N') {
potentialNaN = true;
} else {
}
// compare Input string to "NaN" or "Infinity"
int j = 0;
while(i < l && j < targetChars.length) {
i++; j++;
}
else // something is amiss, throw exception
break parseNumber;
}
// For the candidate string to be a NaN or infinity,
// all characters in input string and target char[]
// must be matched ==> j must equal targetChars.length
// and i must equal l
: new FloatingDecimal(isNegative?
}
else { // something went wrong, throw exception
break parseNumber;
}
} else if (c == '0') { // check for hexadecimal floating-point number
if (l > i+1 ) {
return parseHexString(in);
}
} // look for and process decimal floating-point string
char[] digits = new char[ l ];
int nDigits= 0;
boolean decSeen = false;
int decPt = 0;
int nLeadZero = 0;
int nTrailZero= 0;
while ( i < l ){
case '0':
if ( nDigits > 0 ){
nTrailZero += 1;
} else {
nLeadZero += 1;
}
break; // out of switch.
case '1':
case '2':
case '3':
case '4':
case '5':
case '6':
case '7':
case '8':
case '9':
while ( nTrailZero > 0 ){
nTrailZero -= 1;
}
break; // out of switch.
case '.':
if ( decSeen ){
// already saw one ., this is the 2nd.
throw new NumberFormatException("multiple points");
}
decPt = i;
if ( signSeen ){
decPt -= 1;
}
decSeen = true;
break; // out of switch.
default:
break digitLoop;
}
i++;
}
/*
* At this point, we've scanned all the digits and decimal
* point we're going to see. Trim off leading and trailing
* zeros, which will just confuse us later, and adjust
* our initial decimal exponent accordingly.
* To review:
* we have seen i total characters.
* nLeadZero of them were zeros before any other digits.
* nTrailZero of them were zeros after any other digits.
* if ( decSeen ), then a . was seen after decPt characters
* ( including leading zeros which have been discarded )
* nDigits characters were neither lead nor trailing
* zeros, nor point
*/
/*
* special hack: if we saw no non-zero digits, then the
* answer is zero!
* Unfortunately, we feel honor-bound to keep parsing!
*/
if ( nDigits == 0 ){
nDigits = 1;
if ( nLeadZero == 0 ){
// we saw NO DIGITS AT ALL,
// not even a crummy 0!
// this is not allowed.
break parseNumber; // go throw exception
}
}
/* Our initial exponent is decPt, adjusted by the number of
* discarded zeros. Or, if there was no decPt,
* then its just nDigits adjusted by discarded trailing zeros.
*/
if ( decSeen ){
} else {
}
/*
* Look for 'e' or 'E' and an optionally signed integer.
*/
int expSign = 1;
int expVal = 0;
boolean expOverflow = false;
case '-':
expSign = -1;
//FALLTHROUGH
case '+':
i++;
}
int expAt = i;
while ( i < l ){
// the next character will cause integer
// overflow.
expOverflow = true;
}
case '0':
case '1':
case '2':
case '3':
case '4':
case '5':
case '6':
case '7':
case '8':
case '9':
continue;
default:
i--; // back up.
break expLoop; // stop parsing exponent.
}
}
//
// The intent here is to end up with
// infinity or zero, as appropriate.
// The reason for yielding such a small decExponent,
// rather than something intuitive such as
// expSign*Integer.MAX_VALUE, is that this value
// is subject to further manipulation in
// doubleValue() and floatValue(), and I don't want
// it to be able to cause overflow there!
// (The only way we can get into trouble here is for
// really outrageous nDigits+nTrailZero, such as 2 billion. )
//
} else {
// this should not overflow, since we tested
// for expVal > (MAX+N), where N >= abs(decExp)
}
// if we saw something not a digit ( or end of string )
// after the [Ee][+-], without seeing any digits at all
// this is certainly an error. If we saw some digits,
// but then some trailing garbage, that might be ok.
// so we just fall through in that case.
// HUMBUG
if ( i == expAt )
break parseNumber; // certainly bad
}
/*
* We parsed everything we could.
* If there are leftovers, then this is not good input!
*/
if ( i < l &&
((i != l - 1) ||
break parseNumber; // go throw exception
}
} catch ( StringIndexOutOfBoundsException e ){ }
}
/*
* Take a FloatingDecimal, which we presumably just scanned in,
* and find out what its value is, as a double.
*
* AS A SIDE EFFECT, SET roundDir TO INDICATE PREFERRED
* ROUNDING DIRECTION in case the result is really destined
* for a single-precision float.
*/
public strictfp double doubleValue(){
long lValue;
double dValue;
// First, check for NaN and Infinity values
if(digits == notANumber)
else
}
else {
if (mustSetRoundDir) {
roundDir = 0;
}
/*
* convert the lead kDigits to a long integer.
*/
// (special performance hack: start to do it using int)
for ( int i=1; i < iDigits; i++ ){
}
}
/*
* lValue now contains a long integer with the value of
* the first kDigits digits of the number.
* dValue contains the (double) of the same.
*/
if ( nDigits <= maxDecimalDigits ){
/*
* possibly an easy case.
* We know that the digits can be represented
* exactly. And if the exponent isn't too outrageous,
* the whole thing can be done with one operation,
* thus one rounding error.
* Note that all our constructors trim all leading and
* trailing zeros, so simple values (including zero)
* will always end up here
*/
else if ( exp >= 0 ){
if ( exp <= maxSmallTen ){
/*
* Can get the answer with one operation,
* thus one roundoff.
*/
if ( mustSetRoundDir ){
: -1;
}
}
/*
* We can multiply dValue by 10^(slop)
* and it is still "small" and exact.
* Then we can multiply by 10^(exp-slop)
* with one rounding.
*/
if ( mustSetRoundDir ){
: -1;
}
}
/*
* Else we have a hard case with a positive exp.
*/
} else {
if ( exp >= -maxSmallTen ){
/*
* Can get the answer in one division.
*/
if ( mustSetRoundDir ){
: -1;
}
}
/*
* Else we have a hard case with a negative exp.
*/
}
}
/*
* Harder cases:
* The sum of digits plus exponent is greater than
* what we think we can do with one error.
*
* Start by approximating the right answer by,
* naively, scaling by powers of 10.
*/
if ( exp > 0 ){
/*
* Lets face it. This is going to be
* Infinity. Cut to the chase.
*/
}
}
int j;
}
/*
* The reason for the weird exp > 1 condition
* in the above loop was so that the last multiply
* would get unrolled. We handle it here.
* It could overflow.
*/
if ( Double.isInfinite( t ) ){
/*
* It did overflow.
* Look more closely at the result.
* If the exponent is just one too large,
* then use the maximum finite as our estimate
* value. Else call the result infinity
* and punt it.
* ( I presume this could happen because
* rounding forces the result here to be
* an ULP or two larger than
* Double.MAX_VALUE ).
*/
t = dValue / 2.0;
t *= big10pow[j];
if ( Double.isInfinite( t ) ){
}
}
dValue = t;
}
} else if ( exp < 0 ){
/*
* Lets face it. This is going to be
* zero. Cut to the chase.
*/
}
}
int j;
}
/*
* The reason for the weird exp > 1 condition
* in the above loop was so that the last multiply
* would get unrolled. We handle it here.
* It could underflow.
*/
if ( t == 0.0 ){
/*
* It did underflow.
* Look more closely at the result.
* If the exponent is just one too small,
* then use the minimum finite as our estimate
* value. Else call the result 0.0
* and punt it.
* ( I presume this could happen because
* rounding forces the result here to be
* an ULP or two less than
* Double.MIN_VALUE ).
*/
t = dValue * 2.0;
t *= tiny10pow[j];
if ( t == 0.0 ){
}
}
dValue = t;
}
}
/*
* dValue is now approximately the result.
* The hard part is adjusting it, by comparison
* with FDBigInt arithmetic.
* Formulate the EXACT big-number result as
* bigD0 * 10^exp
*/
while(true){
/* AS A SIDE EFFECT, THIS METHOD WILL SET THE INSTANCE VARIABLES
* bigIntExp and bigIntNBits
*/
/*
* Scale bigD, bigB appropriately for
* big-integer operations.
* Naively, we multiply by powers of ten
* and powers of two. What we actually do
* is keep track of the powers of 5 and
* powers of 2 we would use, then factor out
* common divisors before doing the work.
*/
int Ulp2; // powers of 2 in halfUlp.
if ( exp >= 0 ){
} else {
}
if ( bigIntExp >= 0 ){
} else {
}
// shift bigB and bigD left by a number s. t.
// halfUlp is still an integer.
int hulpbias;
// This is going to be a denormalized number
// (if not actually zero).
// half an ULP is at 2^-(expBias+expShift+1)
} else {
}
// if there are common factors of 2, we might just as well
// factor them out, as they add nothing useful.
// do multiplications by powers of 5 and 2
//
// to recap:
// bigB is the scaled-big-int version of our floating-point
// candidate.
// bigD is the scaled-big-int version of the exact value
// as we understand it.
// halfUlp is 1/2 an ulp of bigB, except for special cases
// of exact powers of 2
//
// the plan is to compare bigB with bigD, and if the difference
// is less than halfUlp, then we're satisfied. Otherwise,
// use the ratio of difference to halfUlp to calculate a fudge
// factor to add to the floating value, then go 'round again.
//
int cmpResult;
boolean overvalue;
overvalue = true; // our candidate is too big.
// candidate is a normalized exact power of 2 and
// is too big. We will be subtracting.
// For our purposes, ulp is the ulp of the
// next smaller range.
Ulp2 -= 1;
if ( Ulp2 < 0 ){
// rats. Cannot de-scale ulp this far.
// must scale diff in other direction.
Ulp2 = 0;
}
}
} else if ( cmpResult < 0 ){
overvalue = false; // our candidate is too small.
} else {
// the candidate is exactly right!
// this happens with surprising frequency
break correctionLoop;
}
// difference is small.
// this is close enough
if (mustSetRoundDir) {
}
break correctionLoop;
} else if ( cmpResult == 0 ){
// difference is exactly half an ULP
// round to some other value maybe, then finish
// should check for bigIntNBits == 1 here??
if (mustSetRoundDir) {
}
break correctionLoop;
} else {
// difference is non-trivial.
// could scale addend by ratio of difference to
// halfUlp here, if we bothered to compute that difference.
// Most of the time ( I hope ) it is about 1 anyway.
break correctionLoop; // oops. Fell off end of range.
continue; // try again.
}
}
}
}
/*
* Take a FloatingDecimal, which we presumably just scanned in,
* and find out what its value is, as a float.
* This is distinct from doubleValue() to avoid the extremely
* unlikely case of a double rounding error, wherein the conversion
* to double has one rounding error, and the conversion of that double
* to a float has another rounding error, IN THE WRONG DIRECTION,
* ( because of the preference to a zero low-order bit ).
*/
public strictfp float floatValue(){
int iValue;
float fValue;
// First, check for NaN and Infinity values
if(digits == notANumber)
else
}
else {
/*
* convert the lead kDigits to an integer.
*/
for ( int i=1; i < kDigits; i++ ){
}
/*
* iValue now contains an integer with the value of
* the first kDigits digits of the number.
* fValue contains the (float) of the same.
*/
if ( nDigits <= singleMaxDecimalDigits ){
/*
* possibly an easy case.
* We know that the digits can be represented
* exactly. And if the exponent isn't too outrageous,
* the whole thing can be done with one operation,
* thus one rounding error.
* Note that all our constructors trim all leading and
* trailing zeros, so simple values (including zero)
* will always end up here.
*/
else if ( exp >= 0 ){
if ( exp <= singleMaxSmallTen ){
/*
* Can get the answer with one operation,
* thus one roundoff.
*/
}
/*
* We can multiply dValue by 10^(slop)
* and it is still "small" and exact.
* Then we can multiply by 10^(exp-slop)
* with one rounding.
*/
}
/*
* Else we have a hard case with a positive exp.
*/
} else {
if ( exp >= -singleMaxSmallTen ){
/*
* Can get the answer in one division.
*/
}
/*
* Else we have a hard case with a negative exp.
*/
}
/*
* In double-precision, this is an exact floating integer.
* So we can compute to double, then shorten to float
* with one round, and get the right answer.
*
* First, finish accumulating digits.
* Then convert that integer to a double, multiply
* by the appropriate power of ten, and convert to float.
*/
}
}
/*
* Harder cases:
* The sum of digits plus exponent is greater than
* what we think we can do with one error.
*
* Start by weeding out obviously out-of-range
* results, then convert to double and go to
* common hard-case code.
*/
/*
* Lets face it. This is going to be
* Infinity. Cut to the chase.
*/
/*
* Lets face it. This is going to be
* zero. Cut to the chase.
*/
}
/*
* Here, we do 'way too much work, but throwing away
* our partial results, and going and doing the whole
* thing as double, then throwing away half the bits that computes
* when we convert back to float.
*
* The alternative is to reproduce the whole multiple-precision
* algorithm for float precision, or to try to parameterize it
* for common usage. The former will take about 400 lines of code,
* and the latter I tried without success. Thus the semi-hack
* answer here.
*/
double dValue = doubleValue();
return stickyRound( dValue );
}
}
/*
* All the positive powers of 10 that can be
*/
private static final double small10pow[] = {
1.0e0,
1.0e1, 1.0e2, 1.0e3, 1.0e4, 1.0e5,
1.0e6, 1.0e7, 1.0e8, 1.0e9, 1.0e10,
1.0e11, 1.0e12, 1.0e13, 1.0e14, 1.0e15,
1.0e16, 1.0e17, 1.0e18, 1.0e19, 1.0e20,
1.0e21, 1.0e22
};
private static final float singleSmall10pow[] = {
1.0e0f,
1.0e1f, 1.0e2f, 1.0e3f, 1.0e4f, 1.0e5f,
1.0e6f, 1.0e7f, 1.0e8f, 1.0e9f, 1.0e10f
};
private static final double big10pow[] = {
1e16, 1e32, 1e64, 1e128, 1e256 };
private static final double tiny10pow[] = {
1e-16, 1e-32, 1e-64, 1e-128, 1e-256 };
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,
};
/*
* Grammar is compatible with hexadecimal floating-point constants
* described in section 6.4.4.2 of the C99 specification.
*/
if (hexFloatPattern == null) {
//1 234 56 7 8 9
"([-+])?0[xX](((\\p{XDigit}+)\\.?)|((\\p{XDigit}*)\\.(\\p{XDigit}+)))[pP]([-+])?(\\p{Digit}+)[fFdD]?"
);
}
return hexFloatPattern;
}
/*
* Convert string s to a suitable floating decimal; uses the
* double constructor and set the roundDir variable appropriately
* in case the value is later converted to a float.
*/
// Verify string is a member of the hexadecimal floating-point
// string language.
boolean validInput = m.matches();
if (!validInput) {
// Input does not match pattern
} else { // validInput
/*
* We must isolate the sign, significand, and exponent
* fields. The sign value is straightforward. Since
* floating-point numbers are stored with a normalized
* representation, the significand and exponent are
* interrelated.
*
* After extracting the sign, we normalized the
* significand as a hexadecimal value, calculating an
* exponent adjust for any shifts made during
* normalization. If the significand is zero, the
* exponent doesn't need to be examined since the output
* will be zero.
*
* Next the exponent in the input string is extracted.
* Afterwards, the significand is normalized as a *binary*
* value and the input value's normalized exponent can be
* computed. The significand bits are copied into a
* double significand; if the string has more logical bits
* than can fit in a double, the extra bits affect the
* round and sticky bits which are used to round the final
* value.
*/
// Extract significand sign
// Extract Significand magnitude
/*
* Based on the form of the significand, calculate how the
* binary exponent needs to be adjusted to create a
* normalized *hexadecimal* floating-point number; that
* is, a number where there is one nonzero hex digit to
* the left of the (hexa)decimal point. Since we are
* adjusting a binary, not hexadecimal exponent, the
* exponent is adjusted by a multiple of 4.
*
* There are a number of significand scenarios to consider;
* letters are used in indicate nonzero digits:
*
* 1. 000xxxx => x.xxx normalized
* increase exponent by (number of x's - 1)*4
*
* 2. 000xxx.yyyy => x.xxyyyy normalized
* increase exponent by (number of x's - 1)*4
*
* 3. .000yyy => y.yy normalized
* decrease exponent by (number of zeros + 1)*4
*
* 4. 000.00000yyy => y.yy normalized
* decrease exponent by (number of zeros to right of point + 1)*4
*
* If the significand is exactly zero, return a properly
* signed zero.
*/
int signifLength = 0;
int exponentAdjust = 0;
{
// left of "decimal" point
// (leading zeros stripped)
// "decimal" point; leading zeros
// must always be accounted for
/*
* The significand is made up of either
*
* 1. group 4 entirely (integer portion only)
*
* OR
*
* 2. the fractional portion from group 7 plus any
* (optional) integer portions from group 6.
*/
// Leading zeros never matter on the integer portion
}
else {
// Group 6 is the optional integer; leading zeros
// never matter on the integer portion
// fraction
// Turn "integer.fraction" into "integer"+"fraction"
// check necessary?
}
/*
* Adjust exponent as described above
*/
} else { // Cases 3 and 4
}
// If the significand is zero, the exponent doesn't
// matter; return a properly signed zero.
}
}
// Extract Exponent
/*
* Use an int to read in the exponent value; this should
* provide more than sufficient range for non-contrived
* inputs. If reading the exponent in as an int does
* overflow, examine the sign of the exponent and
* significand to determine what to do.
*/
long unsignedRawExponent;
try {
}
catch (NumberFormatException e) {
// At this point, we know the exponent is
// syntactically well-formed as a sequence of
// digits. Therefore, if an NumberFormatException
// is thrown, it must be due to overflowing int's
// range. Also, at this point, we have already
// checked for a zero significand. Thus the signs
// of the exponent and significand determine the
// final result:
//
// significand
// + -
// exponent + +infinity -infinity
// - +0.0 -0.0
}
long rawExponent =
unsignedRawExponent; // exponent magnitude
// Calculate partially adjusted exponent
// Starting copying non-zero bits into proper position in
// a long; copy explicit bit too; this will be masked
// later for normal values.
boolean round = false;
boolean sticky = false;
int bitsCopied=0;
int nextShift=0;
long significand=0L;
// First iteration is different, since we only copy
// from the leading significand bit; one more exponent
// adjust will be needed...
// IMPORTANT: make leadingDigit a long to avoid
// surprising shift semantics!
/*
* Left shift the leading digit (53 - (bit position of
* leading 1 in digit)); this sets the top bit of the
* significand to 1. The nextShift value is adjusted
* to take into account the number of bit positions of
* the leadingDigit actually used. Finally, the
* exponent is adjusted to normalize the significand
* as a binary value, not just a hex value.
*/
if (leadingDigit == 1) {
/* exponent += 0 */ }
exponent += 1;
}
exponent += 2;
}
exponent += 3;
} else {
throw new AssertionError("Result from digit conversion too large!");
}
// The preceding if-else could be replaced by a single
// code block based on the high-order bit set in
// leadingDigit. Given leadingOnePosition,
// significand |= leadingDigit << (SIGNIFICAND_WIDTH - leadingOnePosition);
// nextShift = 52 - (3 + leadingOnePosition);
// exponent += (leadingOnePosition-1);
/*
* Now the exponent variable is equal to the normalized
* binary exponent. Code below will make representation
* adjustments if the exponent is incremented after
* rounding (includes overflows to infinity) or if the
* result is subnormal.
*/
// Copy digit into significand until the significand can't
// hold another full hex digit or there are no more input
// hex digits.
int i = 0;
for(i = 1;
i++) {
nextShift-=4;
}
// After the above loop, the bulk of the string is copied.
// Now, we must copy any partial hex digits into the
// significand AND compute the round bit and start computing
// sticky bit.
if ( i < signifLength ) { // at least one hex input digit exists
// from nextShift, figure out how many bits need
// to be copied, if any
switch(nextShift) { // must be negative
case -1:
// three bits need to be copied in; can
// set round bit
break;
case -2:
// two bits need to be copied in; can
// set round and start sticky
break;
case -3:
// one bit needs to be copied in
// Now set round and start sticky, if possible
break;
case -4:
// all bits copied into significand; set
// round and start sticky
// nonzeros in three low order bits?
break;
default:
throw new AssertionError("Unexpected shift distance remainder.");
// break;
}
// Round is set; sticky might be set.
// For the sticky bit, it suffices to check the
// current digit and test for any nonzero digits in
// the remaining unprocessed input.
i++;
while(i < signifLength && !sticky) {
i++;
}
}
// else all of string was seen, round and sticky are
// correct as false.
// Check for overflow and update exponent accordingly.
// overflow to properly signed infinity
} else { // Finite return value
// The result returned in this block cannot be a
// zero or subnormal; however after the
// significand is adjusted from rounding, we could
// still overflow in infinity.
// AND exponent bits into significand; if the
// significand is incremented and overflows from
// rounding, this combination will update the
// exponent correctly, even in the case of
// Double.MAX_VALUE overflowing to infinity.
significand = (( ((long)exponent +
(long)DoubleConsts.EXP_BIAS) <<
& DoubleConsts.EXP_BIT_MASK) |
} else { // Subnormal or zero
// (exponent < DoubleConsts.MIN_EXPONENT)
// No way to round back to nonzero value
// regardless of significand if the exponent is
// less than -1075.
} else { // -1075 <= exponent <= MIN_EXPONENT -1 = -1023
/*
* Find bit position to round to; recompute
* round and sticky bits, and shift
* significand right appropriately.
*/
round = false;
// Number of bits of significand to preserve is
// exponent - abs_min_exp +1
// check:
// -1075 +1074 + 1 = 0
// -1023 +1074 + 1 = 52
int bitsDiscarded = 53 -
// What to do here:
// First, isolate the new round bit
if (bitsDiscarded > 1) {
// create mask to update sticky bits; low
// order bitsDiscarded bits should be 1
}
// Now, discard the bits
(long)DoubleConsts.EXP_BIAS) <<
& DoubleConsts.EXP_BIT_MASK) |
}
}
// The significand variable now contains the currently
// appropriate exponent bits too.
/*
* Determine if significand should be incremented;
* making this determination depends on the least
* significant bit and the round and sticky bits.
*
* Round to nearest even rounding table, adapted from
* table 4.7 in "Computer Arithmetic" by IsraelKoren.
* The digit to the left of the "decimal" point is the
* least significant bit, the digits to the right of
* the point are the round and sticky bits
*
* Number Round(x)
* x0.00 x0.
* x0.01 x0.
* x0.10 x0.
* x0.11 x1. = x0. +1
* x1.00 x1.
* x1.01 x1.
* x1.10 x1. + 1
* x1.11 x1. + 1
*/
boolean incremented = false;
incremented = true;
significand++;
}
sign));
/*
* Set roundingDir variable field of fd properly so
* that the input string can be properly rounded to a
* float value. There are two cases to consider:
*
* 1. rounding to double discards sticky bit
* information that would change the result of a float
* rounding (near halfway case between two floats)
*
* 2. rounding to double rounds up when rounding up
* would not occur when rounding to float.
*
* For former case only needs to be considered when
* the bits rounded away when casting to float are all
* zero; otherwise, float round bit is properly set
* and sticky will already be true.
*
* The lower exponent bound for the code below is the
* minimum (normalized) subnormal exponent - 1 since a
* value with that exponent can round up to the
* minimum subnormal value and the sticky bit
* information must be preserved (i.e. case 1).
*/
// Outside above exponent range, the float value
// will be zero or infinity.
/*
* If the low-order 28 bits of a rounded double
* significand are 0, the double could be a
* half-way case for a rounding to float. If the
* double value is a half-way case, the double
* significand may have to be modified to round
* the the right float value (see the stickyRound
* method). If the rounding to double has lost
* what would be float sticky bit information, the
* double significand must be incremented. If the
* double value's significand was itself
* incremented, the float value may end up too
* large so the increment should be undone.
*/
// For negative values, the sign of the
// roundDir is the same as for positive values
// since adding 1 increasing the significand's
// magnitude and subtracting 1 decreases the
// significand's magnitude. If neither round
// nor sticky is true, the double value is
// exact and no adjustment is required for a
// proper float rounding.
if (leastZero) { // prerounding lsb is 0
// If round and sticky were both true,
// and the least significant
// significand bit were 0, the rounded
// significand would not have its
// low-order bits be zero. Therefore,
// we only need to adjust the
// significand if round XOR sticky is
// true.
}
}
else { // prerounding lsb is 1
// If the prerounding lsb is 1 and the
// resulting significand has its
// low-order bits zero, the significand
// was incremented. Here, we undo the
// increment, which will ensure the
// right guard and sticky bits for the
// float rounding.
if (round)
}
}
}
}
return fd;
}
}
}
/**
* Return <code>s</code> with any leading zeros removed.
*/
}
/**
* Extract a hexadecimal digit from position <code>position</code>
* of string <code>s</code>.
*/
throw new AssertionError("Unexpected failure of digit conversion of " +
}
return value;
}
}
/*
* A really, really simple bigint package
* tailored to the needs of floating base conversion.
*/
class FDBigInt {
public FDBigInt( int v ){
nWords = 1;
data = new int[1];
data[0] = v;
}
public FDBigInt( long v ){
data = new int[2];
data[0] = (int)v;
}
}
private FDBigInt( int [] d, int n ){
data = d;
nWords = n;
}
if ( n < 2 ) n = 2;
data = new int[n]; // allocate enough space
int i = nd0;
int v;
while ( i < limit ){
int ilim = i+5;
v = (int)digit[i++]-(int)'0';
while( i <ilim ){
}
}
int factor = 1;
v = 0;
while ( i < nd ){
factor *= 10;
}
if ( factor != 1 ){
}
}
/*
* Left shift by c bits.
* Shifts this in place.
*/
public void
if ( c <= 0 ){
if ( c == 0 )
return; // silly.
else
throw new IllegalArgumentException("negative shift count");
}
int wordcount = c>>5;
int bitcount = c & 0x1f;
int t[] = data;
int s[] = data;
// reallocate.
}
if ( bitcount == 0 ){
// special hack, since an anticount of 32 won't go!
} else {
while ( src >= 1 ){
}
}
while( target >= 0 ){
t[target--] = 0;
}
data = t;
// may have constructed high-order word of 0.
// if so, trim it
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
int src;
int wordcount = 0;
int bitcount = 0;
int v = 0;
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...
*/
/*
* 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.
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 )
return bitcount;
}
/*
* Multiply a FDBigInt by an int.
* Result is a new FDBigInt.
*/
public FDBigInt
long v = iv;
int r[];
long p;
// guess adequate size of r.
p = 0L;
for( int i=0; i < nWords; i++ ) {
p += v * ((long)data[i]&0xffffffffL);
r[i] = (int)p;
p >>>= 32;
}
if ( p == 0L){
} else {
r[nWords] = (int)p;
}
}
/*
* Multiply a FDBigInt by an int and add another int.
* Result is computed in place.
* Hope it fits!
*/
public void
long v = iv;
long p;
// unroll 0th iteration, doing addition.
data[0] = (int)p;
p >>>= 32;
for( int i=1; i < nWords; i++ ) {
p += v * ((long)data[i]&0xffffffffL);
data[i] = (int)p;
p >>>= 32;
}
if ( p != 0L){
nWords++;
}
}
/*
* Multiply a FDBigInt by another FDBigInt.
* Result is a new FDBigInt.
*/
public FDBigInt
// crudely guess adequate size for r
int i;
// I think I am promised zeros...
for( i = 0; i < this.nWords; i++ ){
long p = 0L;
int 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.
if ( r[i] != 0 )
break;
return new FDBigInt( r, i+1 );
}
/*
* Add one FDBigInt to another. Return a FDBigInt
*/
public FDBigInt
int i;
int a[], b[];
int n, m;
long c = 0L;
// arrange such that a.nWords >= b.nWords;
// n = a.nWords, m = b.nWords
a = this.data;
n = this.nWords;
} else {
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 ];
s[i++] = (int)c;
return new FDBigInt( s, i );
}
return new FDBigInt( r, i );
}
/*
* Subtract one FDBigInt from another. Return a FDBigInt
* Assert that the result is positive.
*/
public FDBigInt
int r[] = new int[ this.nWords ];
int i;
int n = this.nWords;
int nzeros = 0;
long c = 0L;
for ( i = 0; i < n; i++ ){
c += (long)this.data[i] & 0xffffffffL;
if ( i < m ){
}
if ( ( r[i] = (int) c ) == 0 )
nzeros++;
else
nzeros = 0;
c >>= 32; // signed shift
}
assert c == 0L : c; // borrow out of subtract
}
while ( i < m )
return false;
return true;
}
/*
* Compare FDBigInt with another FDBigInt. Return an integer
* >0: this > other
* 0: this == other
* <0: this < other
*/
public int
int i;
// if any of my high-order words is non-zero,
// then the answer is evident
for ( i = this.nWords-1; i > j ; i-- )
// if any of other's high-order words is non-zero,
// then the answer is evident
int j = this.nWords-1;
} else{
i = this.nWords-1;
}
for ( ; i > 0 ; i-- )
break;
// careful! want unsigned compare!
// use brute force here.
int a = this.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
// ensure that this and S have the same number of
// digits. If S is properly normalized and q < 10 then
// this must be so.
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 diff = 0L;
for ( int i = 0; i <= n ; i++ ){
}
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++ ){
}
/*
* 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.
*/
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++ ){
data[i] = (int)p;
p >>= 32; // SIGNED shift.
}
assert p == 0L : p; // Carry out of *10
return (int)q;
}
public long
longValue(){
// if this can be represented as a long, return the value
if (this.nWords == 1)
}
public String
toString() {
r.append('[');
}
for( ; i> 0 ; i-- ){
r.append(' ');
}
r.append(']');
return new String( r );
}
}