1、背景
工作中遇到过需要进行极大数据的存储和运算的场景,当时使用Python解决了这个问题,在Python中,整数没有位数限制,使用起来很方便。但是当程序主体使用C/C++实现时,就比较麻烦。所以考虑实现一个大数类,用于大数的存储和运算,后面生成静态库,需要的时候直接调用。
2、算法设计
(1)存储
大数一般都是以整数队列的形式存储,取一个较大的进制S,队列中的每一个值,相当于一位。假设数组长度为N(0为低位,n-1为高位),则数组表示的值 value = A0 + A1*S + A2*S2 + ... + An-1*Sn-1。例如 “1234567890000” 在数组中存储形式为 567890000, 1234。以下大数的运算全部基于此方式。
(2)比较运算
两个大数比较的方法是:优先比较位数,位数大者为大,位数小者为小;如果位数相等,则依次从高位向低位比较,当前位大的为大,当前位小的为小;如果每一位都相等,则两者相等。
对于A(a0a1a2...an-1) 和 B(b0b1b2...bm-1),比较流程图如下:
(3)加法运算
大数的加减乘运算都与十进制运算的手算方法相同,区别在于后者是10进制,而前者是S进制。加法计算时,A和B由低位到高位依次求和,需考虑进位。
对于A(a0a1a2...an-1) + B(b0b1b2...bm-1),并假设N≥M,则计算过程如下(注意取值范围):
注意:
- A、B位数可能不相等,因此需要分两部分进行加和;
- 两次循环结束,需要考虑进位是否为零;
(4)减法运算
对于A(a0a1a2...an-1) - B(b0b1b2...bm-1),并假设A≥B(不考虑负数),一定有N≥M,则计算过程如下:
注意:
- A、B位数可能不相等,因此需要分两部分求差值;
- 减法运算可能会导致差值的若干个高位是无意义的零值,因此需额外的去零操作;
(5)乘法运算
在计算十进制数 123 * 45 时,首先计算 123 * 5 = 615,然后计算 123 * 4 = 492,最终结果是 615 + 492 * 10 = 5535。
大数乘法的计算过程与此相同:对于A(a0a1a2...an-1) * B(b0b1b2...bm-1),并假设N≥M(被乘数位数不小于乘数时,效率最高),对于B中的每一位 bi,计算其与A的乘积,并在后面补上 i 个零,将结果累加,即为A与B的乘积。
流程图如下:
(6)除法运算
除法运算是四则运算中最复杂的运算, 常用的方法是使用多次减法来模拟除法。最简单的方法是被除数不断减去除数,直到结果为负,减的次数即为商,当然这种方法效率太低。
网上提供了一些优化方法,例如,在被除数减除数的过程中,被除数每次减去除数的10N倍,以加速被除数衰减过程。
这里采用的方法与上述方法原理相似,以23456除以123为例,在进行竖式计算时:以234除以123,得1,余111,以1115除以123,得9,余8,以86除以123,得0,余86,则23456除以123的结果是190。
上述过程是为了将未知位数的两个大数相除,转换为多次位数相近的两个数相除,即被除数最多比除数多1位。算法的核心是上述过程中的 234除以124、1115除以123、86除以123,如何得到商和余数。
以下提出一种收敛算法:
首先,对于 A(a0a1a2...an-1) / B(b0b1b2...bm-1),一定有 A(a0a1a2...an-1) / B(b0b1b2...bm-1) ≤ A(aiai+1...an-1) / B(bibi+1...bm-1),其中 0 ≤ i ≤ n-1 且 0 ≤ i ≤ m-1。
所以,可以通过取A、B的最高两位或一位来预估A、B相除的结果,且预估值≥实际值,并以预估值更新A,多次迭代,直到A<B。设A0=A,A/B的最终结果为C0,推导过程如下:
第1轮:C0 = A0 / B,求得预估值 V0,有V0 >= C0,设V0 - C0 = C1,B * V0 - A0 = A1,则
第2轮:C1 = A1 / B,求得预估值 V1,有V1 >= C1,设V1 - C1 = C2,B * V1 - A1 = A2,则
第3轮:C2 = A2 / B,求得预估值 V2,有V2 >= C2,设V2 - C2 = C3,B * V2 - A2 = A3,则
......
第x轮:Cx-1 = Ax-1 / B,此时Ax-1 < B,即 Vx-1 = 0
通过以上过程,可以推导出:C0 = V0 - V1 + V2 - ... +/- Vx-1
举例(进制S取10亿):
A0 = 86,517999162,161442275,630671648,031880106,681829550,207222443
B = 92,784489371,679693896,011626721,067864399,494212548
第1轮:V0 = 86,517999162 / 92 = 940413034
A1 = B * V0 - A0 = 737743996,000612337,308902736,003021973,147110937,498328189
第2轮:V1 = 737743996,000612337 / 92,784489371 = 7951156
A2 = B * V1 - A1 = 0(实际值是个负值,直接取0)
此时 A2 < B,因此 C = 940413034 - 7951156 = 932461878,实际值为932461877。
因为误差问题,算法值比实际值大1或与实际值相等,所以算法的最终结果需要向下微调。
以下给出收敛过程以及预估值的计算方法:
3、代码实现
(1)类的设计
计划实现一个 LargeInt 类,其含义是一个大整数(无负数和小数),实现的核心功能是:字符串构造、格式化字符串输出、加减乘除四则运算以及逻辑比较运算。
设计类图如下:
(2)类定义
以C++模板类vector作为大数存储,下标0表示低位,进制S取1000000000(10亿,便于格式化输入输出),大数每一位是一个无符号32位整型,乘法运算时的中间值以无符号64位整型承载。
1 #include <iostream> 2 #include <vector> 3 #include <string> 4 5 typedef unsigned int uint; 6 typedef unsigned long long uint64; 7 8 #define MAX_VAL 1000000000 // 10亿 9 #define VAL_LEN 9 10 #define FORMAT_STR "%09d" 11 12 class LargeInt 13 { 14 private: 15 std::vector<uint> _data; 16 17 bool checkValStr(const std::string &valStr); 18 int compare(const LargeInt &li) const; 19 void arrange(); 20 21 static uint estimateQuotient(const LargeInt &liA, const LargeInt &liB); 22 static uint getMaxCycle(const LargeInt &liA, const LargeInt &liB); 23 24 public: 25 LargeInt(); 26 LargeInt(uint val); 27 LargeInt(const std::string &valStr); 28 29 // 四则运算符重载 30 LargeInt operator+(const LargeInt &li) const; 31 LargeInt operator-(const LargeInt &li) const; 32 LargeInt operator*(const LargeInt &li) const; 33 LargeInt operator/(const LargeInt &li) const; 34 35 // 比较运算符重载 36 bool operator==(const LargeInt &li) const; 37 bool operator!=(const LargeInt &li) const; 38 bool operator<(const LargeInt &li) const; 39 bool operator>(const LargeInt &li) const; 40 bool operator<=(const LargeInt &li) const; 41 bool operator>=(const LargeInt &li) const; 42 43 // 字符串格式化输出 44 std::string toString() const; 45 };
行9、10定义了格式化输入和格式化输出所需的宏,后面会用到。
(3)构造方法
LargeInt 提供了三种构造方法,分别是无符号整型、字符串以及无参数构造。所谓字符串构造,是以十进制字符串为入参,系统对字符串进行分割,转为整数数组,构造LargeInt。
以 ”1234567890123“ 为例,字符串构造的方法是,从字符串结尾向开头不断截取长度为9的子串,并转为整数,即 67890123,12345。
1 inline bool isDigit(const char ch) 2 { 3 return ch >= '0' && ch <= '9'; 4 } 5 6 LargeInt::LargeInt() {} 7 8 LargeInt::LargeInt(uint val) 9 { 10 this->_data.push_back(val % MAX_VAL); 11 if (val >= MAX_VAL) 12 this->_data.push_back(val / MAX_VAL); 13 } 14 15 // 字符串合法性检查,只允许出现字符 0~9 16 bool LargeInt::checkValStr(const string &valStr) 17 { 18 for (auto it = valStr.begin(); it != valStr.end(); ++it) 19 { 20 if (!isDigit(*it)) return false; 21 } 22 return true; 23 } 24 25 // 字符串构造 26 // valStr:十进制字符串 27 LargeInt::LargeInt(const string &valStr) 28 { 29 if (checkValStr(valStr))// 检查valStr合法性 30 { 31 int len = valStr.length(); 32 // 按长度9截取子串 33 while (len >= VAL_LEN) 34 { 35 string s = valStr.substr(len - VAL_LEN, VAL_LEN); 36 this->_data.push_back(stoi(s)); 37 len -= VAL_LEN; 38 } 39 // 残留子串 40 if (len > 0) 41 { 42 string s = valStr.substr(0, len); 43 this->_data.push_back(stoi(s)); 44 } 45 } 46 47 this->arrange(); // 去零 48 } 49 50 // 去零操作,避免整数队列的高位存在多余的零 51 void LargeInt::arrange() 52 { 53 int idx = this->_data.size(); 54 55 // 注意,如果队列中全为0,要保留最低位的0 56 while (--idx >= 1) 57 { 58 if (this->_data[idx] > 0) break; 59 60 this->_data.pop_back(); 61 } 62 }
说明:为避免字符串开头有无意义的零,定义了 ”arrange“ 函数,用于整理队列,去除多余的零值,字符串构造以及后面的减法、除法都需要调用arrange函数。
(4)比较运算
对于所有的比较运算,可以提取出一个公共函数 ”compare“,该函数比较两个LargeInt的大小,返回1、0、-1,其他比较函数通过调用该函数判断返回值即可。
1 // 比较函数,0 相等,1 大于,-1 小于 2 int LargeInt::compare(const LargeInt &li) const 3 { 4 int len1 = this->_data.size(); 5 int len2 = li._data.size(); 6 7 // step1: 比较长度 8 if (len1 != len2) 9 return (len1 > len2) ? 1 : -1; 10 11 // step2: 由高位到低位比较值 12 for (int idx = len1 - 1; idx >= 0; --idx) 13 { 14 if (this->_data[idx] == li._data[idx]) continue; 15 return this->_data[idx] > li._data[idx] ? 1 : -1; 16 } 17 18 return 0; 19 } 20 21 bool LargeInt::operator==(const LargeInt &li) const 22 { 23 return compare(li) == 0; 24 } 25 26 bool LargeInt::operator!=(const LargeInt &li) const 27 { 28 return compare(li) != 0; 29 } 30 31 bool LargeInt::operator<(const LargeInt &li) const 32 { 33 return compare(li) < 0; 34 } 35 36 bool LargeInt::operator>(const LargeInt &li) const 37 { 38 return compare(li) > 0; 39 } 40 41 bool LargeInt::operator<=(const LargeInt &li) const 42 { 43 return compare(li) <= 0; 44 } 45 46 bool LargeInt::operator>=(const LargeInt &li) const 47 { 48 return compare(li) >= 0; 49 }
(5)加法运算
首先定义两个简单的MAX/MIN宏,后面会用到。
1 #ifndef MIN 2 #define MIN(a,b) ((a) > (b) ? (b) : (a)) 3 #endif 4 5 #ifndef MAX 6 #define MAX(a,b) ((a) > (b) ? (a) : (b)) 7 #endif
按照前面提供的算法,实现代码。考虑加法的进位只会是0或1,所以为了减少除法和求模运算,函数中使用 if 判断替代。
1 LargeInt LargeInt::operator+(const LargeInt &li) const 2 { 3 int len1 = this->_data.size(); 4 int len2 = li._data.size(); 5 6 int minLen = MIN(len1, len2); 7 int maxLen = MAX(len1, len2); 8 const LargeInt &extraLi = (len1 > len2) ? (*this) : li; 9 10 uint value = 0; // 和,不超过最大值的2倍 11 uint carry = 0; // 进位 12 LargeInt retVal; 13 14 for (int idx = 0; idx < minLen; ++idx) 15 { 16 value = this->_data[idx] + li._data[idx] + carry; 17 18 if (value >= MAX_VAL) 19 { 20 retVal._data.push_back(value - MAX_VAL); 21 carry = 1; 22 } 23 else 24 { 25 retVal._data.push_back(value); 26 carry = 0; 27 } 28 } 29 30 for (int idx = minLen; idx < maxLen; ++idx) 31 { 32 value = extraLi._data[idx] + carry; 33 34 if (value >= MAX_VAL) 35 { 36 retVal._data.push_back(value - MAX_VAL); 37 carry = 1; 38 } 39 else 40 { 41 retVal._data.push_back(value); 42 carry = 0; 43 } 44 } 45 46 if (carry > 0) 47 retVal._data.push_back(carry); 48 49 return retVal; 50 }
(6)减法运算
注意:代码中要避免出现负值,否则会导致计算结果错误;减法运算的结果需要做去零操作。
1 LargeInt LargeInt::operator-(const LargeInt &li) const 2 { 3 if (*this <= li) 4 { 5 return LargeInt(0); 6 } 7 8 int len1 = this->_data.size(); 9 int len2 = li._data.size(); 10 11 uint value = 0; // 差 12 uint carry = 0; // 借位 13 LargeInt retVal; 14 15 for (int idx = 0; idx < len2; ++idx) 16 { 17 if (this->_data[idx] < li._data[idx] + carry) // 注意细节,carry放在右侧,避免出现差值为负数的情况 18 { 19 value = this->_data[idx] + MAX_VAL - carry - li._data[idx]; 20 carry = 1; 21 } 22 else 23 { 24 value = this->_data[idx] - carry - li._data[idx]; 25 carry = 0; 26 } 27 28 retVal._data.push_back(value); 29 } 30 31 for (int idx = len2; idx < len1; ++idx) 32 { 33 if (this->_data[idx] < carry) 34 { 35 value = this->_data[idx] + MAX_VAL - carry; 36 carry = 1; 37 } 38 else 39 { 40 value = this->_data[idx] - carry; 41 carry = 0; 42 } 43 retVal._data.push_back(value); 44 } 45 46 retVal.arrange(); 47 return retVal; 48 }
(7)乘法运算
为保证被乘数位数大于乘数,提高计算效率,如果A的位数小于B的位数,则返回B*A。
1 LargeInt LargeInt::operator*(const LargeInt &li) const 2 { 3 int len1 = this->_data.size(); 4 int len2 = li._data.size(); 5 6 if (len1 < len2) return li.operator*(*this); // 优化,保证被乘数位数大于乘数 7 8 uint64 value; // 积 9 uint64 carry = 0; // 进位 10 LargeInt retVal(0); 11 12 for (int idx2 = 0; idx2 < len2; ++idx2) 13 { 14 LargeInt mulTemp; 15 carry = 0; 16 17 // 补零 18 for (int tmpIdx = 0; tmpIdx < idx2; ++tmpIdx) 19 mulTemp._data.push_back(0); 20 21 for (int idx1 = 0; idx1 < len1; ++idx1) 22 { 23 value = (uint64)(li._data[idx2]) * (uint64)(this->_data[idx1]) + carry; 24 25 mulTemp._data.push_back((uint)(value % MAX_VAL)); 26 carry = value / MAX_VAL; 27 } 28 29 if (carry) 30 mulTemp._data.push_back((uint)carry); 31 32 retVal = retVal + mulTemp; 33 } 34 35 return retVal; 36 }
(8)除法
注意:除法计算由高位向低位进行,即求得的第一个商值是结果的最高位,所以要注意插入位置;
算法结果需要向下微调,正常情况下,算法结果与实际结果相等,或比实际结果大1,为保险起见,使用了while循环。
1 LargeInt LargeInt::operator/(const LargeInt &li) const 2 { 3 if (li._data.empty() || li == 0) return LargeInt(""); 4 if (*this < li) return LargeInt(0); 5 6 int len1 = this->_data.size(); 7 int len2 = li._data.size(); 8 9 uint value; 10 LargeInt retVal; 11 LargeInt divTemp; 12 13 for (int idx = len1 - len2 + 1; idx < len1; ++idx) 14 { 15 divTemp._data.push_back(this->_data[idx]); 16 } 17 18 // len1 >= len2 19 for (int idx = len1 - len2; idx >= 0; --idx) 20 { 21 value = 0; 22 divTemp._data.insert(divTemp._data.begin(), this->_data[idx]); 23 divTemp.arrange(); 24 25 value = getMaxCycle(divTemp, li); // 商 26 27 divTemp = divTemp - li * value; // 余数 28 29 retVal._data.insert(retVal._data.begin(), value); // 除法是由高位向低位进行,所以插入位置在begin 30 } 31 32 retVal.arrange(); 33 return retVal; 34 } 35 36 // 计算商值 37 uint LargeInt::getMaxCycle(const LargeInt &liA, const LargeInt &liB) 38 { 39 LargeInt tempA = liA; 40 const LargeInt& tempB = liB; 41 uint tempC; 42 uint res = 0; 43 uint counter = 0; // 调试用 44 bool flag = true; 45 46 47 while (tempA >= tempB) 48 { 49 counter++; 50 51 tempC = estimateQuotient(tempA, tempB); 52 tempA = tempB * tempC - tempA; 53 54 res = flag ? (res + tempC) : (res - tempC); 55 flag = !flag; 56 } 57 58 // 微调 59 while (res > 0 && liB * res > liA) res--; 60 61 return res; 62 } 63 64 // 估值 65 uint LargeInt::estimateQuotient(const LargeInt &liA, const LargeInt &liB) 66 { 67 int lenA = liA._data.size(); 68 int lenB = liB._data.size(); 69 uint64 valA, valB; 70 71 if (lenA == lenB) 72 { 73 if (lenA > 1) 74 { 75 valA = (uint64)liA._data[lenA - 1] * MAX_VAL + liA._data[lenA - 2]; 76 valB = (uint64)liB._data[lenB - 1] * MAX_VAL + liB._data[lenB - 2]; 77 } 78 else 79 { 80 valA = (uint64)liA._data[lenA - 1]; 81 valB = (uint64)liB._data[lenB - 1]; 82 } 83 } 84 else 85 { 86 valA = (uint64)liA._data[lenA - 1] * MAX_VAL + liA._data[lenA - 2]; 87 valB = (uint64)liB._data[lenB - 1]; 88 } 89 90 return (uint)(valA / valB); 91 }
(9)格式化字符串输出
为便于测试和打印,需要得到 LargeInt 的字符串表现形式,类似于十进制字符串构造,这里将 LargeInt 转为十进制字符串。
1 string LargeInt::toString() const 2 { 3 int len = this->_data.size(); 4 int shift = 0; 5 char *buff = new char[len * VAL_LEN + 1]; 6 7 if (len > 0) 8 shift += sprintf(buff + shift, "%d", this->_data[len - 1]); 9 10 for (int idx = len - 2; idx >= 0; --idx) 11 { 12 shift += sprintf(buff + shift, FORMAT_STR, this->_data[idx]); 13 } 14 buff[shift] = '