高精度运算

之前早就想把学过的算法记录下来，但是一直没有时间。最近在给五年级的小学生上OI的算法课，所以正好可以把所思所想留存下来。

一、为什么需要高精度运算

当要处理的数据超过了任何一种数据类型所能容纳的范围，这种数称为高精度数，必须自定义类型来存储。同时还要自行编写高精度数的运算函数（加\减\乘\除等）。

基本整数类型的取值范围

类型	数值范围	占字节数
char	-128 .. 127	1
unsigned char	0 .. 255	1
int (long)	-2147483648 .. 2147483647	4
unsigned int (long)	0 .. 4294967295	4
long long	-9223372036854775808 .. 9223372036854775807	8
unsigned long long	0 .. 18446744073709551615	8

二、高精度数的输入和存储

在运算对象的数值范围为任何数据类型所无法容纳的情况下，采用数组（每一个元素对应一位十进制数，由其下标顺序指明位序号）来表示一个数，就称为高精度数。

1. 采用字符串形式输入，并将其转化为整数数组。
2. 用一个整型变量记录数据的实际长度（即数组的元素个数）

字符串到整数数组的转换

字符串存储时，数的高位被存放在字符串的低位。转换成整数数组时，要把高精度数的低位“还原”到数组的低位。这样才便于后续计算。a[1]存放最低位，a[6]存放最高位。高精度数的位数可存放在a[0]中。也可以另用一个变量存放。

例如
字符串s: 763218
转化成整数数组：6812367
其中数组第一个元素6是数组的位数

高精度数类型的定义与读入

const int MAXLEN = 241;     //最大长度为240位
typedef int hp[MAXLEN];     //定义hp为高精度类型
void Str2hp(string s,hp &a)   //s转换后存入a
{
    int i, len;
    memset(a, 0, sizeof(a));    //clear a
    len = s.size();
    a[0] = len;     //save the length
    for (i=0; i<len; i++)   //convert
        a[len-i] = s[i] - '0';
}

输出

void Print(hp a)
{
    int i, len;
    len = a[0];     //get the length
    for (i=len; i >=1; i--)
        cout << a[i];
    cout << endl;
}

三、高精度加法

问题表述：输入a，b（<10^240）两个数，输出a+b的值。
样例输入：
99999999999999999999
12345678901234567890
样例输出：
112345678901234567889

算法思想类似于竖式加法运算，注意进位处理。把计算结果存到c中：

1. 先计算：直接将a[i]和b[i]相加，结果放在c[i]中。
2. 处理进位：从最低位开始，逐位整理：把本位模10，向高一位进位。

c[i+1] += c[i] / 10;
c[i] = c[i] % 10;

3. 去掉最高位的0：因为两数相加,也可能不产生进位，因此要把这种情况下最高位的0去掉，其实就是减小和的长度len的值。

while ( (len>1) && (c[len] ==0)) 
    len--;

void add(hp a, hp b, hp &c) //Add a, b to c
{
    hp d;
    int i, len;
    memset(d, 0, sizeof(d));    //d清0
    if (a[0] > b[0]) len = a[0];//求和的位数
    else  len = b[0];
      len++;                    //位数加1

    for (i=1; i<=len; i++)    //逐位相加
        d[i] = a[i] + b[i];

    for (i=1; i<len; i++)     //处理进位
    {
        d[i+1] += d[i] / 10;
        d[i] %= 10;
    }
    while ( (len > 1) && (d[len] == 0)) //处理最高位
        len--;
    d[0] = len;
    memcpy(c, d, sizeof(d));    //保存结果
}

四、高精度减法

问题表述：输入a，b（<10^240）两个数，输出a-b的值。
样例1输入：
456
409
样例1输出：
47
样例2输入：
999
1000
样例2输出：
-1

1. 比较a和b的大小，从而确定结果的正负号
2. a、依照由低位至高位的顺序进行减法运算；b、在每一次位运算中，若出现不够减的情况(a[i]<b[i])，则向高位借位；c、在进行了减运算后，若高位为0，则要减少结果的长度。在具体计算过程中，仍然用abc三步走的办法。

void sub(hp a, hp b, hp &c) //a must be greater than b
{
    int i, len;
    hp d;
    memset(d, 0, sizeof(d));   //clear d
    len = a[0];
    for (i=1; i<=len; i++)
        d[i] = a[i] - b[i];
  for (i=1; i<=len; i++)
        if (d[i] < 0)   //处理借位
        {
            d[i] += 10;
            d[i+1] -= 1;
        }

    while ( (len > 1) && (d[len] == 0) ) //整理差的长度
        len--;
    d[0] = len;
    memcpy(c, d, sizeof(d));
}

注意：保存结果的数组使用前一般都要清零；循环控制变量一般用 i 、 j 、 k，一般不要用m、n等，长度用len表示，这样容易找错误。

五、比较两个高精度数大小

当a>b，返回1； a<b,返回-1；a==b，返回0

int compare(hp a, hp b)
{
    int i;
    if (a[0] > b[0]) return 1;
    if (a[0] < b[0]) return -1;
    for (i=a[0]; i>=1; i--)
        if (a[i] > b[i]) return 1;
        else if (a[i] < b[i]) return -1;
    return 0;
}

例题：求Fibonacci数列的第n项
Fibonacci数列的代表问题是由意大利著名数学家Fibonacci于1202年提出的“兔子繁殖问题”(又称“Fibonacci问题”)。问题：有雌雄一对兔子，假定过两个月后便每个月可繁殖雌雄各一的一对小兔子。问过n个月后共有多少对兔子？已知：N <= 1000
当N<=93时，可以用基本类型long long解决，但是大于93就需要使用到高精度算法。

void fibo (int n, hp & ans)
{
    hp f0, f1, t;
    int i;
    memset(ans, 0, sizeof(ans));
    f0[0] = 1; f0[1] = 1;
    f1[0] = 1; f1[1] = 1;
    if ( (n == 1) || (n == 2) )
    {
        ans[0] = 1;
        ans[1] = 1;
        return ;
    }
    for (i=3; i<=n; i++)
    {
        add(f0, f1, t);
        memcpy(f0, f1, sizeof(f1));
        memcpy(f1, t, sizeof(t));
    }
    memcpy(ans, f1, sizeof(f1));
}

六、高精度数乘以整形数运算

问题表述：精确计算n的阶乘n!（7<n<80）
样例输入：
20
样例输出：
2432902008176640000

1. 估算n!所需的高精度数组长度
2. 被乘数（高精度）从低位向高位逐位乘以乘数（整数）
   因为80!<80^80<10080=(10²⁾80=10^160，所以80！可以用160个数组元素a[1]a[2]…a[160]来存放，一个数组元素存放一个数位上的数字。同样的方法，可以估算出100！可以用200位的数组来存放。

void multiply(hp a, int b,hp &c)
{
    hp d;
    int i, len, t;

    memset(d, 0, sizeof(d));
    len = a[0];
    for (i=1; i<=len; i++)
    {
        t = a[i] * b + d[i];
        d[i] = t % 10;
        d[i+1] = t / 10;
    }
    len++;
    while (d[len])
    {
        d[len+1] = d[len] / 10;
        d[len] %= 10;
        len++;
    }
    while ( (d[len] == 0) && ( len > 1))
        len--;
    d[0] = len;
    memcpy(c, d, sizeof(d));
}
int main()
{
    int n, i;
    hp ans;
    cin >> n;
    ans[0] = 1;
    ans[1] = 1;
    for (i=2; i<=n; i++)
        multiply(ans, i, ans);
    for (i=ans[0]; i>=1; i--)
        cout << ans[i];
    cout << endl;
    return 0;
}

七、高精度数乘以高精度数

问题表述：输入a，b（<10^100）两个数，输出a*b的值。
样例输入：
12345678900
98765432100
样例输出：
1219326311126352690000

1. 积的位数为lena+lenb-1或者lena+lenb；
2. 如果暂且不考虑进位关系，则aibj应该累加在积的第j+i-1位上：c[i+j-1] := a[i]b[j] + c[i+j-1];
3. 可以先乘、后处理进位

void high_multiply(hp a, hp b, hp &c)
{
    hp d;
    int i, j, len;

    memset(d, 0, sizeof(d));
    for (i=1; i<= a[0]; i++)
        for (j=1; j<= b[0]; j++)
            d[i+j-1] += a[i] * b[j];

    len = a[0] + b[0] + 1;
    for (i=1; i<= len; i++)
    {
        d[i+1] += d[i] / 10;
        d[i] = d[i] % 10;
    }
    while ((d[len] == 0) && (len > 1))
        len--;
    d[0] = len;
    memcpy(c, d, sizeof(d));
}
int main()
{
    int i;
    hp a, b, ans;
    string s1, s2;
    cin >> s1;
    cin >> s2;
    str2hp(s1, a);
    str2hp(s2, b);
    high_multiply(a, b, ans);
    for (i=ans[0]; i>=1; i--)
        cout << ans[i];
    cout << endl;
    return 0;
}

八、高精度数除以整形数

问题表述：输入a（<10^240） ，b（<10^9）两个数，输出a/b的值和余数。
样例输入：
99887766554433221100
2001920
样例输出：
49895983133408
1077740
基本方法是从被除数的最高位开始，逐位计算商和余数。存储商。余数则转移到下一次的被除数中。注意在下一次的计算中，余数要乘以10，然后再加上当前位的被除数。

void divide(hp a, int b, hp &c, int &r)
{
    hp d;
    int i, len;

    memset(d, 0, sizeof(d));
    r = 0;
    len = a[0];

    for (i = len; i>=1; i--)
    {
        r = r * 10 + a[i];
        d[i] = r / b;
        r = r % b;
    }
    while ( (len > 1) && (d[len] == 0))
        len--;
    d[0] = len;
    memcpy(c, d, sizeof(d));
}

九、练习题

阶乘之和，求1！+2！+3！+…+n！
输入数据(sumfac.in)：
一行，一个整数n （n<=100）
输出数据(sumfac.out):
一行，1~n的阶乘之和。
输入样例：
5
输出样例：
153