8014319: Faster division of large integers
Summary: Implement Burnickel-Ziegler division algorithm in BigInteger
Reviewed-by: bpb, martin
Contributed-by: Tim Buktu <tbuktu@hotmail.com>
diff --git a/test/java/math/BigInteger/BigIntegerTest.java b/test/java/math/BigInteger/BigIntegerTest.java
index 17c58bc..e02c3e5 100644
--- a/test/java/math/BigInteger/BigIntegerTest.java
+++ b/test/java/math/BigInteger/BigIntegerTest.java
@@ -43,12 +43,12 @@
* this test is a strong assurance that the BigInteger operations
* are working correctly.
*
- * Three arguments may be specified which give the number of
- * decimal digits you desire in the three batches of test numbers.
+ * Four arguments may be specified which give the number of
+ * decimal digits you desire in the four batches of test numbers.
*
* The tests are performed on arrays of random numbers which are
* generated by a Random class as well as special cases which
- * throw in boundary numbers such as 0, 1, maximum sized, etc.
+ * throw in boundary numbers such as 0, 1, maximum SIZEd, etc.
*
*/
public class BigIntegerTest {
@@ -63,11 +63,14 @@
//
// SCHOENHAGE_BASE_CONVERSION_THRESHOLD = 8 ints = 256 bits
//
+ // BURNIKEL_ZIEGLER_THRESHOLD = 50 ints = 1600 bits
+ //
static final int BITS_KARATSUBA = 1600;
static final int BITS_TOOM_COOK = 2400;
static final int BITS_KARATSUBA_SQUARE = 2880;
static final int BITS_TOOM_COOK_SQUARE = 4480;
static final int BITS_SCHOENHAGE_BASE = 256;
+ static final int BITS_BURNIKEL_ZIEGLER = 1600;
static final int ORDER_SMALL = 60;
static final int ORDER_MEDIUM = 100;
@@ -80,14 +83,15 @@
// #bits for testing Toom-Cook squaring
static final int ORDER_TOOM_COOK_SQUARE = 4600;
+ static final int SIZE = 1000; // numbers per batch
+
static Random rnd = new Random();
- static int size = 1000; // numbers per batch
static boolean failure = false;
public static void pow(int order) {
int failCount1 = 0;
- for (int i=0; i<size; i++) {
+ for (int i=0; i<SIZE; i++) {
// Test identity x^power == x*x*x ... *x
int power = rnd.nextInt(6) + 2;
BigInteger x = fetchNumber(order);
@@ -106,7 +110,7 @@
public static void square(int order) {
int failCount1 = 0;
- for (int i=0; i<size; i++) {
+ for (int i=0; i<SIZE; i++) {
// Test identity x^2 == x*x
BigInteger x = fetchNumber(order);
BigInteger xx = x.multiply(x);
@@ -121,7 +125,7 @@
public static void arithmetic(int order) {
int failCount = 0;
- for (int i=0; i<size; i++) {
+ for (int i=0; i<SIZE; i++) {
BigInteger x = fetchNumber(order);
while(x.compareTo(BigInteger.ZERO) != 1)
x = fetchNumber(order);
@@ -187,7 +191,7 @@
int failCount = 0;
BigInteger base = BigInteger.ONE.shiftLeft(BITS_KARATSUBA - 32 - 1);
- for (int i=0; i<size; i++) {
+ for (int i=0; i<SIZE; i++) {
BigInteger x = fetchNumber(BITS_KARATSUBA - 32 - 1);
BigInteger u = base.add(x);
int a = 1 + rnd.nextInt(31);
@@ -210,7 +214,7 @@
failCount = 0;
base = base.shiftLeft(BITS_TOOM_COOK - BITS_KARATSUBA);
- for (int i=0; i<size; i++) {
+ for (int i=0; i<SIZE; i++) {
BigInteger x = fetchNumber(BITS_TOOM_COOK - 32 - 1);
BigInteger u = base.add(x);
BigInteger u2 = u.shiftLeft(1);
@@ -241,7 +245,7 @@
int failCount = 0;
BigInteger base = BigInteger.ONE.shiftLeft(BITS_KARATSUBA_SQUARE - 32 - 1);
- for (int i=0; i<size; i++) {
+ for (int i=0; i<SIZE; i++) {
BigInteger x = fetchNumber(BITS_KARATSUBA_SQUARE - 32 - 1);
BigInteger u = base.add(x);
int a = 1 + rnd.nextInt(31);
@@ -259,7 +263,7 @@
failCount = 0;
base = base.shiftLeft(BITS_TOOM_COOK_SQUARE - BITS_KARATSUBA_SQUARE);
- for (int i=0; i<size; i++) {
+ for (int i=0; i<SIZE; i++) {
BigInteger x = fetchNumber(BITS_TOOM_COOK_SQUARE - 32 - 1);
BigInteger u = base.add(x);
int a = 1 + rnd.nextInt(31);
@@ -276,10 +280,61 @@
report("squareLarge Toom-Cook", failCount);
}
+ /**
+ * Sanity test for Burnikel-Ziegler division. The Burnikel-Ziegler division
+ * algorithm is used when each of the dividend and the divisor has at least
+ * a specified number of ints in its representation. This test is based on
+ * the observation that if {@code w = u*pow(2,a)} and {@code z = v*pow(2,b)}
+ * where {@code abs(u) > abs(v)} and {@code a > b && b > 0}, then if
+ * {@code w/z = q1*z + r1} and {@code u/v = q2*v + r2}, then
+ * {@code q1 = q2*pow(2,a-b)} and {@code r1 = r2*pow(2,b)}. The test
+ * ensures that {@code v} is just under the B-Z threshold and that {@code w}
+ * and {@code z} are both over the threshold. This implies that {@code u/v}
+ * uses the standard division algorithm and {@code w/z} uses the B-Z
+ * algorithm. The results of the two algorithms are then compared using the
+ * observation described in the foregoing and if they are not equal a
+ * failure is logged.
+ */
+ public static void divideLarge() {
+ int failCount = 0;
+
+ BigInteger base = BigInteger.ONE.shiftLeft(BITS_BURNIKEL_ZIEGLER - 33);
+ for (int i=0; i<SIZE; i++) {
+ BigInteger addend = new BigInteger(BITS_BURNIKEL_ZIEGLER - 34, rnd);
+ BigInteger v = base.add(addend);
+
+ BigInteger u = v.multiply(BigInteger.valueOf(2 + rnd.nextInt(Short.MAX_VALUE - 1)));
+
+ if(rnd.nextBoolean()) {
+ u = u.negate();
+ }
+ if(rnd.nextBoolean()) {
+ v = v.negate();
+ }
+
+ int a = 17 + rnd.nextInt(16);
+ int b = 1 + rnd.nextInt(16);
+ BigInteger w = u.multiply(BigInteger.valueOf(1L << a));
+ BigInteger z = v.multiply(BigInteger.valueOf(1L << b));
+
+ BigInteger[] divideResult = u.divideAndRemainder(v);
+ divideResult[0] = divideResult[0].multiply(BigInteger.valueOf(1L << (a - b)));
+ divideResult[1] = divideResult[1].multiply(BigInteger.valueOf(1L << b));
+ BigInteger[] bzResult = w.divideAndRemainder(z);
+
+ if (divideResult[0].compareTo(bzResult[0]) != 0 ||
+ divideResult[1].compareTo(bzResult[1]) != 0) {
+ failCount++;
+ }
+ }
+
+ report("divideLarge", failCount);
+ }
+
public static void bitCount() {
int failCount = 0;
- for (int i=0; i<size*10; i++) {
+ for (int i=0; i<SIZE*10; i++) {
int x = rnd.nextInt();
BigInteger bigX = BigInteger.valueOf((long)x);
int bit = (x < 0 ? 0 : 1);
@@ -300,7 +355,7 @@
public static void bitLength() {
int failCount = 0;
- for (int i=0; i<size*10; i++) {
+ for (int i=0; i<SIZE*10; i++) {
int x = rnd.nextInt();
BigInteger bigX = BigInteger.valueOf((long)x);
int signBit = (x < 0 ? 0x80000000 : 0);
@@ -321,7 +376,7 @@
public static void bitOps(int order) {
int failCount1 = 0, failCount2 = 0, failCount3 = 0;
- for (int i=0; i<size*5; i++) {
+ for (int i=0; i<SIZE*5; i++) {
BigInteger x = fetchNumber(order);
BigInteger y;
@@ -351,7 +406,7 @@
report("clearBit/testBit for " + order + " bits", failCount1);
report("flipBit/testBit for " + order + " bits", failCount2);
- for (int i=0; i<size*5; i++) {
+ for (int i=0; i<SIZE*5; i++) {
BigInteger x = fetchNumber(order);
// Test getLowestSetBit()
@@ -375,7 +430,7 @@
// Test identity x^y == x|y &~ x&y
int failCount = 0;
- for (int i=0; i<size; i++) {
+ for (int i=0; i<SIZE; i++) {
BigInteger x = fetchNumber(order);
BigInteger y = fetchNumber(order);
BigInteger z = x.xor(y);
@@ -387,7 +442,7 @@
// Test identity x &~ y == ~(~x | y)
failCount = 0;
- for (int i=0; i<size; i++) {
+ for (int i=0; i<SIZE; i++) {
BigInteger x = fetchNumber(order);
BigInteger y = fetchNumber(order);
BigInteger z = x.andNot(y);
@@ -440,7 +495,7 @@
public static void divideAndRemainder(int order) {
int failCount1 = 0;
- for (int i=0; i<size; i++) {
+ for (int i=0; i<SIZE; i++) {
BigInteger x = fetchNumber(order).abs();
while(x.compareTo(BigInteger.valueOf(3L)) != 1)
x = fetchNumber(order).abs();
@@ -519,7 +574,7 @@
public static void byteArrayConv(int order) {
int failCount = 0;
- for (int i=0; i<size; i++) {
+ for (int i=0; i<SIZE; i++) {
BigInteger x = fetchNumber(order);
while (x.equals(BigInteger.ZERO))
x = fetchNumber(order);
@@ -536,7 +591,7 @@
public static void modInv(int order) {
int failCount = 0, successCount = 0, nonInvCount = 0;
- for (int i=0; i<size; i++) {
+ for (int i=0; i<SIZE; i++) {
BigInteger x = fetchNumber(order);
while(x.equals(BigInteger.ZERO))
x = fetchNumber(order);
@@ -565,7 +620,7 @@
public static void modExp(int order1, int order2) {
int failCount = 0;
- for (int i=0; i<size/10; i++) {
+ for (int i=0; i<SIZE/10; i++) {
BigInteger m = fetchNumber(order1).abs();
while(m.compareTo(BigInteger.ONE) != 1)
m = fetchNumber(order1).abs();
@@ -883,8 +938,8 @@
/**
* Main to interpret arguments and run several tests.
*
- * Up to three arguments may be given to specify the size of BigIntegers
- * used for call parameters 1, 2, and 3. The size is interpreted as
+ * Up to three arguments may be given to specify the SIZE of BigIntegers
+ * used for call parameters 1, 2, and 3. The SIZE is interpreted as
* the maximum number of decimal digits that the parameters will have.
*
*/
@@ -945,6 +1000,7 @@
multiplyLarge();
squareLarge();
+ divideLarge();
if (failure)
throw new RuntimeException("Failure in BigIntegerTest.");
@@ -952,7 +1008,7 @@
/*
* Get a random or boundary-case number. This is designed to provide
- * a lot of numbers that will find failure points, such as max sized
+ * a lot of numbers that will find failure points, such as max SIZEd
* numbers, empty BigIntegers, etc.
*
* If order is less than 2, order is changed to 2.
@@ -987,13 +1043,13 @@
break;
case 4: // Random bit density
- int iterations = rnd.nextInt(order-1);
- result = BigInteger.ONE.shiftLeft(rnd.nextInt(order));
- for(int i=0; i<iterations; i++) {
- BigInteger temp = BigInteger.ONE.shiftLeft(
- rnd.nextInt(order));
- result = result.or(temp);
+ byte[] val = new byte[(order+7)/8];
+ int iterations = rnd.nextInt(order);
+ for (int i=0; i<iterations; i++) {
+ int bitIdx = rnd.nextInt(order);
+ val[bitIdx/8] |= 1 << (bitIdx%8);
}
+ result = new BigInteger(1, val);
break;
case 5: // Runs of consecutive ones and zeros
result = ZERO;