BoundedRational.java revision 82e5a2f60afd1cc330ccfa0b81b8824964e976c7
1682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm/* 2682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm * Copyright (C) 2015 The Android Open Source Project 3682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm * 4682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm * Licensed under the Apache License, Version 2.0 (the "License"); 5682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm * you may not use this file except in compliance with the License. 6682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm * You may obtain a copy of the License at 7682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm * 8682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm * http://www.apache.org/licenses/LICENSE-2.0 9682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm * 10682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm * Unless required by applicable law or agreed to in writing, software 11682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm * distributed under the License is distributed on an "AS IS" BASIS, 12682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 13682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm * See the License for the specific language governing permissions and 14682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm * limitations under the License. 15682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm */ 16682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 17682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehmpackage com.android.calculator2; 18682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 19682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm// We implement rational numbers of bounded size. 20682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm// If the length of the nuumerator plus the length of the denominator 21682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm// exceeds a maximum size, we simply return null, and rely on our caller 22682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm// do something else. 23682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm// We currently never return null for a pure integer. 24682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm// TODO: Reconsider that. With some care, large factorials might 25682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm// become much faster. 26682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm// 27682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm// We also implement a number of irrational functions. These return 28682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm// a non-null result only when the result is known to be rational. 29682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 30682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehmimport java.math.BigInteger; 31682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehmimport com.hp.creals.CR; 32682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 33682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehmpublic class BoundedRational { 34682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // TODO: Maybe eventually make this extend Number? 3550ed320c1cdd1b3624f73956a80eb2f2c2f5a01dHans Boehm private static final int MAX_SIZE = 800; // total, in bits 36682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 37682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm private final BigInteger mNum; 38682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm private final BigInteger mDen; 39682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 40682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public BoundedRational(BigInteger n, BigInteger d) { 41682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm mNum = n; 42682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm mDen = d; 43682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 44682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 45682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public BoundedRational(BigInteger n) { 46682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm mNum = n; 47682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm mDen = BigInteger.ONE; 48682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 49682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 50682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public BoundedRational(long n, long d) { 51682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm mNum = BigInteger.valueOf(n); 52682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm mDen = BigInteger.valueOf(d); 53682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 54682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 55682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public BoundedRational(long n) { 56682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm mNum = BigInteger.valueOf(n); 57682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm mDen = BigInteger.valueOf(1); 58682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 59682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 60682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Debug or log messages only, not pretty. 6175ca21c698808b61238a3aff3e0a3dfd5ba95d0eHans Boehm public String toString() { 6275ca21c698808b61238a3aff3e0a3dfd5ba95d0eHans Boehm return mNum.toString() + "/" + mDen.toString(); 6375ca21c698808b61238a3aff3e0a3dfd5ba95d0eHans Boehm } 6475ca21c698808b61238a3aff3e0a3dfd5ba95d0eHans Boehm 654a6b7cb235c305761af5d7f40e74d4704e5058c8Hans Boehm // Output to user, more expensive, less useful for debugging 66013969e98ce9e3eb4f87ec6159b06a74d07b2592Hans Boehm // Not internationalized. 674a6b7cb235c305761af5d7f40e74d4704e5058c8Hans Boehm public String toNiceString() { 689e855e82dc86b3038b3e048a56be77af114e24f4Hans Boehm BoundedRational nicer = reduce().positiveDen(); 694a6b7cb235c305761af5d7f40e74d4704e5058c8Hans Boehm String result = nicer.mNum.toString(); 704a6b7cb235c305761af5d7f40e74d4704e5058c8Hans Boehm if (!nicer.mDen.equals(BigInteger.ONE)) { 714a6b7cb235c305761af5d7f40e74d4704e5058c8Hans Boehm result += "/" + nicer.mDen; 724a6b7cb235c305761af5d7f40e74d4704e5058c8Hans Boehm } 734a6b7cb235c305761af5d7f40e74d4704e5058c8Hans Boehm return result; 744a6b7cb235c305761af5d7f40e74d4704e5058c8Hans Boehm } 754a6b7cb235c305761af5d7f40e74d4704e5058c8Hans Boehm 76682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static String toString(BoundedRational r) { 77682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r == null) return "not a small rational"; 7875ca21c698808b61238a3aff3e0a3dfd5ba95d0eHans Boehm return r.toString(); 79682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 80682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 81682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Primarily for debugging; clearly not exact 82682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public double doubleValue() { 83682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return mNum.doubleValue() / mDen.doubleValue(); 84682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 85682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 86682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public CR CRValue() { 87682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return CR.valueOf(mNum).divide(CR.valueOf(mDen)); 88682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 89682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 9082e5a2f60afd1cc330ccfa0b81b8824964e976c7Hans Boehm // Approximate number of bits to left of binary point. 9182e5a2f60afd1cc330ccfa0b81b8824964e976c7Hans Boehm public int wholeNumberBits() { 9282e5a2f60afd1cc330ccfa0b81b8824964e976c7Hans Boehm return mNum.bitLength() - mDen.bitLength(); 9382e5a2f60afd1cc330ccfa0b81b8824964e976c7Hans Boehm } 9482e5a2f60afd1cc330ccfa0b81b8824964e976c7Hans Boehm 95682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm private boolean tooBig() { 96682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (mDen.equals(BigInteger.ONE)) return false; 97682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return (mNum.bitLength() + mDen.bitLength() > MAX_SIZE); 98682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 99682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 100682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // return an equivalent fraction with a positive denominator. 1019e855e82dc86b3038b3e048a56be77af114e24f4Hans Boehm private BoundedRational positiveDen() { 102682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (mDen.compareTo(BigInteger.ZERO) > 0) return this; 103682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return new BoundedRational(mNum.negate(), mDen.negate()); 104682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 105682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 106682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Return an equivalent fraction in lowest terms. 107682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm private BoundedRational reduce() { 108682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (mDen.equals(BigInteger.ONE)) return this; // Optimization only 109682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm BigInteger divisor = mNum.gcd(mDen); 110682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return new BoundedRational(mNum.divide(divisor), mDen.divide(divisor)); 111682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 112682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 113682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Return a possibly reduced version of this that's not tooBig. 114682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Return null if none exists. 115682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm private BoundedRational maybeReduce() { 116682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (!tooBig()) return this; 1179e855e82dc86b3038b3e048a56be77af114e24f4Hans Boehm BoundedRational result = positiveDen(); 118682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (!result.tooBig()) return this; 119682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm result = result.reduce(); 120682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (!result.tooBig()) return this; 121682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return null; 122682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 123682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 124682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public int compareTo(BoundedRational r) { 125682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Compare by multiplying both sides by denominators, 126682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // invert result if denominator product was negative. 127682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return mNum.multiply(r.mDen).compareTo(r.mNum.multiply(mDen)) 128682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm * mDen.signum() * r.mDen.signum(); 129682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 130682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 131682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public int signum() { 13275ca21c698808b61238a3aff3e0a3dfd5ba95d0eHans Boehm return mNum.signum() * mDen.signum(); 133682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 134682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 135682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public boolean equals(BoundedRational r) { 136682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return compareTo(r) == 0; 137682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 138682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 139682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // We use static methods for arithmetic, so that we can 140682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // easily handle the null case. 141682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // We try to catch domain errors whenever possible, sometimes even when 142682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // one of the arguments is null, but not relevant. 143682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 144682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Returns equivalent BigInteger result if it exists, null if not. 145682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BigInteger asBigInteger(BoundedRational r) { 146682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r == null) return null; 147682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (!r.mDen.equals(BigInteger.ONE)) r = r.reduce(); 148682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (!r.mDen.equals(BigInteger.ONE)) return null; 149682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return r.mNum; 150682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 151682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational add(BoundedRational r1, BoundedRational r2) { 152682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r1 == null || r2 == null) return null; 153682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm final BigInteger den = r1.mDen.multiply(r2.mDen); 154682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm final BigInteger num = r1.mNum.multiply(r2.mDen) 155682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm .add(r2.mNum.multiply(r1.mDen)); 156682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return new BoundedRational(num,den).maybeReduce(); 157682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 158682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 159682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational negate(BoundedRational r) { 160682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r == null) return null; 161682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return new BoundedRational(r.mNum.negate(), r.mDen); 162682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 163682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 164682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm static BoundedRational subtract(BoundedRational r1, BoundedRational r2) { 165682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return add(r1, negate(r2)); 166682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 167682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 168682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm static BoundedRational multiply(BoundedRational r1, BoundedRational r2) { 169682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // It's tempting but marginally unsound to reduce 0 * null to zero. 170682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // The null could represent an infinite value, for which we 171682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // failed to throw an exception because it was too big. 172682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r1 == null || r2 == null) return null; 173682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm final BigInteger num = r1.mNum.multiply(r2.mNum); 174682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm final BigInteger den = r1.mDen.multiply(r2.mDen); 175682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return new BoundedRational(num,den).maybeReduce(); 176682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 177682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 178fbcef7005de4436682072927f83000b502928d25Hans Boehm public static class ZeroDivisionException extends ArithmeticException { 179fbcef7005de4436682072927f83000b502928d25Hans Boehm public ZeroDivisionException() { 180fbcef7005de4436682072927f83000b502928d25Hans Boehm super("Division by zero"); 181fbcef7005de4436682072927f83000b502928d25Hans Boehm } 182fbcef7005de4436682072927f83000b502928d25Hans Boehm } 183fbcef7005de4436682072927f83000b502928d25Hans Boehm 184682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm static BoundedRational inverse(BoundedRational r) { 185682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r == null) return null; 186682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r.mNum.equals(BigInteger.ZERO)) { 187fbcef7005de4436682072927f83000b502928d25Hans Boehm throw new ZeroDivisionException(); 188682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 189682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return new BoundedRational(r.mDen, r.mNum); 190682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 191682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 192682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm static BoundedRational divide(BoundedRational r1, BoundedRational r2) { 193682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return multiply(r1, inverse(r2)); 194682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 195682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 196682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm static BoundedRational sqrt(BoundedRational r) { 197682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Return non-null if numerator and denominator are small perfect 198682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // squares. 199682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r == null) return null; 2009e855e82dc86b3038b3e048a56be77af114e24f4Hans Boehm r = r.positiveDen().reduce(); 201682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r.mNum.compareTo(BigInteger.ZERO) < 0) { 202682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm throw new ArithmeticException("sqrt(negative)"); 203682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 204682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm final BigInteger num_sqrt = BigInteger.valueOf(Math.round(Math.sqrt( 205682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm r.mNum.doubleValue()))); 206682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (!num_sqrt.multiply(num_sqrt).equals(r.mNum)) return null; 207682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm final BigInteger den_sqrt = BigInteger.valueOf(Math.round(Math.sqrt( 208682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm r.mDen.doubleValue()))); 20975ca21c698808b61238a3aff3e0a3dfd5ba95d0eHans Boehm if (!den_sqrt.multiply(den_sqrt).equals(r.mDen)) return null; 210682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return new BoundedRational(num_sqrt, den_sqrt); 211682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 212682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 213682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public final static BoundedRational ZERO = new BoundedRational(0); 214682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public final static BoundedRational HALF = new BoundedRational(1,2); 215682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public final static BoundedRational MINUS_HALF = new BoundedRational(-1,2); 216682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public final static BoundedRational ONE = new BoundedRational(1); 217682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public final static BoundedRational MINUS_ONE = new BoundedRational(-1); 218682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public final static BoundedRational TWO = new BoundedRational(2); 219682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public final static BoundedRational MINUS_TWO = new BoundedRational(-2); 220682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public final static BoundedRational THIRTY = new BoundedRational(30); 221682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public final static BoundedRational MINUS_THIRTY = new BoundedRational(-30); 222682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public final static BoundedRational FORTY_FIVE = new BoundedRational(45); 223682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public final static BoundedRational MINUS_FORTY_FIVE = 224682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm new BoundedRational(-45); 225682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public final static BoundedRational NINETY = new BoundedRational(90); 226682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public final static BoundedRational MINUS_NINETY = new BoundedRational(-90); 227682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 228682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm private static BoundedRational map0to0(BoundedRational r) { 229682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r == null) return null; 230682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r.mNum.equals(BigInteger.ZERO)) { 231682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return ZERO; 232682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 233682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return null; 234682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 235682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 2364db31b490443e4454d98a5ae2bc44b87149accfeHans Boehm private static BoundedRational map0to1(BoundedRational r) { 2374db31b490443e4454d98a5ae2bc44b87149accfeHans Boehm if (r == null) return null; 2384db31b490443e4454d98a5ae2bc44b87149accfeHans Boehm if (r.mNum.equals(BigInteger.ZERO)) { 2394db31b490443e4454d98a5ae2bc44b87149accfeHans Boehm return ONE; 2404db31b490443e4454d98a5ae2bc44b87149accfeHans Boehm } 2414db31b490443e4454d98a5ae2bc44b87149accfeHans Boehm return null; 2424db31b490443e4454d98a5ae2bc44b87149accfeHans Boehm } 2434db31b490443e4454d98a5ae2bc44b87149accfeHans Boehm 244682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm private static BoundedRational map1to0(BoundedRational r) { 245682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r == null) return null; 246682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r.mNum.equals(r.mDen)) { 247682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return ZERO; 248682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 249682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return null; 250682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 251682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 252682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Throw an exception if the argument is definitely out of bounds for asin 253682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // or acos. 254682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm private static void checkAsinDomain(BoundedRational r) { 255682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r == null) return; 256682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r.mNum.abs().compareTo(r.mDen.abs()) > 0) { 257682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm throw new ArithmeticException("inverse trig argument out of range"); 258682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 259682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 260682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 261682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational sin(BoundedRational r) { 262682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return map0to0(r); 263682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 264682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 265682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm private final static BigInteger BIG360 = BigInteger.valueOf(360); 266682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 267682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational degreeSin(BoundedRational r) { 268682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm final BigInteger r_BI = asBigInteger(r); 269682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r_BI == null) return null; 270682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm final int r_int = r_BI.mod(BIG360).intValue(); 271682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r_int % 30 != 0) return null; 272682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm switch (r_int / 10) { 273682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm case 0: 274682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return ZERO; 275682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm case 3: // 30 degrees 276682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return HALF; 277682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm case 9: 278682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return ONE; 279682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm case 15: 280682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return HALF; 281682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm case 18: // 180 degrees 282682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return ZERO; 283682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm case 21: 284682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return MINUS_HALF; 285682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm case 27: 286682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return MINUS_ONE; 287682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm case 33: 288682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return MINUS_HALF; 289682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm default: 290682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return null; 291682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 292682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 293682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 294682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational asin(BoundedRational r) { 295682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm checkAsinDomain(r); 296682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return map0to0(r); 297682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 298682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 299682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational degreeAsin(BoundedRational r) { 300682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm checkAsinDomain(r); 301682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm final BigInteger r2_BI = asBigInteger(multiply(r, TWO)); 302682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r2_BI == null) return null; 303682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm final int r2_int = r2_BI.intValue(); 304682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Somewhat surprisingly, it seems to be the case that the following 305682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // covers all rational cases: 306682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm switch (r2_int) { 307682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm case -2: // Corresponding to -1 argument 308682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return MINUS_NINETY; 309682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm case -1: // Corresponding to -1/2 argument 310682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return MINUS_THIRTY; 311682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm case 0: 312682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return ZERO; 313682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm case 1: 314682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return THIRTY; 315682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm case 2: 316682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return NINETY; 317682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm default: 318682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm throw new AssertionError("Impossible asin arg"); 319682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 320682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 321682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 322682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational tan(BoundedRational r) { 323682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Unlike the degree case, we cannot check for the singularity, 324682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // since it occurs at an irrational argument. 325682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return map0to0(r); 326682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 327682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 328682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational degreeTan(BoundedRational r) { 329682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm final BoundedRational degree_sin = degreeSin(r); 330682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm final BoundedRational degree_cos = degreeCos(r); 331682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (degree_cos != null && degree_cos.mNum.equals(BigInteger.ZERO)) { 332682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm throw new ArithmeticException("Tangent undefined"); 333682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 334682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return divide(degree_sin, degree_cos); 335682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 336682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 337682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational atan(BoundedRational r) { 338682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return map0to0(r); 339682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 340682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 341682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational degreeAtan(BoundedRational r) { 342682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm final BigInteger r_BI = asBigInteger(r); 343682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r_BI == null) return null; 344682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r_BI.abs().compareTo(BigInteger.ONE) > 0) return null; 345682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm final int r_int = r_BI.intValue(); 346682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Again, these seem to be all rational cases: 347682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm switch (r_int) { 348682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm case -1: 349682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return MINUS_FORTY_FIVE; 350682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm case 0: 351682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return ZERO; 352682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm case 1: 353682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return FORTY_FIVE; 354682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm default: 355682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm throw new AssertionError("Impossible atan arg"); 356682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 357682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 358682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 359682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational cos(BoundedRational r) { 3604db31b490443e4454d98a5ae2bc44b87149accfeHans Boehm return map0to1(r); 361682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 362682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 363682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational degreeCos(BoundedRational r) { 364682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return degreeSin(add(r, NINETY)); 365682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 366682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 367682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational acos(BoundedRational r) { 368682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm checkAsinDomain(r); 369682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return map1to0(r); 370682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 371682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 372682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational degreeAcos(BoundedRational r) { 373682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm final BoundedRational asin_r = degreeAsin(r); 374682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return subtract(NINETY, asin_r); 375682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 376682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 377682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm private static final BigInteger BIG_TWO = BigInteger.valueOf(2); 378682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 379682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Compute an integral power of this 380682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm private BoundedRational pow(BigInteger exp) { 381682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (exp.compareTo(BigInteger.ZERO) < 0) { 382682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return inverse(pow(exp.negate())); 383682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 384682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (exp.equals(BigInteger.ONE)) return this; 385682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (exp.and(BigInteger.ONE).intValue() == 1) { 386682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return multiply(pow(exp.subtract(BigInteger.ONE)), this); 387682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 388682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (exp.equals(BigInteger.ZERO)) { 389682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return ONE; 390682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 391682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm BoundedRational tmp = pow(exp.shiftRight(1)); 392c1ea091ae1a4d5145069bfc6248f83f5ca8f0b31Hans Boehm if (Thread.interrupted()) { 39319e93c9806c245084d81df601d5ebc38c73dd5a7Hans Boehm throw new CR.AbortedException(); 394c1ea091ae1a4d5145069bfc6248f83f5ca8f0b31Hans Boehm } 395682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return multiply(tmp, tmp); 396682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 397682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 398682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational pow(BoundedRational base, BoundedRational exp) { 399682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (exp == null) return null; 400682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (exp.mNum.equals(BigInteger.ZERO)) { 401682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return new BoundedRational(1); 402682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 403682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (base == null) return null; 4049e855e82dc86b3038b3e048a56be77af114e24f4Hans Boehm exp = exp.reduce().positiveDen(); 405682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (!exp.mDen.equals(BigInteger.ONE)) return null; 406682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return base.pow(exp.mNum); 407682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 408682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 409682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational ln(BoundedRational r) { 4109e855e82dc86b3038b3e048a56be77af114e24f4Hans Boehm if (r != null && r.signum() <= 0) { 4119e855e82dc86b3038b3e048a56be77af114e24f4Hans Boehm throw new ArithmeticException("log(non-positive)"); 412682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 413682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return map1to0(r); 414682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 415682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 4164db31b490443e4454d98a5ae2bc44b87149accfeHans Boehm public static BoundedRational exp(BoundedRational r) { 4174db31b490443e4454d98a5ae2bc44b87149accfeHans Boehm return map0to1(r); 4184db31b490443e4454d98a5ae2bc44b87149accfeHans Boehm } 4194db31b490443e4454d98a5ae2bc44b87149accfeHans Boehm 420682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Return the base 10 log of n, if n is a power of 10, -1 otherwise. 4219e855e82dc86b3038b3e048a56be77af114e24f4Hans Boehm // n must be positive. 422682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm private static long b10Log(BigInteger n) { 423682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // This algorithm is very naive, but we doubt it matters. 424682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm long count = 0; 425682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm while (n.mod(BigInteger.TEN).equals(BigInteger.ZERO)) { 42682e5a2f60afd1cc330ccfa0b81b8824964e976c7Hans Boehm if (Thread.interrupted()) { 42782e5a2f60afd1cc330ccfa0b81b8824964e976c7Hans Boehm throw new CR.AbortedException(); 42882e5a2f60afd1cc330ccfa0b81b8824964e976c7Hans Boehm } 429682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm n = n.divide(BigInteger.TEN); 430682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm ++count; 431682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 432682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (n.equals(BigInteger.ONE)) { 433682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return count; 434682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 435682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return -1; 436682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 437682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 438682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational log(BoundedRational r) { 439682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r == null) return null; 4409e855e82dc86b3038b3e048a56be77af114e24f4Hans Boehm if (r.signum() <= 0) { 4419e855e82dc86b3038b3e048a56be77af114e24f4Hans Boehm throw new ArithmeticException("log(non-positive)"); 442682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 4439e855e82dc86b3038b3e048a56be77af114e24f4Hans Boehm r = r.reduce().positiveDen(); 444682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r == null) return null; 445682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r.mDen.equals(BigInteger.ONE)) { 446682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm long log = b10Log(r.mNum); 447682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (log != -1) return new BoundedRational(log); 448682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } else if (r.mNum.equals(BigInteger.ONE)) { 449682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm long log = b10Log(r.mDen); 450682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (log != -1) return new BoundedRational(-log); 451682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 452682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return null; 453682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 454682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 455682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Generalized factorial. 456682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Compute n * (n - step) * (n - 2 * step) * ... 457682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // This can be used to compute factorial a bit faster, especially 458682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // if BigInteger uses sub-quadratic multiplication. 459682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm private static BigInteger genFactorial(long n, long step) { 460682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (n > 4 * step) { 461682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm BigInteger prod1 = genFactorial(n, 2 * step); 462c023b734310f231244f17189788ae9c17c90b9a8Hans Boehm if (Thread.interrupted()) { 46319e93c9806c245084d81df601d5ebc38c73dd5a7Hans Boehm throw new CR.AbortedException(); 464c023b734310f231244f17189788ae9c17c90b9a8Hans Boehm } 465682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm BigInteger prod2 = genFactorial(n - step, 2 * step); 466c023b734310f231244f17189788ae9c17c90b9a8Hans Boehm if (Thread.interrupted()) { 46719e93c9806c245084d81df601d5ebc38c73dd5a7Hans Boehm throw new CR.AbortedException(); 468c023b734310f231244f17189788ae9c17c90b9a8Hans Boehm } 469682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return prod1.multiply(prod2); 470682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } else { 471682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm BigInteger res = BigInteger.valueOf(n); 472682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm for (long i = n - step; i > 1; i -= step) { 473682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm res = res.multiply(BigInteger.valueOf(i)); 474682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 475682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return res; 476682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 477682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 478682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 479682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Factorial; 480682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // always produces non-null (or exception) when called on non-null r. 481682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm public static BoundedRational fact(BoundedRational r) { 482682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r == null) return null; // Caller should probably preclude this case. 483682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm final BigInteger r_BI = asBigInteger(r); 484682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r_BI == null) { 485682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm throw new ArithmeticException("Non-integral factorial argument"); 486682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 487682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r_BI.signum() < 0) { 488682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm throw new ArithmeticException("Negative factorial argument"); 489682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 490682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r_BI.bitLength() > 30) { 491682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Will fail. LongValue() may not work. Punt now. 492682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm throw new ArithmeticException("Factorial argument too big"); 493682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 494682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return new BoundedRational(genFactorial(r_BI.longValue(), 1)); 495682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 496682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 497682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm private static final BigInteger BIG_FIVE = BigInteger.valueOf(5); 498cd740596f6f8bdacf02f735f116ffc544f922893Hans Boehm private static final BigInteger BIG_MINUS_ONE = BigInteger.valueOf(-1); 499682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm 500682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Return the number of decimal digits to the right of the 501682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // decimal point required to represent the argument exactly, 502682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // or Integer.MAX_VALUE if it's not possible. 503682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Never returns a value les than zero, even if r is 504682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // a power of ten. 505682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm static int digitsRequired(BoundedRational r) { 506682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r == null) return Integer.MAX_VALUE; 507682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm int powers_of_two = 0; // Max power of 2 that divides denominator 508682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm int powers_of_five = 0; // Max power of 5 that divides denominator 509682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Try the easy case first to speed things up. 510682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm if (r.mDen.equals(BigInteger.ONE)) return 0; 511682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm r = r.reduce(); 512682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm BigInteger den = r.mDen; 51382e5a2f60afd1cc330ccfa0b81b8824964e976c7Hans Boehm if (den.bitLength() > MAX_SIZE) { 51482e5a2f60afd1cc330ccfa0b81b8824964e976c7Hans Boehm return Integer.MAX_VALUE; 51582e5a2f60afd1cc330ccfa0b81b8824964e976c7Hans Boehm } 516682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm while (!den.testBit(0)) { 517682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm ++powers_of_two; 518682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm den = den.shiftRight(1); 519682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 520682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm while (den.mod(BIG_FIVE).equals(BigInteger.ZERO)) { 521682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm ++powers_of_five; 522682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm den = den.divide(BIG_FIVE); 523682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 524682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // If the denominator has a factor of other than 2 or 5 525682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // (the divisors of 10), the decimal expansion does not 526682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // terminate. Multiplying the fraction by any number of 527682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // powers of 10 will not cancel the demoniator. 528682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // (Recall the fraction was in lowest terms to start with.) 529682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // Otherwise the powers of 10 we need to cancel the denominator 530682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm // is the larger of powers_of_two and powers_of_five. 531cd740596f6f8bdacf02f735f116ffc544f922893Hans Boehm if (!den.equals(BigInteger.ONE) && !den.equals(BIG_MINUS_ONE)) { 532cd740596f6f8bdacf02f735f116ffc544f922893Hans Boehm return Integer.MAX_VALUE; 533cd740596f6f8bdacf02f735f116ffc544f922893Hans Boehm } 534682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm return Math.max(powers_of_two, powers_of_five); 535682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm } 536682ff5e8ad465d74b289590e5c88e0cf129ca90bHans Boehm} 537