diff options
-rw-r--r-- | src/com/android/calculator2/BoundedRational.java | 164 | ||||
-rw-r--r-- | src/com/android/calculator2/UnifiedReal.java | 138 |
2 files changed, 259 insertions, 43 deletions
diff --git a/src/com/android/calculator2/BoundedRational.java b/src/com/android/calculator2/BoundedRational.java index f3452a2..16fa581 100644 --- a/src/com/android/calculator2/BoundedRational.java +++ b/src/com/android/calculator2/BoundedRational.java @@ -19,6 +19,7 @@ package com.android.calculator2; import com.hp.creals.CR; import java.math.BigInteger; +import java.util.Objects; import java.util.Random; /** @@ -62,6 +63,60 @@ public class BoundedRational { } /** + * Produce BoundedRational equal to the given double. + */ + public static BoundedRational valueOf(double x) { + final long l = Math.round(x); + if ((double) l == x && Math.abs(l) <= 1000) { + return valueOf(l); + } + final long allBits = Double.doubleToRawLongBits(Math.abs(x)); + long mantissa = (allBits & ((1L << 52) - 1)); + final int biased_exp = (int)(allBits >>> 52); + if ((biased_exp & 0x7ff) == 0x7ff) { + throw new ArithmeticException("Infinity or NaN not convertible to BoundedRational"); + } + final long sign = x < 0.0 ? -1 : 1; + int exp = biased_exp - 1075; // 1023 + 52; we treat mantissa as integer. + if (biased_exp == 0) { + exp += 1; // Denormal exponent is 1 greater. + } else { + mantissa += (1L << 52); // Implied leading one. + } + BigInteger num = BigInteger.valueOf(sign * mantissa); + BigInteger den = BigInteger.ONE; + if (exp >= 0) { + num = num.shiftLeft(exp); + } else { + den = den.shiftLeft(-exp); + } + return new BoundedRational(num, den); + } + + /** + * Produce BoundedRational equal to the given long. + */ + public static BoundedRational valueOf(long x) { + if (x >= -2 && x <= 10) { + switch((int) x) { + case -2: + return MINUS_TWO; + case -1: + return MINUS_ONE; + case 0: + return ZERO; + case 1: + return ONE; + case 2: + return TWO; + case 10: + return TEN; + } + } + return new BoundedRational(x); + } + + /** * Convert to String reflecting raw representation. * Debug or log messages only, not pretty. */ @@ -108,12 +163,50 @@ public class BoundedRational { /** * Return a double approximation. - * The result is correctly rounded if numerator and denominator are - * exactly representable as a double. - * TODO: This should always be correctly rounded. + * The result is correctly rounded to nearest, with ties rounded away from zero. + * TODO: Should round ties to even. */ public double doubleValue() { - return mNum.doubleValue() / mDen.doubleValue(); + final int sign = signum(); + if (sign < 0) { + return -BoundedRational.negate(this).doubleValue(); + } + // We get the mantissa by dividing the numerator by denominator, after + // suitably prescaling them so that the integral part of the result contains + // enough bits. We do the prescaling to avoid any precision loss, so the division result + // is correctly truncated towards zero. + final int apprExp = mNum.bitLength() - mDen.bitLength(); + if (apprExp < -1100 || sign == 0) { + // Bail fast for clearly zero result. + return 0.0; + } + final int neededPrec = apprExp - 80; + final BigInteger dividend = neededPrec < 0 ? mNum.shiftLeft(-neededPrec) : mNum; + final BigInteger divisor = neededPrec > 0 ? mDen.shiftLeft(neededPrec) : mDen; + final BigInteger quotient = dividend.divide(divisor); + final int qLength = quotient.bitLength(); + int extraBits = qLength - 53; + int exponent = neededPrec + qLength; // Exponent assuming leading binary point. + if (exponent >= -1021) { + // Binary point is actually to right of leading bit. + --exponent; + } else { + // We're in the gradual underflow range. Drop more bits. + extraBits += (-1022 - exponent) + 1; + exponent = -1023; + } + final BigInteger bigMantissa = + quotient.add(BigInteger.ONE.shiftLeft(extraBits - 1)).shiftRight(extraBits); + if (exponent > 1024) { + return Double.POSITIVE_INFINITY; + } + if (exponent > -1023 && bigMantissa.bitLength() != 53 + || exponent <= -1023 && bigMantissa.bitLength() >= 53) { + throw new AssertionError("doubleValue internal error"); + } + final long mantissa = bigMantissa.longValue(); + final long bits = (mantissa & ((1l << 52) - 1)) | (((long) exponent + 1023) << 52); + return Double.longBitsToDouble(bits); } public CR crValue() { @@ -138,6 +231,10 @@ public class BoundedRational { } } + /** + * Is this number too big for us to continue with rational arithmetic? + * We return fals for integers on the assumption that we have no better fallback. + */ private boolean tooBig() { if (mDen.equals(BigInteger.ONE)) { return false; @@ -198,8 +295,16 @@ public class BoundedRational { return mNum.signum() * mDen.signum(); } - public boolean equals(BoundedRational r) { - return compareTo(r) == 0; + @Override + public int hashCode() { + // Note that this may be too expensive to be useful. + BoundedRational reduced = reduce().positiveDen(); + return Objects.hash(reduced.mNum, reduced.mDen); + } + + @Override + public boolean equals(Object r) { + return r != null && r instanceof BoundedRational && compareTo((BoundedRational) r) == 0; } // We use static methods for arithmetic, so that we can easily handle the null case. We try @@ -332,6 +437,7 @@ public class BoundedRational { public final static BoundedRational MINUS_NINETY = new BoundedRational(-90); private static final BigInteger BIG_TWO = BigInteger.valueOf(2); + private static final BigInteger BIG_MINUS_ONE = BigInteger.valueOf(-1); /** * Compute integral power of this, assuming this has been reduced and exp is >= 0. @@ -350,32 +456,59 @@ public class BoundedRational { if (Thread.interrupted()) { throw new CR.AbortedException(); } - return rawMultiply(tmp, tmp); + BoundedRational result = rawMultiply(tmp, tmp); + if (result == null || result.tooBig()) { + return null; + } + return result; } /** * Compute an integral power of this. */ public BoundedRational pow(BigInteger exp) { - if (exp.signum() < 0) { - return inverse(pow(exp.negate())); + int expSign = exp.signum(); + if (expSign == 0) { + // Questionable if base has undefined or zero value. + // java.lang.Math.pow() returns 1 anyway, so we do the same. + return BoundedRational.ONE; } if (exp.equals(BigInteger.ONE)) { return this; } // Reducing once at the beginning means there's no point in reducing later. - return reduce().rawPow(exp); + BoundedRational reduced = reduce().positiveDen(); + // First handle cases in which huge exponents could give compact results. + if (reduced.mDen.equals(BigInteger.ONE)) { + if (reduced.mNum.equals(BigInteger.ZERO)) { + return ZERO; + } + if (reduced.mNum.equals(BigInteger.ONE)) { + return ONE; + } + if (reduced.mNum.equals(BIG_MINUS_ONE)) { + if (exp.testBit(0)) { + return MINUS_ONE; + } else { + return ONE; + } + } + } + if (exp.bitLength() > 1000) { + // Stack overflow is likely; a useful rational result is not. + return null; + } + if (expSign < 0) { + return inverse(reduced).rawPow(exp.negate()); + } else { + return reduced.rawPow(exp); + } } public static BoundedRational pow(BoundedRational base, BoundedRational exp) { if (exp == null) { return null; } - if (exp.mNum.signum() == 0) { - // Questionable if base has undefined value. Java.lang.Math.pow() returns 1 anyway, - // so we do the same. - return new BoundedRational(1); - } if (base == null) { return null; } @@ -388,7 +521,6 @@ public class BoundedRational { private static final BigInteger BIG_FIVE = BigInteger.valueOf(5); - private static final BigInteger BIG_MINUS_ONE = BigInteger.valueOf(-1); /** * Return the number of decimal digits to the right of the decimal point required to represent diff --git a/src/com/android/calculator2/UnifiedReal.java b/src/com/android/calculator2/UnifiedReal.java index 61cac29..f85dd3e 100644 --- a/src/com/android/calculator2/UnifiedReal.java +++ b/src/com/android/calculator2/UnifiedReal.java @@ -84,6 +84,23 @@ public class UnifiedReal { this(new BoundedRational(n)); } + public static UnifiedReal valueOf(double x) { + if (x == 0.0 || x == 1.0) { + return valueOf((long) x); + } + return new UnifiedReal(BoundedRational.valueOf(x)); + } + + public static UnifiedReal valueOf(long x) { + if (x == 0) { + return UnifiedReal.ZERO; + } else if (x == 1) { + return UnifiedReal.ONE; + } else { + return new UnifiedReal(BoundedRational.valueOf(x)); + } + } + // Various helpful constants private final static BigInteger BIG_24 = BigInteger.valueOf(24); private final static int DEFAULT_COMPARE_TOLERANCE = -1000; @@ -173,8 +190,8 @@ public class UnifiedReal { } /** - * Given a constructive real cr, try to determine whether cr is the square root of - * a small integer. If so, return its square as a BoundedRational. Otherwise return null. + * Given a constructive real cr, try to determine whether cr is the logarithm of a small + * integer. If so, return exp(cr) as a BoundedRational. Otherwise return null. * We make this determination by simple table lookup, so spurious null returns are * entirely possible, or even likely. */ @@ -419,7 +436,8 @@ public class UnifiedReal { /** * Return a double approximation. - * TODO: Result is correctly rounded if known to be rational. + * Rational arguments are currently rounded to nearest, with ties away from zero. + * TODO: Improve rounding. */ public double doubleValue() { if (mCrFactor == CR_ONE) { @@ -506,11 +524,27 @@ public class UnifiedReal { /** * Returns true if values are definitely known to be equal, false in all other cases. + * This does not satisfy the contract for Object.equals(). */ public boolean definitelyEquals(UnifiedReal u) { return isComparable(u) && compareTo(u) == 0; } + @Override + public int hashCode() { + // Better useless than wrong. Probably. + return 0; + } + + @Override + public boolean equals(Object r) { + if (r == null || !(r instanceof UnifiedReal)) { + return false; + } + // This is almost certainly a programming error. Don't even try. + throw new AssertionError("Can't compare UnifiedReals for exact equality"); + } + /** * Returns true if values are definitely known not to be equal, false in all other cases. * Performs no approximate evaluation. @@ -669,7 +703,14 @@ public class UnifiedReal { return multiply(u.inverse()); } + /** + * Return the square root. + * This may fail to return a known rational value, even when the result is rational. + */ public UnifiedReal sqrt() { + if (definitelyZero()) { + return ZERO; + } if (mCrFactor == CR_ONE) { BoundedRational ratSqrt; // Check for all arguments of the form <perfect rational square> * small_int, @@ -866,15 +907,23 @@ public class UnifiedReal { private static final BigInteger BIG_TWO = BigInteger.valueOf(2); + // The (in abs value) integral exponent for which we attempt to use a recursive + // algorithm for evaluating pow(). The recursive algorithm works independent of the sign of the + // base, and can produce rational results. But it can become slow for very large exponents. + private static final BigInteger RECURSIVE_POW_LIMIT = BigInteger.valueOf(1000); + // The corresponding limit when we're using rational arithmetic. This should fail fast + // anyway, but we avoid ridiculously deep recursion. + private static final BigInteger HARD_RECURSIVE_POW_LIMIT = BigInteger.ONE.shiftLeft(1000); + /** - * Compute an integral power of a constrive real, using the standard recursive algorithm. + * Compute an integral power of a constructive real, using the standard recursive algorithm. * exp is known to be positive. */ private static CR recursivePow(CR base, BigInteger exp) { if (exp.equals(BigInteger.ONE)) { return base; } - if (exp.and(BigInteger.ONE).intValue() == 1) { + if (exp.testBit(0)) { return base.multiply(recursivePow(base, exp.subtract(BigInteger.ONE))); } CR tmp = recursivePow(base, exp.shiftRight(1)); @@ -885,28 +934,62 @@ public class UnifiedReal { } /** + * Compute an integral power of a constructive real, using the exp function when + * we safely can. Use recursivePow when we can't. exp is known to be nozero. + */ + private UnifiedReal expLnPow(BigInteger exp) { + int sign = signum(DEFAULT_COMPARE_TOLERANCE); + if (sign > 0) { + // Safe to take the log. This avoids deep recursion for huge exponents, which + // may actually make sense here. + return new UnifiedReal(crValue().ln().multiply(CR.valueOf(exp)).exp()); + } else if (sign < 0) { + CR result = crValue().negate().ln().multiply(CR.valueOf(exp)).exp(); + if (exp.testBit(0) /* odd exponent */) { + result = result.negate(); + } + return new UnifiedReal(result); + } else { + // Base of unknown sign with integer exponent. Use a recursive computation. + // (Another possible option would be to use the absolute value of the base, and then + // adjust the sign at the end. But that would have to be done in the CR + // implementation.) + if (exp.signum() < 0) { + // This may be very expensive if exp.negate() is large. + return new UnifiedReal(recursivePow(crValue(), exp.negate()).inverse()); + } else { + return new UnifiedReal(recursivePow(crValue(), exp)); + } + } + } + + + /** * Compute an integral power of this. * This recurses roughly as deeply as the number of bits in the exponent, and can, in * ridiculous cases, result in a stack overflow. */ private UnifiedReal pow(BigInteger exp) { - if (exp.signum() < 0) { - return pow(exp.negate()).inverse(); - } if (exp.equals(BigInteger.ONE)) { return this; } if (exp.signum() == 0) { - // Questionable if base has undefined value. Java.lang.Math.pow() returns 1 anyway, - // so we do the same. + // Questionable if base has undefined value or is 0. + // Java.lang.Math.pow() returns 1 anyway, so we do the same. return ONE; } - if (mCrFactor == CR_ONE) { + BigInteger absExp = exp.abs(); + if (mCrFactor == CR_ONE && absExp.compareTo(HARD_RECURSIVE_POW_LIMIT) <= 0) { final BoundedRational ratPow = mRatFactor.pow(exp); + // We count on this to fail, e.g. for very large exponents, when it would + // otherwise be too expensive. if (ratPow != null) { - return new UnifiedReal(mRatFactor.pow(exp)); + return new UnifiedReal(ratPow); } } + if (absExp.compareTo(RECURSIVE_POW_LIMIT) > 0) { + return expLnPow(exp); + } BoundedRational square = getSquare(mCrFactor); if (square != null) { final BoundedRational nRatFactor = @@ -920,19 +1003,15 @@ public class UnifiedReal { } } } - if (signum(DEFAULT_COMPARE_TOLERANCE) > 0) { - // Safe to take the log. This avoids deep recursion for huge exponents, which - // may actually make sense here. - return new UnifiedReal(crValue().ln().multiply(CR.valueOf(exp)).exp()); - } else { - // Possibly negative base with integer exponent. Use a recursive computation. - // (Another possible option would be to use the absolute value of the base, and then - // adjust the sign at the end. But that would have to be done in the CR - // implementation.) - return new UnifiedReal(recursivePow(crValue(), exp)); - } + return expLnPow(exp); } + /** + * Return this ^ expon. + * This is really only well-defined for a positive base, particularly since + * 0^x is not continuous at zero. (0^0 = 1 (as is epsilon^0), but 0^epsilon is 0. + * We nonetheless try to do reasonable things at zero, when we recognize that case. + */ public UnifiedReal pow(UnifiedReal expon) { if (mCrFactor == CR_E) { if (mRatFactor.equals(BoundedRational.ONE)) { @@ -956,6 +1035,14 @@ public class UnifiedReal { } } } + // If the exponent were known zero, we would have handled it above. + if (definitelyZero()) { + return ZERO; + } + int sign = signum(DEFAULT_COMPARE_TOLERANCE); + if (sign < 0) { + throw new ArithmeticException("Negative base for pow() with non-integer exponent"); + } return new UnifiedReal(crValue().ln().multiply(expon.crValue()).exp()); } @@ -964,7 +1051,7 @@ public class UnifiedReal { */ private static long pow16(int n) { if (n > 10) { - throw new AssertionError("Unexpexted pow16 argument"); + throw new AssertionError("Unexpected pow16 argument"); } long result = n*n; result *= result; @@ -1072,9 +1159,6 @@ public class UnifiedReal { } final BoundedRational crExp = getExp(mCrFactor); if (crExp != null) { - if (mRatFactor.signum() < 0) { - return negate().exp().inverse(); - } boolean needSqrt = false; BoundedRational ratExponent = mRatFactor; BigInteger asBI = BoundedRational.asBigInteger(ratExponent); |