使用向量容器存储 32位int值,每一段int存储8位数字( 按10^8 为进制处理)
减法运算效率测试(分为 vec_minus 和 s_minus 函数)
// bignum minus
// 523066680/vicyang
// 2019-01
#include <iostream>
#include <chrono>
#include <vector>
#include <cstdio>
using namespace std;
vector<int> vec_minus(const vector<int> &a, const vector<int> &b);
string s_minus(const string &a, const string &b);
int s_cmp( const string &a, const string &b );
string vec2str( const vector<int> &vec );
void check(const vector<int> &va, const vector<int> &vb)
{
string a = vec2str(va);
string b = vec2str(vb);
string cmd = "perl -Mbignum -e \"print " + a + "-" + b + "\"";
system( cmd.c_str() );
cout << " <- check " << endl;
}
int main(int argc, char *argv[] )
{
auto stage0 = chrono::system_clock::now();
chrono::duration<double> diff;
string a(10000, '8');
string b(10000, '9');
string c;
// s_minus 耗时测试 10000位数,2W次
for (int i = 0; i < 20000; i++) s_minus(a, b);
auto stage1 = chrono::system_clock::now();
diff = stage1-stage0;
cout << " s_minus, Time used: " << diff.count() << endl;
// 小范围的"向量"减法测试
vector<int> va{552, 443, 12345678, 87654321};
vector<int> vb{ 93456, 92432100};
//check(va, vb);
vector<int> vc = vec_minus(va,vb);
cout << vec2str(vc) << endl;
// 结果应为1的测试
va = {1, 0, 0};
vb = {99999999,99999999};
//check(va, vb);
vc = vec_minus(va, vb);
cout << vec2str(vc) << endl;
// vec_minus 耗时测试 10000位数,2W次
stage1 = chrono::system_clock::now();
va.assign( 1250, 99999999 );
vb.assign( 1250, 12345678 );
for (int i = 0; i < 20000; i++) vec_minus(va, vb);
auto stage2 = chrono::system_clock::now();
diff = stage2-stage1;
cout << "vec_minus, Time used: " << diff.count() << endl;
return 0;
}
// 测试vector方案
vector<int> vec_minus(const vector<int> &a, const vector<int> &b)
{
static int ia; // iter
const int base = 100000000;
register const int *pa = a.data();
register const int *pb = b.data();
vector<int> c( a.size() );
register int *pc = c.data();
int t, cut=0, ib=b.size()-1, zero=0;
for (ia = a.size()-1; ia >= 0; ia-- )
{
t = ib >= 0 ? (pa[ia]) - (pb[ib--]) + cut
: (pa[ia]) + cut;
t < 0 ? t += base, cut = -1 : cut = 0;
zero = t == 0 ? zero+1 : 0; // 此判断须独立,t有可能+base后才为0
pc[ia] = t;
}
c.erase(c.begin(), c.begin()+zero);
return c;
}
// 大数减法 字符串操作, 暂时假设 a >= b
string s_minus(const string &a, const string &b)
{
static int ia;
// 获取指针以绕过[]重载,更快地取址
register const char* pa=a.data();
register const char* pb=b.data();
int cmp = s_cmp(a, b);
if (cmp == 0) return "0";
string s( a.length(), '0');
register const char* ps=s.data();
int t, cut=0, ib=b.length()-1, zero=0;
for (ia = a.length()-1; ia >= 0; ia-- )
{
t = ib >= 0 ? (pa[ia]-'0') - (pb[ib--]-'0') - cut
: (pa[ia]-'0') - cut;
cut = t < 0 ? 1 : 0;
s[ia] = t + '0' + cut*10;
ps[ia] == '0' ? zero++ : zero=0;
}
if (zero > 0) s.erase(0, zero);
return s;
}
string vec2str( const vector<int> &vec )
{
const int base = 100000000;
string s("");
s += to_string( vec[0] );
for ( int it = 1; it < vec.size(); it++ )
s += to_string(vec[it]+base).substr(1,8);
return s;
}
int s_cmp(const string &a, const string &b )
{
if ( a.length() > b.length() ) return 1;
else if ( a.length() < b.length() ) return -1;
else return a.compare(b);
}
1W位数字减法,2W次循环,向量(10^8 进制) 和 string 逐位操作,效率对比
代码: 全选
s_minus, Time used: 0.860001
vec_minus, Time used: 0.12
《Modern Computer Arithmetic》(后面简称MCA)里面关于采用高进制计算可能产生溢出问题的探讨
对于64位平台 unsigned long long int 支持的最大数字是 2^64-1=18446744073709551615,20位,如果我们充分利用,使用2^64,或者10^19作为进制。
很容易会遇到溢出问题,MCA给出了几种方案:
Let T be the number of different values taken by the data type representing the coefficients ai, bi.
(Clearly, β ≤ T, but equality does not necessarily hold, for example β = 10^9 and T = 2^32.) At step 3, the value of s can be as large as 2β − 1, which is not representable if β = T. Several workarounds are possible:
either use a machine instruction that gives the possible carry of ai + bi,
or use the fact that, if a carry occurs in ai + bi, then the computed sum – if performed modulo T – equals t := ai +bi −T < ai; thus, comparing t and ai will determine if a carry occurred.
A
third solution is to keep a bit in reserve, taking β ≤ T/2.
用 T 表示一个存储单元所能表示的数的量,则有 β <= T
(这里 β 表示采用的进制,以及β不一定等于T,例如T=2^32,但采用的进制为10^9)。
考虑加法操作 s=a+b+d (d为1或0,是上一次加法补进的数值),s的最大可能值为 2*β-1,当 β = T 时该公式无法正确计算。
考虑以下方案:
1. 内部编码实现 (水平有限,暂时这么翻译)
2. 通过-T取余数判断,t := a + b - T < a ,对比 t 和 a 断定是否进 1
3. 采用一个保守的进制数 β,令 β ≤ T/2.
偷个懒,暂时采用10^8而不是 10^9, 2^32
如果对相同的进制进行乘法运算,可以使用 unsigned long long int 类型获得最高19位数的运算空间,足以应付 8*2+1 位数的情况