大数乘法巧解

下面这个表是用这个算法计算 981　X　1234　=　1210554的实例（蓝色字体标记的是乘数和被乘数，红色背景部分是中间过程，红色字体是结果）。

（这个图用OpenOffice和dia做出来，有点麻烦哦，不知道有没有什么更方便的工具？）

　　很容易看出规律出来：

　　１、粉红色格式里头的数为对应列和行位置上数的乘积。如第一行，9=9 X 1，这样就可以求出所有粉红色格子里头的数。
　　２、把粉红色格子里头的数的不同位分别放在红色斜线的两侧，如果十位没有，那么用零补充。
　　３、从右下角开始，我们把同处于斜线一侧的数相加，结果存放到对应的左下部位的格子里头。如4=4, 5=2+0+3, 5=6+3+4+0+2 mod 10[进位为一], 0=1[进位]+3+7+2+6+0+1 mod 10,以此类推,直到左上角。

　　实际上，仔细分析以下，我们就会这种乘法和我们小学就学习过的“竖式乘法”原理是一样的，所以肯定是正确的。但是，这种描述方式给我们带来了算法设计和实现上的突破。

　　下面我们就需要把它实现了。
　　首先，我们定义数据的存储方式，乘数和被乘数直接存放到两个一維数组里头，中间过程（即粉红色背景部分）用一个三维数据来存放，存放形式如下：

　　

　　假设该三维数组的名称为Rect，则有如下关系：

Rect[3][2][1] = 4 Rect[3][2][0] = 0 Rect[3][1][1] = 2 Rect[2][2][1] = 3 ...

　　很容易发现这样的规律:
　　１、计算结果的最后一位的下标最大，其下标之和为6；倒数第二位的下标之和为5，以此类推，第一位的下标之和是0。因此，下标之和的顺序刚好是结果位的位序，我们可以据此求出结果的各个位。
　　２、最大的下标之和刚好是乘数位数和被乘数位数之和减一，结果位数是乘数位数和被乘数位数之和。据此，我们可以定义结果数组和一维数组，大小是乘数和被乘数位数之和的一维数组来存放结果。

　　下面我们就根据这种思路写出我们的大数乘法的function，然后再通过调用该function，解决上面的题目——即求出大整数的幂。

Code

[1] lnmp.c -- count the product of two large numbers' multiplication

/* set ts=4 */ /*  * FILE: lnmp.c  * Funciton: Count the product of two large numbers' multiplication  */ /*  * Author: Wu Zhangjin, <zhangjinw@gmail.com>  * (c) 2006-12-31 dslab.lzu.edu.cn Lanzhou University  * License GPL Version 2  *  */ #include <stdio.h>        /* standard input/output functions */ #include <string.h>        /* needed by strlen function */ #define MAX_DIGIT 2000        /* define the max digit of number here*/ /**  * mulof_ln - count the product of two large numbers' multiplication  * @ln_arr: the address of two large numbers  * @pofm: the address of the result array  *  * Description:  * the function use the algorithm: multiplication using a rectangle  * like this: 981 X 1234 = 1210554  * 9 8 1  * 1 0/9 0/8 0/1 1  * 2 1/8 1/6 0/2 2  * 1 2/7 2/4 0/3 3  * 0 3/6 3/2 0/4 4  * 5 5 4  * you can find the the rule of it easily.  *  * Returns:  * the address of the product of multiplication.  *  */ char* mulof_ln(char ln_arr[2][MAX_DIGIT], char pofm[2*MAX_DIGIT-1]) {     int digit_arr[2], digit_pofm, i, j, k, m, tmp, carry;     char td_arr[MAX_DIGIT][MAX_DIGIT][2];          carry = 0;     digit_arr[0] = strlen(ln_arr[0]);     digit_arr[1] = strlen(ln_arr[1]);     digit_pofm = digit_arr[0] + digit_arr[1] - 1;     /*count the intermediate result, save them to a three-dimensional array*/     for (i = digit_arr[0] - 1; i >= 0; i --)     for (j = digit_arr[1] - 1; j >= 0; j --)     {         tmp = (ln_arr[0][i] - '0') * (ln_arr[1][j] - '0');         td_arr[j][i][0] = tmp / 10;         td_arr[j][i][1] = tmp % 10;     }     /*count the result, save it to the result array*/     for (m = digit_pofm; m >= 0; m--)     {         tmp = carry;         for (i = digit_arr[0] - 1; i >= 0; i --)         for (j = digit_arr[1] - 1; j >= 0; j --)             for (k = 1; k >= 0; k --)         if (m == i+j+k) tmp = tmp + td_arr[j][i][k];         carry = tmp / 10;         pofm[m] = tmp%10 + '0';     }     /*take out the ZERO before the result*/     tmp = 0;     for (i = 0; i < digit_pofm; i ++) {         if (pofm[i] != '0') break;         else tmp ++;     }     for (j = tmp, i = 0; j <= digit_pofm + 1; j ++, i ++)         pofm[i] = pofm[j];     return pofm;         } /**  * main - the main function for testing the function: mulof_ln  */ int main() {     char ln_arr[2][MAX_DIGIT];     char pofm[2*MAX_DIGIT-1];    /* the product of the multiplication */     printf("Please input two large numbers: \n");     scanf("%s %s", ln_arr[0], ln_arr[1]);              printf("The product of the muliplication is:\n%s\n", mulof_ln(ln_arr, pofm));          return 0; }

Demo

shell> make lnmp cc     lnmp.c   -o lnmp shell> ./lnmp Please input two large numbers: 2 5 The product of the muliplication is: 10 shell> ./lnmp Please input two large numbers: 22222222222222222222222222222222222222222222222222222222222222222222222222 5 The product of the muliplication is: 111111111111111111111111111111111111111111111111111111111111111111111111110 shell> ./lnmp Please input two large numbers: 0000000000000000000000000000000000000000000000002 000000000000000000000000000000000000000000000055 The product of the muliplication is: 110

[2] rton.c -- count the result of r to the nth(r ^ n)
　　说明：这里并不是直接调用上面的大数乘积的function，为了更方便地操作和节省内存，进行了一定的调整。

/* set ts=4 */ /*  * FILE: rton.c  * Funciton: count the r to the nth  *  */ /*  * Author: Wu Zhangjin, <zhangjinw@gmail.com>  * (c) 2006-12-31 dslab.lzu.edu.cn Lanzhou University  * License GPL Version 2  *  */ #include <stdio.h>        /* standard input/output functions */ #include <string.h>        /* needed by strlen function */ #define MAX_DIGIT 10        /* define the max digit of number here*/ #define MAX_N 25 /**  * mulof_ln - count the product of two large numbers' multiplication  * @m0: the address of two large numbers  * @m1: the address of two large numbers  * @pofm: the address of the result array  *  * Description:  * the function use the algorithm: multiplication using a rectangle  * like this: 981 X 1234 = 1210554  * 9 8 1  * 1 0/9 0/8 0/1 1  * 2 1/8 1/6 0/2 2  * 1 2/7 2/4 0/3 3  * 0 3/6 3/2 0/4 4  * 5 5 4  * you can find the the rule of it easily.  *  * Returns:  * the address of the product of multiplication.  *  * Notes:  * the index of td_arr must be engough large, when you use it to count the power like R ^ N  */ char* mulof_ln(char* m0, char* m1, char pofm[2*MAX_DIGIT-1]) {     int digit_arr[2], digit_pofm, i, j, k, m, tmp, carry;     char td_arr[MAX_DIGIT*MAX_N][MAX_DIGIT*MAX_N][2];          carry = 0;     digit_arr[0] = strlen(m0);     digit_arr[1] = strlen(m1);     digit_pofm = digit_arr[0] + digit_arr[1] - 1;     /*count the intermediate result, save them to a three-dimensional array*/     for (i = digit_arr[0] - 1; i >= 0; i --)     for (j = digit_arr[1] - 1; j >= 0; j --)     {             tmp = (m0[i] - '0') * (m1[j] - '0');         td_arr[j][i][0] = tmp / 10;         td_arr[j][i][1] = tmp % 10;     }     /*count the result, save it to the result array*/     for (m = digit_pofm; m >= 0; m--)     {         tmp = carry;         for (i = digit_arr[0] - 1; i >= 0; i --)         for (j = digit_arr[1] - 1; j >= 0; j --)             for (k = 1; k >= 0; k --)         if (m == i+j+k) tmp = tmp + td_arr[j][i][k];         carry = tmp / 10;         pofm[m] = tmp%10 + '0';     }     /*take out the ZERO before the result*/     tmp = 0;     for (i = 0; i < digit_pofm; i ++) {         if (pofm[i] != '0') break;         else tmp ++;     }     if (tmp > 0)     for (j = tmp, i = 0; j <= digit_pofm + 1; j ++, i ++)         pofm[i] = pofm[j];     return pofm;         } /**  * mypower - count the R to the nth  * @r: the base of the power, it is a float number.  * @n: the exponent of the power  * @p: the address of the result  *  * Descript:  * call the mulof_ln function to count the large float number R to the nth  */ char* mypower(char* r, int n, char p[MAX_DIGIT*MAX_N]) {     int digit_r, digit_p, i, j, pos_point, pos_point_result, tmp;     memset(p, 0, MAX_DIGIT*MAX_N);     digit_r = strlen(r);     for (i = 0; i < digit_r; i ++)         if (r[i] == '.') break;     pos_point = i;     /*if R is a decimal, change it to an integer firstly*/     if (pos_point != digit_r) {             pos_point_result = (digit_r - i - 1) * n;             for (i = pos_point; i <= digit_r; i++)                 r[i] = r[i+1];     }     /*take out the ZERO before*/     tmp = 0;     for (i = 0; i < digit_r -1; i ++) {         if (r[i] != '0') break;         else tmp ++;     }     if (tmp > 0)     for (j = tmp, i = 0; j <= digit_r + 1; j ++, i ++)         r[i] = r[j];     strcpy(p, r);     for (i = 1; i < n; i++)         mulof_ln(r, p, p);          //normalize the output     digit_p = strlen(p);     if (pos_point != digit_r)     {         tmp = digit_p - pos_point_result;         if (tmp >= 0) {             j = 1;             for (i = digit_p; i > tmp; i --) {                 if (j && p[i-1] == '0') p[i-1] = '\0';                 else {                     p[i] = p[i-1];                     j = 0;                 }             }             if (j==0) p[tmp] = '.';         } else {             tmp = -tmp;             j = 1;             for (i = digit_p; i >= 0; i --) {                 if (j && p[i-1] == '0') p[i-1] = '\0';                 else {                     p[i+tmp] = p[i-1];                     j = 0;                 }             }             for (i = 1; i <= tmp; i ++)                 p[i] = '0';             p[0] = '.';         }     }     return p; } /**  * main - the main function for testing the function: mulof_ln  */ int main() {     char p[MAX_DIGIT*MAX_N];     char r[MAX_DIGIT];     int n;     while(scanf("%s %d", r, &n)!=EOF){         printf("%s\n", mypower(r, n, p));         memset(p, 0, MAX_DIGIT*MAX_N);         memset(r, 0, MAX_DIGIT);     }          return 0; }

Demo
　　测试数据：文件名test.txt（说明：下面的测试结果我增加了两条测试数据，非常重要的两条）

95.123 12 0.4321 20 5.1234 15 6.7592 9 98.999 10 1.0100 12