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

目录

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

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

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

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

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

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

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

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

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

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

    1.2、Karatsuba平方算法（二分展开式算法）
    /**
     * 使用Karatsuba平方算法平方一个大整数。
     * 当两个数字都大于某一阈值时(实验发现)，应该使用它。
     * 它是一种递归分治算法，其渐近性优于squareToLen算法。
     */
    private BigInteger squareKaratsuba() {
        int half = (mag.length+1) / 2;
        BigInteger xl = getLower(half);
        BigInteger xh = getUpper(half);
        BigInteger xhs = xh.square();  // xhs = xh^2
        BigInteger xls = xl.square();  // xls = xl^2
        // (xh^2 << 64)  +  (((xl+xh)^2 - (xh^2 + xl^2)) << 32) + xl^2
        /**
         * this=(xh<<(32*half))+xl,
         * this*this=(xhs<<(64*half))+xls+2*(xh<<(32*half))*xl=(xhs<<(64*half))+xls+(((xh+xl)*(xh+xl)-xhs-xhl)<<(32*half)),
         * 故返回：(((xhs<<(32*half))+(((xh+xl)*(xh+xl)-(xhs+xhl)))<<(32*half))+xls.
         */
        return xhs.shiftLeft(half*32).add(xl.add(xh).square().subtract(xhs.add(xls))).shiftLeft(half*32).add(xls);
    }
    /**
     * 返回一个新的BigInteger，表示该大型整数的末n位的大型整数。
     * 这用于Karatsuba乘法和Karatsuba平方。
     */
    private BigInteger getLower(int n) {
        int len = mag.length;
        if (len <= n) {
            return abs();
        }
        int lowerInts[] = new int[n];
        System.arraycopy(mag, len-n, lowerInts, 0, n);
        return new BigInteger(trustedStripLeadingZeroInts(lowerInts), 1);
    }
    /**
     * 返回一个新的BigInteger，表示该大型整数的末n位以外的大型整数。
     * 这用于Karatsuba乘法和Karatsuba平方。
     */
    private BigInteger getUpper(int n) {
        int len = mag.length;
        if (len <= n) {
            return ZERO;
        }
        int upperLen = len - n;
        int upperInts[] = new int[upperLen];
        System.arraycopy(mag, 0, upperInts, 0, upperLen);
        return new BigInteger(trustedStripLeadingZeroInts(upperInts), 1);
    }
    /**
     * 返回值为{@code (this << n)}的大型整数。
     * 左移距离{@code n}，可能为负, 那种情况下，要做正确的右移运算。
     */
    public BigInteger shiftLeft(int n) {
        if (signum == 0)
            return ZERO;
        if (n > 0) {
            return new BigInteger(shiftLeft(mag, n), signum);
        } else if (n == 0) {
            return this;
        } else {
            // (-n)中可能的int溢出不是问题，
            // 因为shiftRightImpl认为它的参数是无符号的。
            return shiftRightImpl(-n);
        }
    }
    /**
     * 返回值为{@code (mag << n)}的mag数组。
     * 左移距离{@code n}，被认为是无符号的。
     */
    private static int[] shiftLeft(int[] mag, int n) {
        int nInts = n >>> 5;//需左移int单位（32位）的次数
        int nBits = n & 0x1f;//还需左移的位数（小于32）
        int magLen = mag.length;
        int newMag[] = null;
        if (nBits == 0) {//只需左移int单位的倍数的位数，可以用数组复制操作
            newMag = new int[magLen + nInts];
            System.arraycopy(mag, 0, newMag, 0, magLen);
        } else {
            int i = 0;
            int nBits2 = 32 - nBits;
            int highBits = mag[0] >>> nBits2;//计算最高位
            if (highBits != 0) {//最高位非零
                newMag = new int[magLen + nInts + 1];
                newMag[i++] = highBits;
            } else {//最高位值零
                newMag = new int[magLen + nInts];
            }
            int j=0;
            while (j < magLen-1)//移位操作
                newMag[i++] = mag[j++] << nBits | mag[j] >>> nBits2;
            newMag[i] = mag[j] << nBits;
        }
        return newMag;
    }
    /**
     * 返回值为{@code (this >> n)}的mag数组。
     * 左移距离{@code n}，被认为是无符号的。
     */
    private BigInteger shiftRightImpl(int n) {
        int nInts = n >>> 5;
        int nBits = n & 0x1f;
        int magLen = mag.length;
        int newMag[] = null;
        // 特殊情况:整个内容移动超过末尾
        if (nInts >= magLen)
            return (signum >= 0 ? ZERO : negConst[1]);//negConst[1]值为-1.
        if (nBits == 0) {//只需右移int单位的倍数的位数，可以用数组复制操作
            int newMagLen = magLen - nInts;
            newMag = Arrays.copyOf(mag, newMagLen);
        } else {//移位操作
            int i = 0;
            int highBits = mag[0] >>> nBits;
            if (highBits != 0) {
                newMag = new int[magLen - nInts];
                newMag[i++] = highBits;
            } else {
                newMag = new int[magLen - nInts -1];
            }
            int nBits2 = 32 - nBits;
            int j=0;
            while (j < magLen - nInts - 1)
                newMag[i++] = (mag[j++] << nBits2) | (mag[j] >>> nBits);
        }
        if (signum < 0) {//符号处理
            // 找出是否有one-bit移出了端点。也就是移出的位是否有1.
            boolean onesLost = false;
            for (int i=magLen-1, j=magLen-nInts; i >= j && !onesLost; i--)
                onesLost = (mag[i] != 0);
            if (!onesLost && nBits != 0)
                onesLost = (mag[magLen - nInts - 1] << (32 - nBits) != 0);
            if (onesLost)
                newMag = javaIncrement(newMag);//对newMag数组进行自增操作
        }
        return new BigInteger(newMag, signum);
    }
    int[] javaIncrement(int[] val) {
        int lastSum = 0;
        for (int i=val.length-1;  i >= 0 && lastSum == 0; i--)
            lastSum = (val[i] += 1);
        if (lastSum == 0) {
            val = new int[val.length+1];//没有数组复制操作？
            val[0] = 1;
        }
        return val;
    }
    Karatsuba平方算法，将mag数组分成高位数组和低位数组，再进行乘法运算。因而，又称为二分展开式算法。