大数据乘法

大整数乘法比较容易想到的是做多位数乘法时列竖式进行计算的方法，只要写出模拟这一过程的程序，就能实现任意大整数的乘法运算。 

下面介绍一种更便于编程的方法：“列表法”:

例如当计算8765*234时，把乘数和被乘数照如下列出，见表1：

8	7	6	5	*
16	14	12	10	2
24	21	18	15	3
32	28	24	20	4

 

16	14	12	10
	24	21	18	15
		32	28	24	20
16	38	65	56	39	20

  

	16	38	65	56	39	20
2	16+4=20	38+7=45	65+6=71	56+4=60	39+2=41
留2	留0进2	留5进4	留1进7	留0进6	留1进4	留0进2
2	0	5	1	0	1	0

根据以上思路 就可以编写C程序了，再经分析可得：

1，一个m位的整数与一个n位的整数相乘，乘积为m+n-1位或m+n位。

2，程序中，用三个字符数组分别存储乘数，被乘数与乘积。由第1点分析知，存放乘积的字符数组饿长度应不小于存放乘数与被乘数的两个数组的长度之和。

3，可以把第二步“计算填表”与第三四步“累加进位”放在一起完成，可以节省存储表格2所需的空间。

4，程序关键部分是两层循环，内层循环累计一数组的和，外层循环处理保留的数字和进位。

<pre name="code" class="cpp"><span style="font-size:18px;">#define MAXLENGTH 1000
#include <stdio.h>
#include <string.h>
void compute(char *a,char *b,char *c)
{
	int i,j,m,n;
	long sum,carry;
	m = strlen(a)-1;
	n = strlen(b)-1;
	for(i=m;i>=0;i--)
	{
		a[i] -= '0';
	}
	for(i=n;i>=0;i--)
	{
		b[i] -= '0';
	}
	c[m+n+2] = '\0';
	carry = 0;
	for(i=m+n;i>=0;i--)//倒斜对角线之和的项相加的该位结果
	{
		sum=carry;
		if((j=(i-m))<0)
			j=0;
		for(;j<=i&&j<=n;j++)
			sum += a[i-j]*b[j];
		c[i+1] = sum%10 + '0';
		carry = sum/10;
	}
	if((c[0] = carry+'0')=='0')
		c[0] = ' ';
}
void main()
{
	char a[MAXLENGTH],b[MAXLENGTH],c[MAXLENGTH*2];
	puts("a:");
	gets(a);
	puts("b:");
	gets(b);
	compute(a,b,c);
	puts("Answer:");
	puts(c);
	getchar();
}</span>

效率分析：用以上算法计算m位整数乘以n位整数，需要先进行m*n次乘法，再进行约m+n次加法运算和m+n次取模运算（实为整数除法）。  经过改进，此算法效率可以提高约9倍。

注意到以下事实：8216547*96785 将两数从个位起，每3位分为节，列出乘法表，将斜线间的数字相加：

    8  216  547

         96   785

8	216	547	*
768	20736	52512	96
6250	169560	429395	785

 

768	20736	52512
	6250	169560	429395
768	27016	222072	429395

 

将表中最后一行进行如下处理：从个位数开始，每一个方格里只保留三个数字，超出1000的部分进位到前一个方格里：

	768	27016	222072	429395
	768+27=795	27016+222=27238	222072+429=222501	留395进429
	795	238	501	395

  

所以8216547*96785 = 795238501395

 

也就是说我们在计算生成这个二维表时，不必一位一位的乘，而可以三位三位的乘；在累加时也是满1000进位。这样，我们计算m位整数乘以n位整数，只需要进行m*n/9次乘法运算，再进行约（m+n）/3次加法运算和（m+n）/3次去摸运算。总体看来，效率是前一种算法的9倍。

有人可能会想：既然能用三位三位的乘，为什么不能4位4位的乘，甚至5位。本算法在累加表中斜线间的数字时，如果用无符号长整数（范围0至~4294967295）作为累加变量，在最不利的情况下（两个乘数的所有数字均为9），能够累加约4294967295/（999*999）=4300次，也就是能够准确计算任意两个约不超过12900（每次累加的结果“值”三位，故4300*3=12900）位的整数相乘。如果4位4位地乘，在最不利的情况下，能过累加月4294967295/(9999*9999)=43次，仅能够确保任意两个不超过172位的整数相乘，没什么实用价值，更不要说5位了。

<span style="font-size:18px;">
</span><pre name="code" class="cpp"><span style="font-size:18px;">#include <stdio.h>
#include <string.h>
#include <conio.h>
#include <stdlib.h>
#include <time.h>
#define N 3
FILE *fp;
int max(int a,int b,int c)
{
	int d = (a>b)?a:b;
	return (d>c)?d:c;
}
int initarray(int a[])
{
	int q,p,i;
	q = N + rand()%100;//在linux下stdlib.h包含srandom和random ,但在VC下stdlib.h包含的是srand和rand
	if(q%3 == 0)
		p = q/3;
	else
		p = q/3 + 1;
	for(i=0;i<p;i++)
		a[i] = rand()%1000;
	if(q%3 == 0)
		a[0] = 100 + rand()%900;
	if(q%3 == 2)
		a[0] = 10 + rand()%90;
	if(q%3 == 1)
		a[0] = 1 +rand()%9;
	return p;
}
void write(int a[],int l)
{
	int i;
	char string[10];
	for(i=1;i<l;i++)
	{
		itoa(a[i],string,10);
		if(strlen(string) == 1)
			fprintf(fp,"00");
		if(strlen(string) == 2)
			fprintf(fp,"0");
		fprintf(fp,"%s",string);
		if((i+1)%25 == 0)
			fprintf(fp,"\n");
	}
	fprintf(fp,"\n");
	fprintf(fp,"\n");
}

void main()
{
	int a[5000]={0},b[5000]={0},k[10001]={0};
	unsigned long c,d,e;
	int i,j,la,lb,ma,mi,p,q,t;
	//randomize();
	srand(unsigned(time(NULL)));//randomize()不是随机数生成的函数,而是初始化随机数生成器的函数,而且,他不是一个C语言标准库函数.改用srand(unsigned(time(NULL)));

	la = initarray(a);
	lb = initarray(b);

	if(la<lb)//如果被乘数长度小于乘数，则交换被乘数与乘数
	{
		p = (lb>la)?lb:la;
		for(q=0;q<p;q++)
		{
			t = a[q];
			a[q] = b[q];
			b[q] = t;
		}
		t = la;
		la = lb;
		lb = t;
	}
	c=d=0;
	for(i=la+lb-2;i>=0;i--)//累加斜线间的数，i位横纵坐标之和
	{
		c = d;//将前一位的进位标志存入累加变量C
		ma = max(0,i-la+1,i-lb+1);//求累加的下线
		mi = (i>la)?(la-1):i;//求累加的上线
		for(j=ma;j<mi;j++)
			c+=a[j]*b[i-j];
		d = c/1000;//求进位标志
		if(c>999)
			c%=1000;
		k[i] = c;
	}
	e = k[0] + 1000*d;
	fp = fopen("res.txt","w+");
	fprintf(fp,"%d",a[0]);
	write(a,la);
	fprintf(fp,"%d",b[0]);
	write(b,lb);
	fprintf(fp,"%d",e);
	write(k,la+lb-1);
	fclose(fp);
}</span>

我们将n位的二进制整数X和Y各分为2段，每段的长为n/2位(为简单起见，假设n是2的幂)。由此，X=A2n/2+B ，Y=C2n/2+D。这样，X和Y的乘积为：

XY=(A2n/2+B)(C2n/2+D)=AC2n+(AD+CB)2n/2+BD     (1)

如果按式(1)计算XY，则我们必须进行4次n/2位整数的乘法(AC，AD，BC和BD)，以及3次不超过n位的整数加法(分别对应于式(1)中的加号)，此外还要做2次移位(分别对应于式(1)中乘2n和乘2n/2)。所有这些加法和移位共用O(n)步运算。设T(n)是2个n位整数相乘所需的运算总数，则由式(1)，我们有:

          T(n)=4T(n/2)+θ(n)        (2)

由此可得T(n)=O(n2)。因此，用(1)式来计算X和Y的乘积并不比小学生的方法更有效。要想改进算法的计算复杂性，必须减少乘法次数。为此我们把XY写成另一种形式:

XY=AC2n+[(A-B)(D-C)+AC+BD]2n/2+BD      (3)

虽然，式(3)看起来比式(1)复杂些，但它仅需做3次n/2位整数的乘法(AC，BD和(A-B)(D-C))，6次加、减法和2次移位。由此可得:

                     T(n)=3T(n/2)+θ(n)        (4)

用解递归方程的套用公式法马上可得其解为T(n)=O(nlog3)=O(n1.59)。利用式(3)，并考虑到X和Y的符号对结果的影响，我们给出大整数相乘的完整算法MULT如下:

function MULT(X，Y，n); {X和Y为2个小于2n的整数，返回结果为X和Y的乘积XY}begin  S:=SIGN(X)*SIGN(Y); {S为X和Y的符号乘积}  X:=ABS(X);  Y:=ABS(Y); {X和Y分别取绝对值}  if n=1 then     if (X=1)and(Y=1) then return(S)                      else return(0)         else begin                A:=X的左边n/2位;                B:=X的右边n/2位;                C:=Y的左边n/2位;                D:=Y的右边n/2位;                ml:=MULT(A,C,n/2);                m2:=MULT(A-B,D-C,n/2);                m3:=MULT(B,D,n/2);                S:=S*(m1*2n+(m1+m2+m3)*2n/2+m3);                return(S);              end;end;

上述二进制大整数乘法同样可应用于十进制大整数的乘法以提高乘法的效率减少乘法次数。下面演示算法的计算过程。

设X=314l,Y=5327，用上述算法计算XY的计算过程可列表如下，其中带'号的数值是在计算完成AC，BD，和(A-B)(D-C)之后才填入的。

 

X=3141        A=31       B=41        A-B=-10

Y=5327        C=53       D=27        D-C=-26

           AC=(1643)'

           BD=(1107)'

          (A-B)(D-C)=(260)'

XY=(1643)'104+[(1643)'+(260)'+(1107)']102+(1107)'

  =(16732107)'

 

A=31        A1=3       B1=1        A1-B1=2

C=53        C1=5       D1=3        D1-C1=-2

           A1C1=15     B1D1=3     (A1-B1)(D1-C1)=-4

AC=1500+(15+3-4)10+3=1643

 

B=41        A2=4       B2=1        A2-B2=3

D=27        C2=2       D2=7        D2-C2=5

           A2C2=8     B2D2=7     (A2-B2)(D2-C2)=15

BD=800+(8+7+15)10+7=1107

 

|A-B|=10        A3=1       B3=0        A3-B3=1

|D-C|=26        C3=2       D3=6        D3-C3=4

           A3C3=2     B3D3=0     (A3-B3)(D3-C3)=4

(A-B)(D-C)=200+(2+0+4)10+0=260

 

代码的实现

<span style="font-size:18px;"><span style="font-size:18px;">[cpp] view plaincopy
1. /************************************************************************/  
2. //函数功能：分治法求两个N为的整数的乘积  
3. //输入参数：X,Y分别为两个N为整数  
4. //算法思想：  
5. //时间复杂度为：T(n)=O(nlog3)=O(n1.59)  
6. /************************************************************************/  
7. #define SIGN(A) ((A > 0) ? 1 : -1)  
8. int IntegerMultiply(int X, int Y, int N)  
9. {  
10.     int sign = SIGN(X) * SIGN(Y);  
11.     int x = abs(X);  
12.     int y = abs(Y);  
13.     if((0 == x) || (0 == y))  
14.         return 0;  
15.     if (1 == N)  
16.         return x*y;  
17.     else  
18.     {  
19.         int XL = x / (int)pow(10., (int)N/2);   
20.         int XR = x - XL * (int)pow(10., N/2);  
21.         int YL = y / (int)pow(10., (int)N/2);  
22.         int YR = y - YL * (int)pow(10., N/2);  
23.           
24.         int XLYL = IntegerMultiply(XL, YL, N/2);  
25.         int XRYR = IntegerMultiply(XR, YR, N/2);  
26.         int XLYRXRYL = IntegerMultiply(XL - XR, YR - YL, N/2) + XLYL + XRYR;  
27.         return sign * (XLYL * (int)pow(10., N) + XLYRXRYL * (int)pow(10., N/2) + XRYR);  
28.     }  
29. }  
30. int _tmain(int argc, _TCHAR* argv[])  
31. {  
32.     int x = 1234;  
33.     int y = 4321;  
34.     cout<<"x * y = "<<IntegerMultiply(x, y, 4)<<endl;  
35.     cout<<"x * y = "<<x*y<<endl;  
36.     return 0;  
37. }  </span></span>

<span style="font-size:18px;"><span style="font-size:18px;">#include <iostream>
#include <string>
using namespace std;
//计算两个大整数的加法，结果放到b中(假设b的数组大小比a大)
//由于两个数组的位数都比较多,所以一定不会溢出
void add(int* a,int size1,int*b,int size2)
{
 //存放每一位的进位
 int carry = 0;
 for(int i = size1 - 1;i >= 0;i--)
 {
  int sum = b[i] + a[i] + carry;
 b[i] = sum % 10 ;
  	carry = sum / 10;
 }
 //处理进位
 for(i = size2 - size1 - 1;i >= 0 && carry;i--)
 {
  	int sum = b[i] + carry;
  	b[i] = sum % 10;
  	carry = b[i] / 10;
 }
}
//一个大整数和一个一位数相乘
//将结果放在d中(d的空间一定够)
void tempMult(int *a,int size1,int num,int* d,int sized)
{
 //存放进位信息
 	int carry = 0;
 	for(int i = size1 - 1;i >= 0;i--)
 	{
  	int result = a[i] * num + carry;
  	d[--sized] = result % 10;
  	carry = result / 10;
 	}
 	//将进位放到下一位置上
 	d[--sized] = carry;
}
//两个大整数相乘的算法
void multiply(int* a,int size1,int* b,int size2,int* c)
{
 //存放中间结果
 int* d = new int[size1 + size2];
 	//将d清零
 	memset(d,0,(size1 + size2)*sizeof(int));
 	//标记现在是第几次相乘
 	int times = 0;
 	//执行乘法
 	for(int i = size2 - 1;i >= 0;i--)
 	{
  	tempMult(a,size1,b[i],d,size1 + size2);
  	//如果是第一次的话,直接将c拷贝到d
  	if(!times)
  	{
  	 memcpy(c,d,(size1 + size2)*sizeof(int));
   	times++;
 	 }//否则要先将d移位
  	else
  	{
   	for(int j = times;j < size1 + size2;j++)
    	d[j - times] = d[j];
  	 //将最后几位置零
   	for(int k = j - times;k < size1 + size2;k++)
    	d[k] = 0;
  	 times++;
   	//将c和d相加放到c中
   	add(d,size1 + size2,c,size1 + size2);
  	}
  	//将d清空
  	memset(d,0,(size1 + size2)*sizeof(int));
 	}
}
int main()
{
 string s1,s2;
 cout<<"请输入你要相乘的两个大整数:"<<endl<<"第一个是:"<<endl;
 cin>>s1;
 cout<<"第二个是:"<<endl;
 cin>>s2;
 
 //将大整数放到下面的两个数组中
 int* a,*b;
 int len1 = s1.size();
 int len2 = s2.size();
 a = new int[len1];
 b = new int[len2];
 //将字符串转化为整数数组
 for(int i = 0;i < len1;i++)
  a[i] = s1.at(i) - '0';
 //将字符串转化为整数数组
 for(i = 0;i < len2;i++)
  b[i] = s2.at(i) - '0';
 //将大整数相乘的结果放在c中
 int* c = new int[len1 + len2];
 //先将c数组清零
 memset(c,0,(len1 + len2)*sizeof(int));
 //将两个大整数相乘
 multiply(a,len1,b,len2,c);
 //输出大整数相乘的结果
 //找到最高位
 for(i = 0;i < len1 + len2 - 1;i++)
  if(c[i] != 0)
   break;
 
 cout<<"结果为: "<<endl;
 //输出结果
 for(int j = i;j < len1 + len2;j++)
  cout<<c[j];
 cout<<endl;
 delete []a;
 delete []b;
 delete []c;
 return 0;
}</span></span>