超长整数的基础运算 算法实现之乘除篇

Posted 2020-09-19 jzdwajue

笔算乘法：

对于m位和n位的输入。传统的乘法须要m*n次主要的乘法，也即算法复杂度为O()。我们用纸和笔做乘法运算时，用乘数的每一位乘以被乘数的每一位并加上上一列的进位而产生一行适当移位的中间结果。然后再将各行中间结果相加即得到乘法的终于结果。比如10进制下计算189*34的步骤例如以下表：

 笔算乘法的运算过程

本算法依照上述过程进行计算。但在计算机上最好是把内部的乘法和加法并行的运行。即计算每一行中间结果的同一时候将该行加到终于结果上。这样既可以省去不少步骤，也避免了中间结果的空间开销和内存管理操作，这个简单的优化可以一定程度上提高效率。

本算法的伪代码例如以下表：

输入：分别含有n和t个B进制数位的大整数x、y 输出：B进制的乘积

1. i从0到n+t-1,置Wi = 0 2. i从0到t 2.1 j从0到n 计算wi+j=wi+j + xj * yi 若 wi+j 不小于进位值， 则又一次格式化当前中间值表示 3. 返回

因为两个16位的乘积将超过16位所能表示的最大数，我们知道一个32位的数可以完整地保存两个16位数的乘积，上面的伪码中用一个32位来保存两个16位的数乘积再加上两个单字后的结果，以下讨论一下它的安全性（是否会产生溢出）：上面表达式计算的最大的结果为(2¹⁶ - 1)*(2¹⁶- 1) + (2¹⁶ - 1)。化简后为2³² - 2¹⁶。小于一个32位数的最大值。所以一个32位数可以保存该表达式的结果而不产生溢出。程序不会损失精度。

在乘法中有一类乘法较为特殊--自平方。针对这类乘法使用一般的计算法自然是可行的，可是能够採用更加快捷的方式。使用一次遍历整数每一位做乘法后就可以解决一般乘法须要嵌套循环两次才干处理的结果，由于自平方的两个乘数是“对称”的。以A（a1a2a3a4）为例：每一位的乘法运行后都会有相应的位反复操作一次，比方a1*a4会运行两次，并且相应的乘积结果在同样的结果位置上累加。

所以在此基础上自平方的算法效率要比一般乘法要高。

/*
乘法运算
*/
int mulHBInt(HBigInt *product, HBigInt *biA, HBigInt *biB){
	HBigInt biT;		//计算结果先保存在暂时变量中
	register un_short *pWordA = biA->pBigInt;
	register un_short *pWordB = biB->pBigInt;
	register un_short *pPrd = NULL;
	long result = 0, i = 0, j = 0;

    //初始化暂时大整数变量
	if((result = initHBInt(&biT,biA->length + biB->length)) != RETURN_OK_BINT)
		return result;
	
	biT.length = biA->length + biB->length;
	biT.sign = (biA->sign == biB->sign) ? 1: -1; //同号为正，异号为负
	pPrd = biT.pBigInt;
	
	for(i=0; i<biA->length; ++i) { // 按传统的纸笔乘法的方式计算乘积
		for(j=0; j<biB->length; ++j) {
			pPrd[i+j] = pPrd[i+j] + pWordA[i] * pWordB[j];
			// 假设超出单位长度能表示最大的值（本例中是2<<16）,则进行一次格式化处理
			if(pPrd[i+j] >> BIT_PRE_WORD) formatHBInt(&biT);
		}
	}
	//formatHBInt(&biT);

	trimHBInt(&biT);//去除高位无效的0

	extendHBInt(product,biT.length);
	assignHBInt(product,&biT);
	deleteHBInt(&biT); //清除暂时变量

	return RETURN_OK_BINT;
}

笔算除法：

此算法是模拟笔算除法的形式。变形而来。依据笔算除法的特性。我们知道最重要的也是较难的一步是，对商的估算。

我们做除法运算，是依据被除数和除数最高几位来试商的，此算法也是使用这样的方式。在试商时。为了提快速度。尝试使用的商不可能是从1開始直到i，i满足。

本运算使用二分查找法来试商。从而提高了速度。

       除试商外，另外一个难点就是对除数进行位对齐。多数操作中，除数的位数小于被除数，所以在进行试商前须要对除数进行左移，即人为扩大除数的倍数（Bⁿ）。至于左移详细几位，须要对高位进行比較后确定。规则例如以下：

输入：两个B进制（位数M、N）的大整数x（除数）、y（被除数）       x=m1m2..mM，y=n1n2…nN 输出：左移位数L

1.      默认情况下M < N，故初始化L=N-M 2.      比較两个数的高位部分： 2.1  比較m1 与 n1 若m1 > n1 则L=L-1 若m1 = n1 则置i从0..L           若mi > ni则L=L-1 否则，结束循环 3. 返回L

 

/*
除法运算
a 被除数
b 除数
c 商
d 余数
*/
int divHBInt(HBigInt *a, HBigInt *b, HBigInt *c, HBigInt *d){
	HBigInt ta,tb,tc,temp,tc_1;  //暂时变量
	long result=0,first=1,middle=0,last=0;
	un_short mark=0;
	long len_ta,len_tb,len,i;
	int re;

	if(0 == b->length || (1 == b->length && 0 == b->pBigInt[0])) return -1; //除数不能为0
	//除数为1。则商等于被除数，余数为0
	if(1 == b->length && 1 == b->pBigInt[0]) {
		assignHBInt(c,a);
		setZeroHBInt(d);
		hbi_add_int(d,0);
		return RETURN_OK_BINT;
	}
	
	re = compareHBInt(a,b);
	if(-1 == re) { //除数大于被除数，商为0。余数为被除数
		setZeroHBInt(c);
		assignHBInt(d,a);
		return RETURN_OK_BINT;
	}
	
	if(0 == re) { //除数等于被除数，商为1，余数为0
		setZeroHBInt(c);
		setZeroHBInt(d);
		hbi_add_int(c,1);
		hbi_add_int(d,0);
		return RETURN_OK_BINT;
	}

	//初始化暂时变量
	initHBInt(&ta,INITIAL_BINT);
	initHBInt(&tb,INITIAL_BINT);
	initHBInt(&tc,INITIAL_BINT);
	initHBInt(&temp,INITIAL_BINT);
	initHBInt(&tc_1,INITIAL_BINT);

	len_ta=ta.length;
	len_tb=ta.length;
	len=len_ta-len_tb;
	i=0;

	assignHBInt(&ta,a);				//对暂时变量赋值
	assignHBInt(&tb,b);
	hbi_add_int(&tc,1);				//初始商为1
	extendHBInt(&temp,ta.length);	//扩展暂时变量长度
	extendHBInt(&tc_1,ta.length);

	//二分法试商计算
	while(compareHBInt(&ta,b) > 0 ) {  // 保证ta大于b时循环
		len_ta=ta.length;
		len_tb=tb.length;
		len=len_ta-len_tb;
		i=0;
		mark=0;
		
		if(ta.pBigInt[len_ta-1] <= tb.pBigInt[len_tb-1]) {
			while(i<len_tb){
				if(ta.pBigInt[len_ta-1-i] < tb.pBigInt[len_tb-1-i]) {
					mark=ta.pBigInt[len_ta-1];
					len--;
					i++;
					break;
				} else i++;
			}
		}

		Left_shift_word(&tb,len);	//对除数进行左移位。使其跟被除数有同样数位
		Left_shift_word(&tc,len);	//商同一时候左移位
		
		result=mark*CARRY_RADIX+ta.pBigInt[len_ta-1-i];
		first=1;
		last= result / tb.pBigInt[len_tb-1+len] + 1;	//二分查找法的上限
		middle=(last+first)/2 + 1;						//二分查找的中间值

		for(; (first < last) && (first != middle); ) {	//使用二分查找法来试商
			assignHBInt(&temp,&tb);
			hbi_mul_radix(&temp,middle);
			
			if(-1 == compareHBInt(&ta,&temp)) last=middle;
			else first=middle;
			
			middle=(last+first)/2;
		}

		assignHBInt(&temp,&tb);
		hbi_mul_radix(&temp,middle);
		subHBInt(&ta,&ta,&temp);	//被除数-除数*商
		
		hbi_mul_radix(&tc,middle);  //得到所试的商
		
		addHBInt(&tc_1,&tc_1,&tc);	//对商进行累加

		setZeroHBInt(&tc);          //暂时保存商的值清0
		hbi_add_int(&tc,1);			//对商的暂时变量赋1      

		setZeroHBInt(&temp);	     
		assignHBInt(&tb,b);
	}
	
	if(compareHBInt(&ta,b) == 0) {
		hbi_add_int(c,1);
		setZeroHBInt(d);
	} else {
		if(c!=NULL)	assignHBInt(c,&tc_1);	//得到除法运算的商
		if(d!=NULL) assignHBInt(d,&ta);		//除法运算的余数
	}
	
	// 回收暂时变量空间
	deleteHBInt(&ta); 
	deleteHBInt(&tb);
	deleteHBInt(&tc);
	deleteHBInt(&tc_1);
	deleteHBInt(&temp);

	return RETURN_OK_BINT;
}