【转】高精度(压位存储)
有的时候,数字会大到连long long都不能承受的程度。这时,我们可以自己模拟大数的各种运算。
所谓压位存储,就是在高精度数内部采用10000进制(即每四位放到一个数中)进行存储。它与10进制(即一个数位对应一个数)相比速度要快一些。
高精度数内部也可以采用100000000进制,但是这样就不能计算乘除法了。
(1) 定义
编程时这样做——假设hp是高精度类型。
先用宏定义:#define hp long long,然后集中精力编写算法代码。
最后直接删除这条宏定义,把真正的高精度算法写出来。这样做的好处是无需再修改算法代码,减小了维护代价。
1 const int MAX=100; 2 struct hp 3 { 4 int num[MAX]; 5 6 hp & operator = (const char*); 7 hp & operator = (int); 8 hp(); 9 hp(int); 10 11 // 以下运算符可以根据实际需要来选择。 12 bool operator > (const hp &) const; 13 bool operator < (const hp &) const; 14 bool operator <= (const hp &) const; 15 bool operator >= (const hp &) const; 16 bool operator != (const hp &) const; 17 bool operator == (const hp &) const; 18 19 hp operator + (const hp &) const; 20 hp operator - (const hp &) const; 21 hp operator * (const hp &) const; 22 hp operator / (const hp &) const; 23 hp operator % (const hp &) const; 24 25 hp & operator += (const hp &); 26 hp & operator -= (const hp &); 27 hp & operator *= (const hp &); 28 hp & operator /= (const hp &); 29 hp & operator %= (const hp &); 30 };
不使用宏定义,也可以用typedef long long hp;
实现运算符的代码最好写到结构体的外面。
“<<”和“>>”由于不是hp的成员,所以必须写到结构体的外面。
(2) 赋值和初始化
1 // num[0]用来保存数字位数。利用10000进制可以节省空间和时间。 2 hp & hp::operator = (const char* c) 3 { 4 memset(num,0,sizeof(num)); 5 int n=strlen(c),j=1,k=1; 6 for (int i=1;i<=n;i++) 7 { 8 if (k==10000) j++,k=1;// 10000进制,4个数字才算1位。 9 num[j]+=k*(c[n-i]-'0'); 10 k*=10; 11 } 12 num[0]=j; 13 return *this; 14 } 15 16 hp & hp::operator = (int a) 17 { 18 char s[MAX]; 19 sprintf(s,"%d",a); 20 return *this=s; 21 } 22 23 hp::hp() {memset(num,0,sizeof(num)); num[0]=1;} 24 hp::hp (int n) {*this = n;}// 目的:支持“hp a=1;”之类的代码。
(3) 比较运算符
小学时候学过怎么比较两个数的大小吧?现在,虽为高中生,小学知识却从未过时……
1 // 如果位数不等,大小是可以明显看出来的。如果位数相等,就需要逐位比较。 2 bool hp::operator > (const hp &b) const 3 { 4 if (num[0]!=b.num[0]) return num[0]>b.num[0]; 5 for (int i=num[0];i>=1;i--) 6 if (num[i]!=b.num[i]) 7 return (num[i]>b.num[i]); 8 return false; 9 } 10 bool hp::operator < (const hp &b) const {return b>*this;} 11 bool hp::operator <= (const hp &b) const {return !(*this>b);} 12 bool hp::operator >= (const hp &b) const {return !(b>*this);} 13 bool hp::operator != (const hp &b) const {return (b>*this)||(*this>b);} 14 bool hp::operator == (const hp &b) const {return !(b>*this)&&!(*this>b);}
(4) 四则运算
如果没学过竖式,或者忘了怎么用竖式算数,那么你就悲剧了……
1. 加法和减法
1 // 注意:最高位的位置和位数要匹配。 2 hp hp::operator + (const hp &b) const 3 { 4 hp c; 5 c.num[0] = max(num[0], b.num[0]); 6 for (int i=1;i<=c.num[0];i++) 7 { 8 c.num[i]+=num[i]+b.num[i]; 9 if (c.num[i]>=10000) // 进位 10 { 11 c.num[i]-=10000; 12 c.num[i+1]++; 13 } 14 } 15 if (c.num[c.num[0]+1]>0) c.num[0]++; // 9999+1,计算完成后多了一位 16 17 return c; 18 } 19 20 hp hp::operator - (const hp &b) const 21 { 22 hp c; 23 c.num[0] = num[0]; 24 for (int i=1;i<=c.num[0];i++) 25 { 26 c.num[i]+=num[i]-b.num[i]; 27 if (c.num[i]<0) // 退位 28 { 29 c.num[i]+=10000; 30 c.num[i+1]--; 31 } 32 } 33 while (c.num[c.num[0]]==0&&c.num[0]>1) c.num[0]--; // 100000000-99999999 34 35 return c; 36 } 37 38 hp & hp::operator += (const hp &b) {return *this=*this+b;} 39 hp & hp::operator -= (const hp &b) {return *this=*this-b;}
2. 乘法
1 hp hp::operator * (const hp &b) const 2 { 3 hp c; 4 c.num[0] = num[0]+b.num[0]+1; 5 for (int i=1;i<=num[0];i++) 6 { 7 for (int j=1;j<=b.num[0];j++) 8 { 9 c.num[i+j-1]+=num[i]*b.num[j]; // 和小学竖式的算法一模一样 10 c.num[i+j]+=c.num[i+j-1]/10000; // 进位 11 c.num[i+j-1]%=10000; 12 } 13 } 14 while (c.num[c.num[0]]==0&&c.num[0]>1) c.num[0]--; // 99999999*0 15 16 return c; 17 } 18 19 hp & hp::operator *= (const hp &b) {return *this=*this*b;}
3.除法
高精度除法的使用频率不太高。
以下代码的缺陷是:如果b太大,运算速度会非常慢。该问题可以用二分查找解决。
1 hp hp::operator / (const hp &b) const 2 { 3 hp c, d; 4 c.num[0] = num[0]+b.num[0]+1; 5 d.num[0] = 0; 6 for (int i=num[0];i>=1;i--) 7 { 8 // 以下三行的含义是:d=d*10000+num[i]; 9 memmove(d.num+2, d.num+1, sizeof(d.num)-sizeof(int)*2); 10 d.num[0]++; 11 d.num[1]=num[i]; 12 13 // 以下循环的含义是:c.num[i]=d/b; d%=b; 14 while (d >= b) 15 { 16 d-=b; 17 c.num[i]++; 18 } 19 } 20 while (c.num[c.num[0]]==0&&c.num[0]>1) c.num[0]--; // 99999999/99999999 21 22 return c; 23 } 24 hp hp::operator % (const hp &b) const 25 { 26 …… // 和除法的代码一样。唯一不同的地方是返回值:return d; 27 } 28 29 hp & hp::operator /= (const hp &b) {return *this=*this/b;} 30 hp & hp::operator %= (const hp &b) {return *this=*this%b;}
4. 二分优化的除法
高精度除法速度慢,就慢在上面的while (d>=b)处。如果我们用二分法去猜d/b的值,速度就快了。
1 hp hp::operator / (const hp& b) const 2 { 3 hp c, d; 4 c.num[0] = num[0]+b.num[0]+1; 5 d.num[0] = 0; 6 for (int i=num[0];i>=1;i--) 7 { 8 // 以下三行的含义是:d=d*10000+num[i]; 9 memmove(d.num+2, d.num+1, sizeof(d.num)-sizeof(int)*2); 10 d.num[0]++; 11 d.num[1]=num[i]; 12 13 // 以下循环的含义是:c.num[i]=d/b; d%=b; 利用二分查找求c.num[i]的上界。 14 // 注意,这里是二分优化后除法和朴素除法的区别! 15 int left=0, right=9999, mid; 16 while (left < right) 17 { 18 mid = (left+right)/2; 19 if (b*hp(mid) <= d) left=mid+1; 20 else right=mid; 21 } 22 c.num[i]=right-1; 23 d=d-b*hp(right-1); 24 } 25 while (c.num[c.num[0]]==0&&c.num[0]>1) c.num[0]--; // 99999999/99999999 26 27 return c; // 求余数就改成return d; 28 }
(5) 输入/输出
有了这两段代码,就可以直接用cout和cin输出、输入高精度数了。
1 ostream & operator << (ostream & o, hp &n) 2 { 3 o<<n.num[n.num[0]]; 4 for (int i=n.num[0]-1;i>=1;i--) 5 { 6 o.width(4); 7 o.fill('0'); 8 o<<n.num[i]; 9 } 10 return o; 11 } 12 13 istream & operator >> (istream & in, hp &n) 14 { 15 char s[MAX]; 16 in>>s; 17 n=s; 18 return in; 19 }