java大数运算详解【其八】大数除法

目录

java大数运算详解【其一】大数加减法

java大数运算详解【其二】大数乘法

java大数运算详解【其三】大数乘法之平方算法之按位二次展开式算法

java大数运算详解【其四】大数乘法之平方算法之Karatsuba平方算法

java大数运算详解【其五】大数乘法之平方算法之ToomCook3平方算法

java大数运算详解【其六】大数乘法之单位乘法和经典乘法

java大数运算详解【其七】大数乘法之Karatsuba乘法和ToomCook3乘法

java大数运算详解【其八】大数除法

java大数运算详解【其九】大数除法之试商法（Knuth除法）核心算法

java大数运算详解【其十】大数除法之Burnikel-Ziegler除法算法

四、大数除法
/**
     * 返回值为{@code (this / val)}的大型整数。
     *
     * @param  除数
     * @return {@code this / val}
     * @throws ArithmeticException 如果{@code val}值零。
     */
    public BigInteger divide(BigInteger val) {
        if (val.mag.length < BURNIKEL_ZIEGLER_THRESHOLD ||
                mag.length - val.mag.length < BURNIKEL_ZIEGLER_OFFSET) {//BURNIKEL_ZIEGLER_THRESHOLD值80
            return divideKnuth(val);
        } else {
            return divideBurnikelZiegler(val);
        }
    }
1、Knuth除法算法D（循环带余数试商除法算法，简称试商法）
/**
     * 返回值为{@code (this / val)}的大型整数，使用一个 O(n^2) 的Knuth算法。
     *
     * @param  除数
     * @return {@code this / val}
     * @throws ArithmeticException 如果{@code val}值零。
     * @see MutableBigInteger#divideKnuth(MutableBigInteger, MutableBigInteger, boolean)
     */
    private BigInteger divideKnuth(BigInteger val) {
        MutableBigInteger q = new MutableBigInteger(),
                          a = new MutableBigInteger(this.mag),
                          b = new MutableBigInteger(val.mag);
        a.divideKnuth(b, q, false);
        return q.toBigInteger(this.signum * val.signum);
    }
    /**
     * 计算当前MutableBigInteger除以b的商，并将商放在提供的quotient对象中，然后返回余数对象。
     * 使用Knuth 4.3.1节中的算法D。对该算法的许多优化都是源自Colin Plumb C库。
     * 它的特例，除以一个整型数可以更快。
     * b的内容没有改变。
     */
    MutableBigInteger divideKnuth(MutableBigInteger b, MutableBigInteger quotient, boolean needRemainder) {
        if (b.intLen == 0)
            throw new ArithmeticException("BigInteger divide by zero");
        // 被除数为零
        if (intLen == 0) {
            quotient.intLen = quotient.offset = 0;
            return needRemainder ? new MutableBigInteger() : null;
        }
        int cmp = compare(b);
        // 被除数小于除数
        if (cmp < 0) {
            quotient.intLen = quotient.offset = 0;
            return needRemainder ? new MutableBigInteger(this) : null;
        }
        // 被除数等于除数
        if (cmp == 0) {
            quotient.value[0] = quotient.intLen = 1;
            quotient.offset = 0;
            return needRemainder ? new MutableBigInteger() : null;
        }
        quotient.clear();
        // 除以整型数的特殊情况
        if (b.intLen == 1) {
            int r = divideOneWord(b.value[b.offset], quotient);
            if(needRemainder) {
                if (r == 0)
                    return new MutableBigInteger();
                return new MutableBigInteger(r);
            } else {
                return null;
            }
        }
        // 如果超过KNUTH_POW2_*thresholds（？意义未明），则取消2的公共幂
        if (intLen >= KNUTH_POW2_THRESH_LEN) {//KNUTH_POW2_THRESH_LEN值6
            int trailingZeroBits = Math.min(getLowestSetBit(), b.getLowestSetBit());
            if (trailingZeroBits >= KNUTH_POW2_THRESH_ZEROS*32) {//KNUTH_POW2_THRESH_ZEROS值3
                MutableBigInteger a = new MutableBigInteger(this);
                b = new MutableBigInteger(b);
                a.rightShift(trailingZeroBits);
                b.rightShift(trailingZeroBits);
                MutableBigInteger r = a.divideKnuth(b, quotient);
                r.leftShift(trailingZeroBits);
                return r;
            }
        }
        return divideMagnitude(b, quotient, needRemainder);
    }
     /**
     * 这种方法用于除以一个整型数的除法。
     * 商被放入quotient中。一个整型数作除数是其特例。
     *
     * @return 余数.
     */
    int divideOneWord(int divisor, MutableBigInteger quotient) {
        long divisorLong = divisor & LONG_MASK;
        // 被除数为整型数的特殊情况
        if (intLen == 1) {
            long dividendValue = value[offset] & LONG_MASK;
            int q = (int) (dividendValue / divisorLong);
            int r = (int) (dividendValue - q * divisorLong);
            quotient.value[0] = q;
            quotient.intLen = (q == 0) ? 0 : 1;
            quotient.offset = 0;
            return r;
        }
        if (quotient.value.length < intLen)
            quotient.value = new int[intLen];
        quotient.offset = 0;
        quotient.intLen = intLen;
        // 规范化除数
        int shift = Integer.numberOfLeadingZeros(divisor);//计算除数的前导零的个数
        int rem = value[offset];
        long remLong = rem & LONG_MASK;
        // 确定商的最高位
        if (remLong < divisorLong) {
            quotient.value[0] = 0;
        } else {
            quotient.value[0] = (int)(remLong / divisorLong);
            rem = (int) (remLong - (quotient.value[0] * divisorLong));
            remLong = rem & LONG_MASK;
        }
        // 循环试商
        int xlen = intLen;
        while (--xlen > 0) {
            long dividendEstimate = (remLong << 32) |
                    (value[offset + intLen - xlen] & LONG_MASK);
            int q;
            if (dividendEstimate >= 0) {
                q = (int) (dividendEstimate / divisorLong);
                rem = (int) (dividendEstimate - q * divisorLong);
            } else {
                long tmp = divWord(dividendEstimate, divisor);
                q = (int) (tmp & LONG_MASK);
                rem = (int) (tmp >>> 32);
            }
            quotient.value[intLen - xlen] = q;
            remLong = rem & LONG_MASK;
        }
        quotient.normalize();
        // 反规范化
        if (shift > 0)//这是多余的，因为rem的无符号值小于divisor
            return rem % divisor;
        else
            return rem;
    }
    /**
     * 该方法将一个长整型除以一个整型（来估计两个多精度数字的qhat？意义未明）。
     * 当n的有符号值小于0时使用。
     * 返回长整型，其中，高32位包含余数，低32位包含商值。
     *
     * 由于Java没有无符号数，故亦无无符号除法。
     */
    static long divWord(long n, int d) {
        long dLong = d & LONG_MASK;
        long r;
        long q;
        if (dLong == 1) {
            q = (int)n;
            r = 0;
            return (r << 32) | (q & LONG_MASK);
        }
        // 计算近似商和余数
        q = (n >>> 1) / (dLong >>> 1);
        r = n - q*dLong;
        // 修正近似值
        while (r < 0) {
            r += dLong;
            q--;
        }
        while (r >= dLong) {
            r -= dLong;
            q++;
        }
        // n - q*dlong == r 且 0 <= r <dLong，循环结束
        return (r << 32) | (q & LONG_MASK);
    }
    /**
     * 确保MutableBigInteger是标准形式的，
     * 特别是确保没有前导零，如果大小为零，则intLen为零。
     */
    final void normalize() {
        if (intLen == 0) {
            offset = 0;
            return;
        }
        int index = offset;
        if (value[index] != 0)
            return;
        int indexBound = index+intLen;
        do {
            index++;
        } while(index < indexBound && value[index] == 0);
        int numZeros = index - offset;
        intLen -= numZeros;
        offset = (intLen == 0 ?  0 : offset+numZeros);
    }
    /**
     * 返回该MutableBigInteger中最低SetBit的索引，即尾随零的个数。
     * 如果该MutableBigInteger值为0，则返回-1。
     */
    private final int getLowestSetBit() {
        if (intLen == 0)
            return -1;
        int j, b;
        for (j=intLen-1; (j > 0) && (value[j+offset] == 0); j--)
            ;
        b = value[j+offset];
        if (b == 0)
            return -1;
        return ((intLen-1-j)<<5) + Integer.numberOfTrailingZeros(b);
    }
    /**
     * 该MutableBigInteger右移n位，以正常形式右移。
     */
    void rightShift(int n) {
        if (intLen == 0)
            return;
        int nInts = n >>> 5;
        int nBits = n & 0x1F;
        this.intLen -= nInts;
        if (nBits == 0)
            return;
        int bitsInHighWord = BigInteger.bitLengthForInt(value[offset]);
        if (nBits >= bitsInHighWord) {
            this.primitiveLeftShift(32 - nBits);
            this.intLen--;
        } else {
            primitiveRightShift(nBits);
        }
    }
    /**
     * 包私有方法（package private method）以返回整数的比特长度。
     */
    static int bitLengthForInt(int n) {
        return 32 - Integer.numberOfLeadingZeros(n);
    }
    /**
     * 该MutableBigInteger左移n位，其中n小于32。为了速度，假设intLen > 0 且 n > 0.
     */
    private final void primitiveLeftShift(int n) {
        int[] val = value;
        int n2 = 32 - n;
        for (int i=offset, c=val[i], m=i+intLen-1; i < m; i++) {
            int b = c;
            c = val[i+1];
            val[i] = (b << n) | (c >>> n2);
        }
        val[offset+intLen-1] <<= n;
    }
    /**
     * 该MutableBigInteger右移n位，其中n小于32。为了速度，假设intLen > 0 且 n > 0.
     */
    private final void primitiveRightShift(int n) {
        int[] val = value;
        int n2 = 32 - n;
        for (int i=offset+intLen-1, c=val[i]; i > offset; i--) {
            int b = c;
            c = val[i-1];
            val[i] = (c << n2) | (b >>> n);
        }
        val[offset] >>>= n;
    }
    /**
     * @see #divideKnuth(MutableBigInteger, MutableBigInteger, boolean)
     */
    MutableBigInteger divideKnuth(MutableBigInteger b, MutableBigInteger quotient) {
        return divideKnuth(b,quotient,true);
    }
    /**
     * 该MutableBigInteger左移n位。
     */
    void leftShift(int n) {
        /*
         * 如果这个MutableBigInteger中有足够的存储空间，那么可用空间将被使用。
         * value数组中使用的int值右侧的空间使用起来更快，
         * 因此如果可能，将从右侧占用额外的空间（即优先使用低位空间）。
         */
        if (intLen == 0)
           return;
        int nInts = n >>> 5;
        int nBits = n&0x1F;
        int bitsInHighWord = BigInteger.bitLengthForInt(value[offset]);
        // 如果不需要移动整型位，那么就这样做
        if (n <= (32-bitsInHighWord)) {//也不须使数组增长
            primitiveLeftShift(nBits);
            return;
        }
        int newLen = intLen + nInts +1;
        if (nBits <= (32-bitsInHighWord))
            newLen--;
        // 同时处理nInts整型位的移动
        if (value.length < newLen) {
            // 数组必须增长
            int[] result = new int[newLen];
            for (int i=0; i < intLen; i++)
                result[i] = value[offset+i];
            setValue(result, newLen);
        } else if (value.length - offset >= newLen) {
            // 使用右侧空间
            for(int i=0; i < newLen - intLen; i++)
                value[offset+intLen+i] = 0;
        } else {
            // 使用左侧空间
            for (int i=0; i < intLen; i++)
                value[i] = value[offset+i];
            for (int i=intLen; i < newLen; i++)
                value[i] = 0;
            offset = 0;
        }
        // 处理nBits位的移动
        intLen = newLen;
        if (nBits == 0)
            return;
        if (nBits <= (32-bitsInHighWord))
            primitiveLeftShift(nBits);
        else
            primitiveRightShift(32 -nBits);
    }

   核心算法：（见下节）