高精度计算

高精度计算

有时候C语言内置数据类型难以处理非常大的整数计算，此时需要自己实现大整数计算。

数字存储方式

一般的高精度计算使用字符串作为输入，每个字符数字表示数字的一个十进制位，平常的高精度计算实现类似我们平常在草稿纸上的计算方式，也就是竖式运算。

一般竖式运算从个位开始计算，所以为了计算方便，一般将大数从低位到高位存取。另外，数字的长度可能发生变化，但我们希望同样权值位始终保持对齐（例如，希望所有的个位都在下标 [0] ，所有的十位都在下标 [1] ……），这都给了「反转存储」以充分的理由。

const int LEN = 1004;

//清空数字
void clear(short a[], int len) {
    for (int i = 0; i <= len; ++i) a[i] = 0;
}

void read(short a[]) {
    char s[LEN + 1];
    cin >> s;
    int len = strlen(s);
    clear(a, len);
    a[0] = len;//首位保存数字长度
    for (int i = 0; i < len; ++i) a[len - i] = s[i] - '0';//逆序存储
}

输出也按照存储的逆序输出。因为一般高位的零要省略，所以这里从最高位开始向下寻找第一个非零位，从此处开始输出；终止条件i > 1而不是i >= 1是因为当整个数字等于0时仍希望输出一个字符0 。

ostream &operator<<(ostream &os, const short a[]) {
    int i;
    for (i = a[0]; i > 1; --i)
        if (a[i] != 0) break;
    for (; i > 0; --i) cout << a[i];
    cout << '\n';
    return os;
}

尝试打包一下这段代码

#include<iostream>

using namespace std;
const int LEN = 1004;

class BigNum {
public:
    short a[LEN + 1];

    //清空数字
    void clear() {
        for (int i = 0; i <= LEN; ++i) a[i] = 0;
    }

    BigNum() {
        clear();
    }

    BigNum(const char *s) {
        int len = strlen(s);
        clear();
        a[0] = len;//首位保存数字长度
        for (int i = 0; i < len; ++i) a[len - i] = s[i] - '0';//逆序存储
    }

    friend ostream &operator<<(ostream &os, const BigNum &num) {
        int i;
        for (i = num.a[0]; i > 1; --i)
            if (num.a[i] != 0) break;
        for (; i > 0; --i) cout << num.a[i];
        cout << '\n';
        return os;
    }

    friend istream &operator>>(istream &is, BigNum &num) {
        char s[LEN + 1];
        cin >> s;
        int len = strlen(s);
        num.clear();
        num.a[0] = len;//首位保存数字长度
        for (int i = 0; i < len; ++i) num.a[len - i] = s[i] - '0';//逆序存储
        return is;
    }
};


int main() {
    BigNum a;
    cin >> a;
    cout << a;
    return 0;
}

1、大数加法

大数加法

从最低位开始，将两个加数对应位置上的数码相加，并判断是否达到或超过10 。如果达到，那么处理进位：将更高一位的结果上增加1，当前位的结果减少10。

    friend BigNum operator+(const BigNum &left, const BigNum &right) {
        BigNum c;
        // 高精度实现中，一般令数组的最大长度 LEN 比可能的输入大一些
        // 然后略去末尾的几次循环，这样一来可以省去不少边界情况的处理
        // 因为实际输入不会超过 1000 位，故在此循环到 LEN - 1 = 1003 已经足够
        int len = max(left.a[0], right.a[0]);
        for (int i = 1; i <= LEN; ++i) {
            // 将相应位上的数码相加
            c.a[i] += left.a[i] + right.a[i];
            if (c.a[i] >= 10) {
                // 进位
                c.a[i + 1] += 1;
                c.a[i] -= 10;
            }
            if (i > len && c.a[i] != 0) len += 1;
        }
        c.a[0] = len;
        return c;
    }

2、大数减法

大数减法相比大数加法有些特殊，同样是利用竖式运算，我们平时做竖式减法的时候也可以发现，一般情况下，都用大数减小数，然后相应的添加负号。


大数减法

简洁起见，下面仅支持两个正数相减，且必须是大数减小数

    friend BigNum operator-(const BigNum &left, const BigNum &right) {
        BigNum c;
        for (int i = 1; i <= LEN; ++i) {
            // 将相应位上的数码相加
            c.a[i] += a->a[i] - b->a[i];
            if (c.a[i] < 0) {
                // 进位
                c.a[i + 1] -= 1;
                c.a[i] += 10;
            }
        }
        for (int i = 1; i <= a->a[0]; ++i) {
            if (c.a[i] != 0) c.a[0] = i;
        }
        return c;
    }

为了让算法支持负数，应该为BigNum类添加符号标志isNegative，指示当前大数的正负号；重新修正BigNum中的加法和减法，让他们支持有符号运算。

#include<iostream>
#include<string>

using namespace std;
const int LEN = 1004;

class BigNum {
public:
    short a[LEN + 1];
    bool isNegative = false;

    //清空数字
    void clear() {
        for (int i = 0; i <= LEN; ++i) a[i] = 0;
    }

    BigNum() {
        clear();
    }

    BigNum(const char *s) {
        this->clear();
        int len = strlen(s);
        this->a[0] = len;
        int i = 0;
        this->isNegative = false;
        if (s[0] == '-') this->isNegative = true;
        if (s[0] == '-' || s[0] == '+') {
            i = 1;
            this->a[0] = len - 1;
        }
        for (; i < len; ++i) this->a[len - i] = s[i] - '0';
    }

    BigNum(const BigNum &num) {
        clear();
        this->isNegative = num.isNegative;
        for (int i = 0; i <= num.a[0]; ++i) this->a[i] = num.a[i];
    }

    BigNum &operator=(const BigNum &num) {
        if (this == &num) return *this;
        this->clear();
        this->isNegative = num.isNegative;
        for (int i = 0; i <= num.a[0]; ++i) this->a[i] = num.a[i];
        return *this;
    }

    BigNum &operator=(const char *s) {
        this->clear();
        int len = strlen(s);
        this->a[0] = len;
        int i = 0;
        this->isNegative = false;
        if (s[0] == '-') this->isNegative = true;
        if (s[0] == '-' || s[0] == '+') {
            i = 1;
            this->a[0] = len - 1;
        }
        for (; i < len; ++i) this->a[len - i] = s[i] - '0';
        return *this;
    }

    BigNum operator-() {
        BigNum new_num(*this);
        new_num.isNegative = !this->isNegative;
        return new_num;
    }

    BigNum operator+() {
        return BigNum(*this);
    }

    friend BigNum operator+(const BigNum &left, const BigNum &right) {
        BigNum c;
        //异号数字相加，也就是两数字相减，符号取决于绝对值大的数字
        if (left.isNegative != right.isNegative) {
            BigNum tmp = left;
            tmp.isNegative = right.isNegative;
            c = tmp - right;
            c.isNegative = left.a[0] > right.a[0] ? left.isNegative : right.isNegative;
            return c;
        }

        // 高精度实现中，一般令数组的最大长度 LEN 比可能的输入大一些
        // 然后略去末尾的几次循环，这样一来可以省去不少边界情况的处理
        // 因为实际输入不会超过 1000 位，故在此循环到 LEN - 1 = 1003 已经足够
        int len = max(left.a[0], right.a[0]);
        for (int i = 1; i <= LEN; ++i) {
            // 将相应位上的数码相加
            c.a[i] += left.a[i] + right.a[i];
            if (c.a[i] >= 10) {
                // 进位
                c.a[i + 1] += 1;
                c.a[i] -= 10;
            }
            if (i > len && c.a[i] != 0) len += 1;
        }
        c.a[0] = len;
        c.isNegative = left.isNegative;
        return c;
    }

    friend BigNum operator-(const BigNum &left, const BigNum &right) {
        BigNum c;
        //异号数字相减，相当于无符号的两个数相加，符号取决于被减数的符号
        if (left.isNegative != right.isNegative) {
            BigNum tmp = left;
            tmp.isNegative = right.isNegative;
            c = tmp + right;
            c.isNegative = left.isNegative;
            return c;
        }

        //同号数字相减，相当于无符号大数减小数，符号取决于绝对值大的数在最终表达式中的符号
        const BigNum *a = &left;
        const BigNum *b = &right;
        c.isNegative = left.isNegative;
        if (left.a[0] < right.a[0]) {
            a = &right;
            b = &left;
            c.isNegative = !right.isNegative;
        }
        for (int i = 1; i <= LEN; ++i) {
            // 将相应位上的数码相减
            c.a[i] += a->a[i] - b->a[i];
            if (c.a[i] < 0) {
                // 进位
                c.a[i + 1] -= 1;
                c.a[i] += 10;
            }
        }

        for (int i = a->a[0]; i > 1; --i) {
            if (c.a[i] != 0) {
                c.a[0] = i;
                break;
            }
        }
        return c;
    }

    friend ostream &operator<<(ostream &os, const BigNum &num) {
        int i;
        if (num.isNegative)
            cout << '-';
        for (i = num.a[0]; i > 1; --i)
            if (num.a[i] != 0) break;
        for (; i > 0; --i) cout << num.a[i];
        cout << '\n';
        return os;
    }

    friend istream &operator>>(istream &is, BigNum &num) {
        char s[LEN + 1];
        cin >> s;
        num = s;
        return is;
    }
};

3、高精度乘法--单精度

高精度乘法的单精度版本其实相当于将a的每一位数字乘以b，然后整合结果，如果b过大，那么此算法就有可能出错（溢出）。
重整的方式，也是从个位开始逐位向上处理进位。但是这里的进位可能非常大，甚至远大于9，因为每一位被乘上之后都可能达到9b的数量级（如果确定9b不会超过可用数值类型的范围，可以使用这种方法）。所以这里的进位不能再简单地进行-10运算，而是要通过除以10的商以及余数计算。
例如:a=1337，b=42

高精度乘法--单精度

    //注意：前面将数组定义为short，这里为了增大单精度乘法的运算范围，可以改为int、long long等。
    friend BigNum operator*(const BigNum &a, int b) {
        BigNum c;
        c.isNegative=a.isNegative != b<0;
        for (int i = 1; i <= LEN; ++i) {
            // 直接把 a 的第 i 位数码乘以乘数，加入结果
            c.a[i] += a.a[i] * b;

            if (c.a[i] >= 10) {
                // 处理进位
                // c[i] / 10 即除法的商数成为进位的增量值
                c.a[i + 1] += c.a[i] / 10;
                // 而 c[i] % 10 即除法的余数成为在当前位留下的值
                c.a[i] %= 10;
            }
        }
        for(int i=LEN;i>1;--i){
            if(c.a[i]!=0){
                c.a[0]=i;
                break;
            }
        }
        return c;
    }

4、高精度乘法--高精度

高精度乘法可以完全按照竖式乘法来计算。
回想竖式乘法的每一步，实际上是计算了若干a \times b_i \times 10^i的和。例如计算1337 \times 42，计算的就是1337 \times 2 \times 10^0 + 1337 \times 4 \times 10^1。

高精度乘法

    friend BigNum operator*(const BigNum &left, const BigNum &right) {
        BigNum c;
        //同号为正、异号为负
        c.isNegative = left.isNegative != right.isNegative;
        for (int i = 1; i <= LEN; ++i) {
            // 这里直接计算结果中的从低到高第 i 位，且一并处理了进位
            // 第 i 次循环为 c[i] 加上了所有满足 p + q = i 的 a[p] 与 b[q] 的乘积之和
            // 这样做的效果和直接进行上图的运算最后求和是一样的，只是更加简短的一种实现方式
            for (int j = 1; j <= i; ++j)
                c.a[i] += left.a[j] * right.a[i - j + 1];

            if (c.a[i] >= 10) {
                c.a[i + 1] += c.a[i] / 10;
                c.a[i] %= 10;
            }
        }
        for (int i = LEN; i > 1; --i) {
            if (c.a[i] != 0) {
                c.a[0] = i;
                break;
            }
        }
        return c;
    }

5、高精度除法

高精度除法

竖式长除法实际上可以看作一个逐次减法的过程。例如上图中商数十位的计算可以这样理解：将45减去三次12后变得小于12，不能再减，故此位为3。

为了减少冗余运算，我们提前得到被除数的长度l_a与除数的长度l_b，从下标l_a-l_b开始，从高位到低位来计算商。这和手工计算时将第一次乘法的最高位与被除数最高位对齐的做法是一样的。

参考程序实现了一个函数 greater_eq() 用于判断被除数以下标 last_dg 为最低位，是否可以再减去除数而保持非负。此后对于商的每一位，不断调用 greater_eq() ，并在成立的时候用高精度减法从余数中减去除数，也即模拟了竖式除法的过程

    // 被除数 a 以下标 last_dg 为最低位，是否可以再减去除数 b 而保持非负
    // len 是除数 b 的长度，避免反复计算
    static bool greater_eq(const BigNum &a, const BigNum &b, int last_dg, int lb) {
        // 有可能被除数剩余的部分比除数长，这个情况下最多多出 1 位，故如此判断即可
        if (a.a[last_dg + lb + 1] != 0) return true;
        // 从高位到低位，逐位比较
        for (int i = lb; i >= 1; --i) {
            if (a.a[last_dg + i] > b.a[i]) return true;
            if (a.a[last_dg + i] < b.a[i]) return false;
        }
        // 相等的情形下也是可行的
        return true;
    }

    friend BigNum operator/(const BigNum &left, const BigNum &right) {
        BigNum c;
        c.isNegative = left.isNegative != right.isNegative;

        int la = left.a[0], lb = right.a[0];

        if (lb == 0) {
            cout << "> <";
            return c;
        }  // 除数不能为零

        // c 是商
        // d 是被除数的剩余部分，算法结束后自然成为余数
        BigNum d = left;
        for (int i = la - lb; i >= 0; --i) {
            // 计算商的第 i 位
            while (BigNum::greater_eq(d, right, i, lb)) {
                // 若可以减，则减
                // 这一段是一个高精度减法
                for (int j = 1; j <= lb; ++j) {
                    d.a[i + j] -= right.a[j];
                    if (d.a[i + j] < 0) {
                        d.a[i + j + 1] -= 1;
                        d.a[i + j] += 10;
                    }
                }
                // 使商的这一位增加 1
                c.a[i + 1] += 1;
                // 返回循环开头，重新检查
            }
        }
        c.a[0] = 1;
        for (int i = LEN; i > 1; --i) {
            if (c.a[i] != 0) {
                c.a[0] = i;
                break;
            }
        }
        return c;
    }

无符号简洁模版

class UBigNum {
public:
    int a[LEN + 1];

    //清空数字
    void clear() {
        for (int i = 0; i <= LEN; ++i) a[i] = 0;
    }

    UBigNum() {
        clear();
    }

    UBigNum(const char *s) {
        this->clear();
        int len = strlen(s);
        this->a[0] = len;
        int i = 0;
        for (; i < len; ++i) this->a[len - i] = s[i] - '0';
    }

    UBigNum(const UBigNum &num) {
        clear();
        for (int i = 0; i <= num.a[0]; ++i) this->a[i] = num.a[i];
    }

    UBigNum &operator=(const UBigNum &num) {
        if (this == &num) return *this;
        this->clear();
        for (int i = 0; i <= num.a[0]; ++i) this->a[i] = num.a[i];
        return *this;
    }

    UBigNum &operator=(const char *s) {
        this->clear();
        int len = strlen(s);
        this->a[0] = len;
        int i = 0;
        for (; i < len; ++i) this->a[len - i] = s[i] - '0';
        return *this;
    }

    static int computeLens(const UBigNum &num){
        int lens = 1;
        for (int i = LEN; i > 1; --i) {
            if (num.a[i] != 0) {
                lens=i;
                break;
            }
        }
        return lens;
    }

    friend UBigNum operator+(const UBigNum &left, const UBigNum &right) {
        UBigNum c;
        // 高精度实现中，一般令数组的最大长度 LEN 比可能的输入大一些
        // 然后略去末尾的几次循环，这样一来可以省去不少边界情况的处理
        // 因为实际输入不会超过 1000 位，故在此循环到 LEN - 1 = 1003 已经足够
        int len = max(left.a[0], right.a[0]);
        for (int i = 1; i <= LEN; ++i) {
            // 将相应位上的数码相加
            c.a[i] += left.a[i] + right.a[i];
            if (c.a[i] >= 10) {
                // 进位
                c.a[i + 1] += 1;
                c.a[i] -= 10;
            }
        }
        c.a[0] = c.a[len + 1] == 0 ? len : len + 1;
        return c;
    }

    //只返回大数减小数的结果，符号需要自己判断
    friend UBigNum operator-(const UBigNum &left, const UBigNum &right) {
        UBigNum c;

        //同号数字相减，相当于无符号大数减小数，符号取决于绝对值大的数在最终表达式中的符号
        const UBigNum *a = &left;
        const UBigNum *b = &right;
        if (left.a[0] < right.a[0]) {
            a = &right;
            b = &left;
        }
        for (int i = 1; i <= LEN; ++i) {
            // 将相应位上的数码相减
            c.a[i] += a->a[i] - b->a[i];
            if (c.a[i] < 0) {
                // 进位
                c.a[i + 1] -= 1;
                c.a[i] += 10;
            }
        }
        c.a[0] = UBigNum::computeLens(c);
        return c;
    }

    friend UBigNum operator*(const UBigNum &a, int b) {
        UBigNum c;
        for (int i = 1; i <= LEN; ++i) {
            // 直接把 a 的第 i 位数码乘以乘数，加入结果
            c.a[i] += a.a[i] * b;

            if (c.a[i] >= 10) {
                // 处理进位
                // c[i] / 10 即除法的商数成为进位的增量值
                c.a[i + 1] += c.a[i] / 10;
                // 而 c[i] % 10 即除法的余数成为在当前位留下的值
                c.a[i] %= 10;
            }
        }
        c.a[0] = UBigNum::computeLens(c);
        return c;
    }

    friend UBigNum operator*(const UBigNum &left, const UBigNum &right) {
        UBigNum c;
        for (int i = 1; i <= LEN; ++i) {
            // 这里直接计算结果中的从低到高第 i 位，且一并处理了进位
            // 第 i 次循环为 c[i] 加上了所有满足 p + q = i 的 a[p] 与 b[q] 的乘积之和
            // 这样做的效果和直接进行上图的运算最后求和是一样的，只是更加简短的一种实现方式
            for (int j = 1; j <= i; ++j)
                c.a[i] += left.a[j] * right.a[i - j + 1];

            if (c.a[i] >= 10) {
                c.a[i + 1] += c.a[i] / 10;
                c.a[i] %= 10;
            }
        }
        c.a[0] = UBigNum::computeLens(c);
        return c;
    }

    // 被除数 a 以下标 last_dg 为最低位，是否可以再减去除数 b 而保持非负
    // len 是除数 b 的长度，避免反复计算
    static bool greater_eq(const UBigNum &a, const UBigNum &b, int last_dg, int lb) {
        // 有可能被除数剩余的部分比除数长，这个情况下最多多出 1 位，故如此判断即可
        if (a.a[last_dg + lb + 1] != 0) return true;
        // 从高位到低位，逐位比较
        for (int i = lb; i >= 1; --i) {
            if (a.a[last_dg + i] > b.a[i]) return true;
            if (a.a[last_dg + i] < b.a[i]) return false;
        }
        // 相等的情形下也是可行的
        return true;
    }

    friend UBigNum operator/(const UBigNum &left, const UBigNum &right) {
        UBigNum c;
        int la = left.a[0], lb = right.a[0];

        if (lb == 0) {
            cout << "> <";
            return c;
        }  // 除数不能为零

        // c 是商
        // d 是被除数的剩余部分，算法结束后自然成为余数
        UBigNum d = left;
        for (int i = la - lb; i >= 0; --i) {
            // 计算商的第 i 位
            while (UBigNum::greater_eq(d, right, i, lb)) {
                // 若可以减，则减
                // 这一段是一个高精度减法
                for (int j = 1; j <= lb; ++j) {
                    d.a[i + j] -= right.a[j];
                    if (d.a[i + j] < 0) {
                        d.a[i + j + 1] -= 1;
                        d.a[i + j] += 10;
                    }
                }
                // 使商的这一位增加 1
                c.a[i + 1] += 1;
                // 返回循环开头，重新检查
            }
        }
        c.a[0] = UBigNum::computeLens(c);
        return c;
    }


    friend ostream &operator<<(ostream &os, const UBigNum &num) {
        int i;
        for (i = num.a[0]; i > 1; --i)
            if (num.a[i] != 0) break;
        for (; i > 0; --i) cout << num.a[i];
        cout << '\n';
        return os;
    }

    friend istream &operator>>(istream &is, UBigNum &num) {
        char s[LEN + 1];
        cin >> s;
        num = s;
        return is;
    }
};

优化

实际上前面所述的方案数组中的每个数字表示的是大数中的某一位数，此方案实际上对int数组的利用效率比较低，能不能使得数组中的每个数字表示连续的 t 位有效数字呢？

网友实现的版本

#include<algorithm>
#include<unordered_map>

using namespace std;
const int power = 4; //压位
const int base = 10000;
const int MAXL = 106;

struct num //高精度模板，很好理解的
{
    int a[MAXL];

    num() { memset(a, 0, sizeof a); }

    num(const char *s, int len) {
        memset(a, 0, sizeof(a));
        a[0] = (len + power - 1) / power;//a[0]保存数字长度，此时一个int保存一个4位长度的数字
        for (int i = 0, t = 0, w; i < len; w *= 10, ++i) {
            if (i % power == 0) { w = 1, ++t; }
            a[t] += w * (s[i] - '0');
        }
    }

    void add(int k) { if (k || a[0]) a[++a[0]] = k; }

    void re() { reverse(a + 1, a + a[0] + 1); }          //STL自带的反转字符串的函数
    void print() {
        printf("%d", a[a[0]]);
        for (int i = a[0] - 1; i > 0; --i)
            printf("%0*d", power, a[i]);
        printf("\n");
    }

};

bool operator<(const num &p, const num &q) {
    if (p.a[0] < q.a[0]) return true;
    if (p.a[0] > q.a[0]) return false;
    for (int i = p.a[0]; i > 0; --i) {
        if (p.a[i] != q.a[i]) return p.a[i] < q.a[i];
    }
    return false;
}

num operator*(const num &p, const num &q) {
    num c;
    c.a[0] = p.a[0] + q.a[0] - 1;
    for (int i = 1; i <= p.a[0]; ++i)
        for (int j = 1; j <= q.a[0]; ++j) {
            c.a[i + j - 1] += p.a[i] * q.a[j];
            c.a[i + j] += c.a[i + j - 1] / base;
            c.a[i + j - 1] %= base;
        }
    if (c.a[c.a[0] + 1]) ++c.a[0];
    return c;
}

根据网友版本更改而来

const int power = 4, base = 10000;//base=10^4
class UNum {
public:
    int a[LEN + 1];

    //清空数字
    void clear() {
        for (int i = 0; i <= LEN; ++i) a[i] = 0;
    }

    UNum() {
        clear();
    }

    UNum(const char *s) {
        this->clear();
        int len = strlen(s);
        //a[0]保存数字长度，此时一个int保存一个4位长度的数字
        a[0] = (len + power - 1) / power;
        for (int t = 1; t <= a[0]; ++t) {
            for (int i = max(len - t * power, 0); i < min(len, len - t * power + power); ++i) {
                a[t] = 10 * a[t] + (s[i] - '0');
            }
        }
    }

    UNum(const UNum &num) {
        clear();
        for (int i = 0; i <= num.a[0]; ++i) this->a[i] = num.a[i];
    }

    UNum &operator=(const UNum &num) {
        if (this == &num) return *this;
        this->clear();
        for (int i = 0; i <= num.a[0]; ++i) this->a[i] = num.a[i];
        return *this;
    }

    UNum &operator=(const char *s) {
        this->clear();
        int len = strlen(s);
        //a[0]保存数字长度，此时一个int保存一个4位长度的数字
        a[0] = (len + power - 1) / power;
        for (int t = 1; t <= a[0]; ++t) {
            for (int i = max(len - t * power, 0); i < min(len, len - t * power + power); ++i) {
                a[t] = 10 * a[t] + (s[i] - '0');
            }
        }
        return *this;
    }

    static int computeLens(const UNum &num) {
        int lens = 1;
        for (int i = LEN; i > 1; --i) {
            if (num.a[i] != 0) {
                lens = i;
                break;
            }
        }
        return lens;
    }

    //STL自带的反转字符串的函数
    void re() { reverse(a + 1, a + a[0] + 1); }

    friend UNum operator+(const UNum &left, const UNum &right) {
        UNum c;
        int len = max(left.a[0], right.a[0]);
        for (int i = 1; i <= LEN; ++i) {
            // 将相应位上的数码相加
            c.a[i] += left.a[i] + right.a[i];
            if (c.a[i] >= base) {
                // 进位
                c.a[i + 1] += 1;
                c.a[i] -= base;
            }
        }
        c.a[0] = c.a[len + 1] == 0 ? len : len + 1;
        return c;
    }

    //只返回大数减小数的结果，符号需要自己判断
    friend UNum operator-(const UNum &left, const UNum &right) {
        UNum c;
        const UNum *a = &left;
        const UNum *b = &right;
        if (left.a[0] < right.a[0]) {
            a = &right;
            b = &left;
        }
        for (int i = 1; i <= LEN; ++i) {
            // 将相应位上的数码相减
            c.a[i] += a->a[i] - b->a[i];
            if (c.a[i] < 0) {
                // 进位
                c.a[i + 1] -= 1;
                c.a[i] += base;
            }
        }
        c.a[0] = UNum::computeLens(c);
        return c;
    }

    friend bool operator<(const UNum &p, const UNum &q) {
        if (p.a[0] < q.a[0]) return true;
        if (p.a[0] > q.a[0]) return false;
        for (int i = p.a[0]; i > 0; --i) {
            if (p.a[i] != q.a[i]) return p.a[i] < q.a[i];
        }
        return false;
    }

    friend UNum operator*(const UNum &p, const UNum &q) {
        UNum c;
        c.a[0] = p.a[0] + q.a[0] - 1;
        for (int i = 1; i <= p.a[0]; ++i)
            for (int j = 1; j <= q.a[0]; ++j) {
                c.a[i + j - 1] += p.a[i] * q.a[j];
                c.a[i + j] += c.a[i + j - 1] / base;
                c.a[i + j - 1] %= base;
            }
        if (c.a[c.a[0] + 1]) ++c.a[0];
        return c;
    }

    friend ostream &operator<<(ostream &os, const UNum &num) {
        int i = num.a[0];
        cout << num.a[i--];
        for (; i > 0; --i) cout << setw(power) << setfill('0') << num.a[i];
        cout << '\n';
        return os;
    }

    friend istream &operator>>(istream &is, UNum &num) {
        char s[LEN + 1];
        cin >> s;
        num = s;
        return is;
    }
};

另外高精度乘法还有更优的Karatsuba算法，时间复杂度更低

Karatsuba 乘法

记高精度数字的位数为 n ，那么高精度—高精度竖式乘法需要花费 O(n^2) 的时间。本节介绍一个时间复杂度更为优秀的算法，由前苏联（俄罗斯）数学家 Anatoly Karatsuba 提出，是一种分治算法。

考虑两个十进制大整数 x 和 y ，均包含 n 个数码（可以有前导零）。任取 0 < m < n ，记

\begin{aligned} x &= x_1 \cdot 10^m + x_0, \\ y &= y_1 \cdot 10^m + y_0, \\ x \cdot y &= z_2 \cdot 10^{2m} + z_1 \cdot 10^m + z_0, \end{aligned}

其中 x_0, y_0, z_0, z_1 < 10^m 。可得

\begin{aligned} z_2 &= x_1 \cdot y_1, \\ z_1 &= x_1 \cdot y_0 + x_0 \cdot y_1, \\ z_0 &= x_0 \cdot y_0. \end{aligned}

观察知

z_1 = (x_1 + x_0) \cdot (y_1 + y_0) - z_2 - z_0,

于是要计算 z_1 ，只需计算 (x_1 + x_0) \cdot (y_1 + y_0) ，再与 z_0 、 z_2 相减即可。

上式实际上是 Karatsuba 算法的核心，它将长度为 n 的乘法问题转化为了 3 个长度更小的子问题。若令 m = \left\lceil \frac n 2 \right\rceil ，记 Karatsuba 算法计算两个 n 位整数乘法的耗时为 T(n) ，则有 T(n) = 3 \cdot T \left(\left\lceil \frac n 2 \right\rceil\right) + O(n) ，由主定理可得 T(n) = \Theta(n^{\log_2 3}) \approx \Theta(n^{1.585}) 。

整个过程可以递归实现。为清晰起见，下面的代码通过 Karatsuba 算法实现了多项式乘法，最后再处理所有的进位问题。

??? "karatsuba_mulc.cpp"

int *karatsuba_polymul(int n, int *a, int *b) {
  if (n <= 32) {
    // 规模较小时直接计算，避免继续递归带来的效率损失
    int *r = new int[n * 2 + 1]();
    for (int i = 0; i <= n; ++i)
      for (int j = 0; j <= n; ++j) r[i + j] += a[i] * b[j];
    return r;
  }
    
  int m = n / 2 + 1;
  int *r = new int[m * 4 + 1]();
  int *z0, *z1, *z2;
    
  z0 = karatsuba_polymul(m - 1, a, b);
  z2 = karatsuba_polymul(n - m, a + m, b + m);
    
  // 计算 z1
  // 临时更改，计算完毕后恢复
  for (int i = 0; i + m <= n; ++i) a[i] += a[i + m];
  for (int i = 0; i + m <= n; ++i) b[i] += b[i + m];
  z1 = karatsuba_polymul(m - 1, a, b);
  for (int i = 0; i + m <= n; ++i) a[i] -= a[i + m];
  for (int i = 0; i + m <= n; ++i) b[i] -= b[i + m];
  for (int i = 0; i <= (m - 1) * 2; ++i) z1[i] -= z0[i];
  for (int i = 0; i <= (n - m) * 2; ++i) z1[i] -= z2[i];
    
  // 由 z0、z1、z2 组合获得结果
  for (int i = 0; i <= (m - 1) * 2; ++i) r[i] += z0[i];
  for (int i = 0; i <= (m - 1) * 2; ++i) r[i + m] += z1[i];
  for (int i = 0; i <= (n - m) * 2; ++i) r[i + m * 2] += z2[i];
    
  delete[] z0;
  delete[] z1;
  delete[] z2;
  return r;
}
    
void karatsuba_mul(int a[], int b[], int c[]) {
  int *r = karatsuba_polymul(LEN - 1, a, b);
  memcpy(c, r, sizeof(int) * LEN);
  for (int i = 0; i < LEN - 1; ++i)
    if (c[i] >= 10) {
      c[i + 1] += c[i] / 10;
      c[i] %= 10;
    }
  delete[] r;
}

但是这样的实现存在一个问题：在 b 进制下，多项式的每一个系数都有可能达到 n \cdot b^2 量级，在压位高精度实现（即 b > 10 ，下文介绍）中可能造成整数溢出；而若在多项式乘法的过程中处理进位问题，则 x_1 + x_0 与 y_1 + y_0 的结果可能达到 2 \cdot b^m ，增加一个位（如果采用 x_1 - x_0 的计算方式，则不得不特殊处理负数的情况）。因此，需要依照实际的应用场景来决定采用何种实现方式。

Reference

https://en.wikipedia.org/wiki/Karatsuba_algorithm
https://oi-wiki.org/math/bignum/

封装类

这里 有一个封装好的高精度整数类。

#define MAXN 9999
// MAXN 是一位中最大的数字
#define MAXSIZE 10024
// MAXSIZE 是位数
#define DLEN 4
// DLEN 记录压几位
struct Big {
  int a[MAXSIZE], len;
  bool flag;  //标记符号'-'
  Big() {
    len = 1;
    memset(a, 0, sizeof a);
    flag = 0;
  }
  Big(const int);
  Big(const char*);
  Big(const Big&);
  Big& operator=(const Big&);  //注意这里operator有&，因为赋值有修改……
  //由于OI中要求效率
  //此处不使用泛型函数
  //故不重载
  // istream& operator>>(istream&,  BigNum&);   //重载输入运算符
  // ostream& operator<<(ostream&,  BigNum&);   //重载输出运算符
  Big operator+(const Big&) const;
  Big operator-(const Big&) const;
  Big operator*(const Big&)const;
  Big operator/(const int&) const;
  // TODO: Big / Big;
  Big operator^(const int&) const;
  // TODO: Big ^ Big;
    
  // TODO: Big 位运算;
    
  int operator%(const int&) const;
  // TODO: Big ^ Big;
  bool operator<(const Big&) const;
  bool operator<(const int& t) const;
  inline void print();
};
// README::不要随随便便把参数都变成引用，那样没办法传值
Big::Big(const int b) {
  int c, d = b;
  len = 0;
  // memset(a,0,sizeof a);
  CLR(a);
  while (d > MAXN) {
    c = d - (d / (MAXN + 1) * (MAXN + 1));
    d = d / (MAXN + 1);
    a[len++] = c;
  }
  a[len++] = d;
}
Big::Big(const char* s) {
  int t, k, index, l;
  CLR(a);
  l = strlen(s);
  len = l / DLEN;
  if (l % DLEN) ++len;
  index = 0;
  for (int i = l - 1; i >= 0; i -= DLEN) {
    t = 0;
    k = i - DLEN + 1;
    if (k < 0) k = 0;
    g(j, k, i) t = t * 10 + s[j] - '0';
    a[index++] = t;
  }
}
Big::Big(const Big& T) : len(T.len) {
  CLR(a);
  f(i, 0, len) a[i] = T.a[i];
  // TODO:重载此处？
}
Big& Big::operator=(const Big& T) {
  CLR(a);
  len = T.len;
  f(i, 0, len) a[i] = T.a[i];
  return *this;
}
Big Big::operator+(const Big& T) const {
  Big t(*this);
  int big = len;
  if (T.len > len) big = T.len;
  f(i, 0, big) {
    t.a[i] += T.a[i];
    if (t.a[i] > MAXN) {
      ++t.a[i + 1];
      t.a[i] -= MAXN + 1;
    }
  }
  if (t.a[big])
    t.len = big + 1;
  else
    t.len = big;
  return t;
}
Big Big::operator-(const Big& T) const {
  int big;
  bool ctf;
  Big t1, t2;
  if (*this < T) {
    t1 = T;
    t2 = *this;
    ctf = 1;
  } else {
    t1 = *this;
    t2 = T;
    ctf = 0;
  }
  big = t1.len;
  int j = 0;
  f(i, 0, big) {
    if (t1.a[i] < t2.a[i]) {
      j = i + 1;
      while (t1.a[j] == 0) ++j;
      --t1.a[j--];
      // WTF?
      while (j > i) t1.a[j--] += MAXN;
      t1.a[i] += MAXN + 1 - t2.a[i];
    } else
      t1.a[i] -= t2.a[i];
  }
  t1.len = big;
  while (t1.len > 1 && t1.a[t1.len - 1] == 0) {
    --t1.len;
    --big;
  }
  if (ctf) t1.a[big - 1] = -t1.a[big - 1];
  return t1;
}
Big Big::operator*(const Big& T) const {
  Big res;
  int up;
  int te, tee;
  f(i, 0, len) {
    up = 0;
    f(j, 0, T.len) {
      te = a[i] * T.a[j] + res.a[i + j] + up;
      if (te > MAXN) {
        tee = te - te / (MAXN + 1) * (MAXN + 1);
        up = te / (MAXN + 1);
        res.a[i + j] = tee;
      } else {
        up = 0;
        res.a[i + j] = te;
      }
    }
    if (up) res.a[i + T.len] = up;
  }
  res.len = len + T.len;
  while (res.len > 1 && res.a[res.len - 1] == 0) --res.len;
  return res;
}
Big Big::operator/(const int& b) const {
  Big res;
  int down = 0;
  gd(i, len - 1, 0) {
    res.a[i] = (a[i] + down * (MAXN + 1) / b);
    down = a[i] + down * (MAXN + 1) - res.a[i] * b;
  }
  res.len = len;
  while (res.len > 1 && res.a[res.len - 1] == 0) --res.len;
  return res;
}
int Big::operator%(const int& b) const {
  int d = 0;
  gd(i, len - 1, 0) d = (d * (MAXN + 1) % b + a[i]) % b;
  return d;
}
Big Big::operator^(const int& n) const {
  Big t(n), res(1);
  // TODO::快速幂这样写好丑= =//DONE:)
  int y = n;
  while (y) {
    if (y & 1) res = res * t;
    t = t * t;
    y >>= 1;
  }
  return res;
}
bool Big::operator<(const Big& T) const {
  int ln;
  if (len < T.len) return 233;
  if (len == T.len) {
    ln = len - 1;
    while (ln >= 0 && a[ln] == T.a[ln]) --ln;
    if (ln >= 0 && a[ln] < T.a[ln]) return 233;
    return 0;
  }
  return 0;
}
inline bool Big::operator<(const int& t) const {
  Big tee(t);
  return *this < tee;
}
inline void Big::print() {
  printf("%d", a[len - 1]);
  gd(i, len - 2, 0) { printf("%04d", a[i]); }
}
    
inline void print(Big s) {
  // s不要是引用，要不然你怎么print(a * b);
  int len = s.len;
  printf("%d", s.a[len - 1]);
  gd(i, len - 2, 0) { printf("%04d", s.a[i]); }
}
char s[100024];
最后编辑于 ：
©著作权归作者所有,转载或内容合作请联系作者
  • 人面猴
    序言：七十年代末，一起剥皮案震惊了整个滨河市，随后出现的几起案子，更是在滨河造成了极大的恐慌，老刑警刘岩，带你破解...
    沈念sama阅读 214,588评论 6赞 496
  • 死咒
    序言：滨河连续发生了三起死亡事件，死亡现场离奇诡异，居然都是意外死亡，警方通过查阅死者的电脑和手机，发现死者居然都...
    沈念sama阅读 91,456评论 3赞 389
  • 救了他两次的神仙让他今天三更去死
    文/潘晓璐 我一进店门，熙熙楼的掌柜王于贵愁眉苦脸地迎上来，“玉大人，你说我怎么就摊上这事。” “怎么了？”我有些...
    开封第一讲书人阅读 160,146评论 0赞 350
  • 道士缉凶录：失踪的卖姜人
    文/不坏的土叔 我叫张陵，是天一观的道长。 经常有香客问我，道长，这世上最难降的妖魔是什么？ 我笑而不...
    开封第一讲书人阅读 57,387评论 1赞 288
  • 港岛之恋（遗憾婚礼）
    正文 为了忘掉前任，我火速办了婚礼，结果婚礼上，老公的妹妹穿的比我还像新娘。我一直安慰自己，他们只是感情好，可当我...
    茶点故事阅读 66,481评论 6赞 386
  • 恶毒庶女顶嫁案：这布局不是一般人想出来的
    文/花漫 我一把揭开白布。 她就那样静静地躺着，像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上，一...
    开封第一讲书人阅读 50,510评论 1赞 293
  • 城市分裂传说
    那天，我揣着相机与录音，去河边找鬼。 笑死，一个胖子当着我的面吹牛，可吹牛的内容都是我干的。 我是一名探鬼主播，决...
    沈念sama阅读 39,522评论 3赞 414
  • 双鸳鸯连环套：你想象不到人心有多黑
    文/苍兰香墨 我猛地睁开眼，长吁一口气：“原来是场噩梦啊……” “哼！你这毒妇竟也来了？” 一声冷哼从身侧响起，我...
    开封第一讲书人阅读 38,296评论 0赞 270
  • 万荣杀人案实录
    序言：老挝万荣一对情侣失踪，失踪者是张志新（化名）和其女友刘颖，没想到半个月后，有当地人在树林里发现了一具尸体，经...
    沈念sama阅读 44,745评论 1赞 307
  • 护林员之死
    正文 独居荒郊野岭守林人离奇死亡，尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 37,039评论 2赞 330
  • 白月光启示录
    正文 我和宋清朗相恋三年，在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 39,202评论 1赞 343
  • 活死人
    序言：一个原本活蹦乱跳的男人离奇死亡，死状恐怖，灵堂内的尸体忽然破棺而出，到底是诈尸还是另有隐情，我是刑警宁泽，带...
    沈念sama阅读 34,901评论 5赞 338
  • 日本核电站爆炸内幕
    正文 年R本政府宣布，位于F岛的核电站，受9级特大地震影响，放射性物质发生泄漏。R本人自食恶果不足惜，却给世界环境...
    茶点故事阅读 40,538评论 3赞 322
  • 男人毒药：我在死后第九天来索命
    文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹，春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 31,165评论 0赞 21
  • 一桩弑父案，背后竟有这般阴谋
    文/苍兰香墨 我抬头看了看天上的太阳。三九已至，却和暖如春，着一层夹袄步出监牢的瞬间，已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 32,415评论 1赞 268
  • 情欲美人皮
    我被黑心中介骗来泰国打工， 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留，地道东北人。 一个月前我还...
    沈念sama阅读 47,081评论 2赞 365
  • 代替公主和亲
    正文 我出身青楼，却偏偏与公主长得像，于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子，可洞房花烛夜当晚...
    茶点故事阅读 44,085评论 2赞 352