当前位置:文档之家› 无符号大整数相乘优化算法及

无符号大整数相乘优化算法及

Win32下无符号大整数相乘优化算法及其C++实现Lightning[0GiNr]1、问题的引出:两个无符号的大整数相乘是一道实践意味很浓的算法题目,这里的“无符号”(unsigned) 指的是相乘的两个数都是正数,不需要考虑符号。

由于32 位计算机没有指令支持128 及以上二进制位数的大整数的运算,所以必须自己设计算法来计算。

传统的优化算法基本上都是理论层面上的优化,即尽可能地从理论上减少乘法次数,但是往往不能达到预想的优化效果。

比方说二分法优化:将待相乘的整数M分成相等的左右两个部分M1和M2,另一个相乘整数N也同样地分成N1和N2,然后按这样的方法递归分割,直到最后的元素大小小到可以利用CPU旨令直接计算为止。

这时利用公式M*N= (M1 + M2) * (N1 + N2) = M1*N1 + M1*N2 + M2*N1 + M2*N2 结合移位运算再逐层返回得出最终结果。

显然这种算法理论性过强,一来只有当M和N为2的P次方(P为正整数)时的优化才会节省时间,而实际情况下应对随机数据时则会出现大量位移操作,速度不会得到提升;二来使用的递归算法由于调用栈和跳转指令的开销,浪费大量CPU时间;三来这种方法实际上并没有真正地减少乘法次数,因为除了最后一层递归中的乘法可以直接用CPU指令实现,其余各层的乘法由于数值较大仍得另想办法。

由此,我们须从实际出发,探索一些实用的优化方法。

本程序的测试环境为:Windows XP SP2 32bit + 512MB SDRAM + P4 1.80Ghz + VC2、朴素的算法思路:为了简易起见,我们先来设计一个朴素的算法。

使用一个DWOR类型的数组m_buffer作为缓冲区,大小为64,同时声明一个int类型的变量m_nUsed 记录当前缓冲中DWOR使用的个数(即后面所提到的“位数”)。

类的声明如下:代码清单:BigNumber.cpp#include <windows.h>#include <stdio.h>class CBigNumber{public:CBigNumber(){memset(this, 0, sizeof(*this)); m_nUsed = 1;}CBigNumber& operator = (DWORD dwData);CBigNumber& operator *= (const CBigNumber& right);int GetCount() const { return m_nUsed; }const DWORD* GetBuffer() const { return m_buffer; }protected:VOID OffsetAdd(DWORD dwData, int nOffset);int m_nUsed;DWORD m_buffer[64];};首先是赋值函数,这个函数将一个DWOR类型的整数转化到CBigNumber中。

CBigNumber& CBigNumber::operator = (DWORD dwData){memset(this, 0, sizeof(*this));m_nUsed = 1;m_buffer[0] = dwData;return *this;}下面,重载*=运算符,进行大整数乘法运算。

本程序使用双重循环,将被乘数的每一个DWOR位与乘数的每一位相乘,将结果错位累加(类似于手算乘法) 。

由于两个DWOR相乘后结果可能溢出DWOR的表示范围,所以可以将一个DWOR分成两个部分,如0x12345678 分成0x1234 和0x5678 ,另一个被乘数假设是0x44445555 ,则分成0x4444 和0x5555 。

这样,由于(A+B)*(C+D) = AC + (BC + AD) + BD; 结果显然是0x1234*0x4444*0x100000000 + (0x5678*0x4444 +0x1234*0x5555)*0x10000 + 0x5678*0x5555 。

考虑上面的式子,0x1234*0x4444*0x100000000必然是溢出DWORI的,正好可以放到下一个DWORD(0x5678*0x4444 + 0x1234*0x5555)*0x10000 的前16位是溢出的,需要累加到下一位,后16位则不溢出,放到本位。

而0x5678*0x5555 绝对不会溢出,也要加到本位。

这种方法比较“正规” ,稍后我们会介绍更好的办法。

CBigNumber& CBigNumber::operator *= (const CBigNumber& right){// 朴素的乘法CBigNumber csTemp;DWORD A = 0, B = 0, C = 0, D = 0;for(int i = 0; i < right.m_nUsed; i++){for(int n = 0; n < m_nUsed; n++){A = right.m_buffer[i] >> 16;B = right.m_buffer[i] & 0x0000FFFF;C = m_buffer[n] >> 16;D = m_buffer[n] & 0x0000FFFF;DWORD AC = A * C;DWORD BD = B * D;DWORD BC_AD = B * C + A * D;csTemp.OffsetAdd(BD + ((BC_AD & 0x0000FFFF) << 16), i + n); // 强行截断csTemp.OffsetAdd(AC + (BC_AD >> 16), i + n + 1);}}*this = csTemp;return *this;}其中调用的函数OffsetAdd 用于实现错位相加。

防止因为全部相加而引起的性能下降VOID CBigNumber::OffsetAdd(DWORD dwData, const int nOffset){// 判断是否溢出m_buffer[nOffset] += dwData;if(m_nUsed < nOffset + 1) m_nUsed = nOffset + 1;if(m_buffer[nOffset] < dwData) // 按照补码运算法则,这里产生进位OffsetAdd(1, nOffset + 1);}CBigNumber 类设计好了,现在我们要实践一下。

先讨论一个大数与一个较小数相乘,然后在此基础上进行优化。

我们计算:(1+2人0)*(1+2人1)*(1+2人2)* ……*(1+2人31)。

这个积有150多位(十进制),十六进制则有63位,换算成DWOR有16位。

为了看出时间差的效果,我们计算10000遍。

VOID CreateList(int* data){for(int i = 0; i < 32; i++){data[i] = 1 + (1 << i);}}int main(void){int nList[32];memset(nList, 0, 32 * sizeof(int));DWORD dwTickCount = ::GetTickCount();CBigNumber number;CBigNumber addtor;CreateList(nList);for(int c = 0; c < 10000; c++){number = 1;for(int i = 0; i < 32; i++){addtor = nList[i];number *= addtor;}}// printconst DWORD* p = number.GetBuffer();int bit = 0;for(int t = number.GetCount() - 1; t >= 0; t--){if(p[t] != 0){printf("%X", p[t]);bit++;}}printf("\n\nCost Time = %d ms BIT = %d", ::GetTickCount() - dwTickCount, bit);getchar();return 0;}上面的这些代码在DEBUG版本中计算10000次表达式,用时约2000ms,平均每次0.5ms, RELEASED 本用时约500ms,平均每次0.02ms。

3、简单的优化朴素的算法效果已经不错,那么我们还有多少优化余地呢?其实我们大有文章可做。

正如你看到的,为了解决溢出问题,我们将好好的两个数拆开,用了 4 次乘法, 4 次移位,3次加法,效率低得可怜。

网上有资料说这种方法计算BC*AD时可用(A + B) * (C + D) - AC - BD; 来等效BC * AD,用减少一次乘法和添加三次加(减)法来平衡时间。

但是我经过测试,这种所谓的优化方法比原方法更慢!虽然乘法比好几次加法都浪费时间,但是现在的CPU早就加入了UV流水线,将相邻的两个相互不干扰的指令进行配对并行执行;更高级的还有Pentium Pro\II\III 就引进的乱序执行技术(这些处理器现在已经基本不见,取而代之的是更先进的处理器) 。

这些技术使得指令用掉的时间不再是简单地将单独执行它们所用时间之和,而是很大程度上取决于关键指令的相互依靠性。

我们来看(A + B) * (C + D) - AC - BD,其中A+B和C+D是可以并行的,但是乘法和后面的减法由于需要依靠A+B和C+D的运行结果,不能并行执行,降低效率。

我们看A*C+B*D用VC编译后的汇编代码(DEBUG:0040105C mov edx,dword ptr [ebp+0Ch]0040105F imul edx,dword ptr [ebp+10h]00401063 mov eax,dword ptr [ebp+8]00401066 imul eax,dword ptr [ebp+14h]0040106A add edx,eax其中而(A + B) * (C + D) -AC - BD:004010EC mov edx,dword ptr [ebp+8]004010EF add edx,dword ptr [ebp+0Ch]004010F2 mov eax,dword ptr [ebp+10h]004010F5 add eax,dword ptr [ebp+14h]004010F8 imul edx,eax004010FB sub edx,dword ptr [ebp-4]004010FE sub edx,dword ptr [ebp-8]00401101 mov dword ptr [ebp-0Ch],edx其实从本质上而言,我们为什么要如此麻烦地处理两个个比DWOR大一倍的数据结构,不就不会溢出了么?我们想到了unsigned __int64 ,它是VC支持的一种DWOR的相乘?溢出问题。

相关主题