我写的另一个大数类,基于二进制(确切的说在32位机下是2^31次进制),乘法是用移位和加法模拟的(如果不是考虑到64位兼容性的话,用long long可以显著提高乘法效率),除法是移位和减法(状况同乘法),效率上比我原来的那个低了很多(有兴趣的可以看我之前的文章,基于10^9次进制的大数类),尤其是输出函数。但毕竟是是我精心设计的,很多方面感觉还是很有价值的
有需要的朋友可以下载我的源码(基于code::blocks的工程)
下面是其中的一个头文件,包括了大数类的主要声明
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
| #ifndef LARGEINT_H_INCLUDED
#define LARGEINT_H_INCLUDED
#include "numarray.h"
const int bits=8*sizeof(numunit)-1; //每一个存储单位的2进制位数
const int intlim=1000000000;//numunit能够容纳的最大的10的幂
const int intlimbits=9;//上面那个intLim的十进制位数(幂)
class largeInt
{
private:
numArray number; //反序保存数剧本身,最低位是number[0]
int length; //number的长度
int sign; //符号1表示正,-1表示负数,特殊的,有两个0,保证0的符号是1
public:
largeInt() {length=1;sign=1;}
largeInt(int in) {length=1;sign=(in>=0)?1:-1;number[0]=(numunit)sign*in;}
largeInt(const largeInt &in) {number=in.number;length=in.length;sign=in.sign;}
largeInt& operator=(const largeInt &in) {number=in.number;length=in.length;sign=in.sign;return *this;}
largeInt operator-() const {largeInt opposite=(*this);opposite.sign=-sign;return opposite;}
void opposite() {sign*=-1;}
//##I/O################################
largeInt& read();
largeInt& read(const char *);
void print() const;
private:
void printAll(largeInt &num) const;
void printDecimalUnit(int unit,int len) const;
void printDecimalUnit(int unit) const;
public:
//##compare############################
bool operator==(const largeInt&) const;
bool operator>(const largeInt&) const;
bool operator<(const largeInt&) const;
bool operator>=(const largeInt &right) const {return !(*this<right);}
bool operator<=(const largeInt &right) const {return !(*this>right);}
bool operator!=(const largeInt &right) const {return !(*this==right);}
//##movement###########################
largeInt operator<<(int n) const {largeInt solution=*this;solution<<=n;return solution;}
largeInt operator>>(int n) const {largeInt solution=*this;solution>>=n;return solution;}
const largeInt& operator<<=(int);
const largeInt& operator>>=(int);
private:
void moveLeft(int); //大移位,移1就是移动bits位
void moveRight(int);
public:
//##+ - * / %##########################
largeInt operator++(int) {*this+=1;return *this-1;}//后缀运算
const largeInt& operator++() {return *this+=1;}//前缀
largeInt operator--(int) {*this-=1;return *this+1;}
const largeInt& operator--() {return *this-=1;}
//*************************************
largeInt operator+(const largeInt &right) const {largeInt solution=*this;solution+=right;return solution;}
const largeInt& operator+=(const largeInt&);
largeInt operator-(const largeInt &right) const {largeInt solution=*this;solution-=right;return solution;}
const largeInt& operator-=(const largeInt&);
largeInt operator*(const largeInt&) const;
largeInt operator*(int) const;
const largeInt operator*=(const largeInt &right) {return *this=*this*right;}
const largeInt operator*=(int right) {return *this=*this*right;}
largeInt operator/(largeInt&) const;
largeInt operator/(int) const;
largeInt operator/(const largeInt &divi) const {largeInt divisor=divi;return *this/divisor;}
const largeInt& operator/=(largeInt &divisor) {return *this=*this/divisor;}
const largeInt& operator/=(const largeInt &divi) {largeInt divisor=divi;return *this/=divisor;}
const largeInt& operator/=(int divisor) {return *this=*this/divisor;}
largeInt operator%(const largeInt &divisor) const {largeInt solution;solution%=divisor;return solution;}
largeInt operator%(largeInt &divisor) const {largeInt solution;solution%=divisor;return solution;}
const largeInt& operator%=(largeInt&);
const largeInt& operator%=(const largeInt &divisor) {largeInt divi=divisor;return *this%=divi;}
const largeInt& operator%=(int);
int operator%(int) const; //特殊的取模函数,用于返回int型的余数
//**************************************
friend largeInt operator+(int a,const largeInt& b) {return b+a;}
friend largeInt operator-(int a,const largeInt& b) {return -(b-a);}
friend largeInt operator*(int a,const largeInt& b) {return b*a;}
friend largeInt operator/(int a,const largeInt& b) {largeInt x=a;return a/b;}
friend largeInt operator%(int a,const largeInt& b) {largeInt x=a;return a%b;}
private:
};
#endif // LARGEINT_H_INCLUDED |
给指针分配的一块内存,都应该是连续的。这必然会导致一个问题:当我分配的空间过大,系统将需要整理出一块有效的地方,这显然是不太合适的(事实上,系统不会支持过大的内存请求的,在我这里不能超过一百万)
我的策略是通过分配较小的单元用另一个指针数组对这些单元索引,类似下图
如果单元的大小为k,指针数组的大小为t,那么就可以处理k*t大小的空间。比如k=100,t=100,就是大小10000空间,只要我建立合适的函数operator[],就可以对这部分空间索引,实现一个模拟的数组
当然这只是双层的索引(基本的数组就是单层的),我在下面给出了一个三层的索引,能够处理lengthunit^3的空间,这其实是我为大数类设计的,用初始3*lengthunit的空间,实现一个开放数组(虽然是有上限的,但它是使用时才分配的)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
| //############numarray.h#################################
//###初始化注意:确保ptr1和ptr2指向的位置不是NULL,并且已经分配了空间######
#ifndef NUMARRAY_H_INCLUDED
#define NUMARRAY_H_INCLUDED
typedef unsigned int numunit;
const int lengthunit=1000;
class numArray
{
private:
numunit ***num;
int ptr1; //顶层的索引,已分配空间和未分配空间的分界线
int ptr2; //第二层索引,类似ptr1,该位置相对于ptr1
public:
numArray();
numArray(const numArray&);
~numArray();
numArray& operator=(const numArray&);
//###########
numunit& operator[](int n);//若空间不存在,会自动创建,并初始化为0
const numunit& operator[](int n) const;//该函数不会分配空间
void resize(int n); //可以额外维护一个长度,因为长度的定义不属于该类
private:
void resize(int,int);//resize(int)就是调用它实现的
};
#endif // NUMARRAY_H_INCLUDED
//#################numarray.cpp###################
#include "numarray.h"
#include <stdio.h>
numArray::numArray() //构造函数,负责最基本的初始化。
{ //注意这里的初始化原则,ptr1和ptr2指向的位置必须分配空间
num=new numunit**[lengthunit];//为顶层的索引数组分配空间
num[0]=new numunit*[lengthunit];//为顶层第一个位置指定空间
num[0][0]=new numunit[lengthunit];//为第二层索引的第一个位置分配空间
num[0][0][0]=0; //初始化归0,下同
for(int i=1;i<lengthunit;i++)
{
num[i]=NULL;
num[0][i]=NULL;
num[0][0][i]=0;
}
ptr1=0;
ptr2=0;
}
numArray::numArray(const numArray &in)
{
num=new numunit**[lengthunit];
for(int i=0;i<lengthunit;i++)
{
if(in.num[i]==NULL)
{
num[i]=NULL;
}
else
{
num[i]=new numunit*[lengthunit];
for(int j=0;j<lengthunit;j++)
{
if(in.num[i][j]==NULL)
{
num[i][j]=NULL;
}
else
{
num[i][j]=new numunit[lengthunit];
for(int t=0;t<lengthunit;t++)
{
num[i][j][t]=in.num[i][j][t];
}
}
}
}
}
ptr1=in.ptr1;
ptr2=in.ptr2;
}
numArray::~numArray()
{
for(int i=0;i<lengthunit;i++)
{
if(num[i]!=NULL)
{
for(int j=0;j<lengthunit;j++)
{
if(num[i][j]!=NULL)
{
delete [] num[i][j];
}
}
delete [] num[i];
}
}
delete [] num;
}
numArray& numArray::operator=(const numArray &in)
{
resize(in.ptr1,in.ptr2);
for(int i=0;i<=ptr1;i++)
{
for(int j=0;j<=ptr2;j++)
{
for(int t=0;t<lengthunit;t++)
{
num[i][j][t]=in.num[i][j][t];
}
}
}
return *this;
}
numunit& numArray::operator[](int n)
{
int p1=(n/lengthunit)/lengthunit;
int p2=(n/lengthunit)%lengthunit;
if(ptr1<p1||(ptr1==p1&&ptr2<p2))
{
resize(n);
}
return num[p1][p2][n%lengthunit];
}
const numunit& numArray::operator[](int n) const
{
int p1=(n/lengthunit)/lengthunit;
int p2=(n/lengthunit)%lengthunit;
if(num[p1][p2]==NULL)
{
printf("error:out of space!\n");
system("pause");
}
return num[p1][p2][n%lengthunit];
}
void numArray::resize(int n)
{
int p1=(n/lengthunit)/lengthunit;
int p2=(n/lengthunit)%lengthunit;
resize(p1,p2);
}
void numArray::resize(int p1,int p2)
{
if(ptr1>p1)
{
for(int i=p1+1;i<=ptr1;i++)
{
for(int j=0;j<=lengthunit;j++)
{
if(num[i][j]!=NULL)
{
delete [] num[i][j];
}
}
delete [] num[i];
num[i]=NULL;
}
}
else if(ptr1<p1)
{
for(int j=1;j<lengthunit;j++)//分配num[ptr1]的空间,注意该空间被初始化过
{
if(num[ptr1][j]==NULL)
{
num[ptr1][j]=new numunit[lengthunit];
for(int t=0;t<lengthunit;t++)
{
num[ptr1][j][t]=0;
}
}
}
for(int i=ptr1+1;i<p1;i++)//这部分空间是完全用0填充的
{
num[i]=new numunit*[lengthunit];
for(int j=0;j<lengthunit;j++)
{
num[i][j]=new numunit[lengthunit];
for(int t=0;t<lengthunit;t++)
{
num[i][j][t]=0;
}
}
}
num[p1]=new numunit*[lengthunit];//num[p2]是要初始化的,不能随便填0
ptr2=0;
num[p1][0]=new numunit[lengthunit];//确保了ptr2指向的位置分配了空间
for(int t=0;t<lengthunit;t++)
{
num[p1][0][t]=0;
}
for(int j=1;j<lengthunit;j++)
{
num[p1][j]=NULL;
}
} //到此为止,p1和ptr1应该是相等的
ptr1=p1;
if(ptr2<p2) //开始处理ptr2
{
for(int j=ptr2+1;j<=p2;j++)
{
if(num[ptr1][j]==NULL)
{
num[ptr1][j]=new numunit[lengthunit];
for(int t=0;t<lengthunit;t++)
{
num[ptr1][j][t]=0;
}
}
}
}
else if(ptr2>p2)
{
for(int j=p2+1;j<=ptr2;j++)
{
if(num[ptr1][j]!=NULL)
{
delete [] num[ptr1][j];
num[ptr1][j]=NULL;
}
}
}
ptr2=p2;
} |
为了编写一个更为有效的大数类,我设计了这个开放数组,就是不存在上限的数组(只要机器不崩溃)。
我简单说下思路:
其实数组使用链表模拟的,每个链表单元存储unitsize个元素(给一个unsigned指针分配的连续内存),然后重载operator[],像数组一样对它索引。unitsize不能太小,因为链表的索引效率不高,同样是1000000个单位,unitsize是1000的话,链表长度就是1000,unitsize是10000的话链表长度就是100;但是unitsize也不能太大,否则当运算位数较低时(比如几百位),初始化效率比较低,也就违背了我的初衷(初始化1000个元素的数组和初始化1000000个元素的数组效率是不同的,尤其是经常有临时对象被创建的时候)。
下面贴一下源代码
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
| //############openarray.h#################
#ifndef OPENARRAY_H_INCLUDED
#define OPENARRAY_H_INCLUDED
const unsigned unitsize=1000;
typedef unsigned numunit;
class openarray
{
private:
class arrayunit
{
public:
numunit *num;
arrayunit *next;
public:
arrayunit();
arrayunit(const arrayunit*);
~arrayunit();
};
arrayunit *root; //链表的根节点
public:
openarray(){root=new arrayunit;}
openarray(const openarray &in){root=new arrayunit(in.root);}
~openarray(){delete root;}
//####################################
numunit& operator[](unsigned offset);
void resize(unsigned length);
};
#endif // OPENARRAY_H_INCLUDED
//##############################################
//***************************************************
//#################openarray.cpp###################
#include “openarray.h”
#include <stdio.h>
//#########Class arrayunit##########################
openarray::arrayunit::arrayunit()
{
num=new numunit[unitsize];
for(unsigned i=0;i<unitsize;i++)
{
num[i]=0;
}
next=NULL;
}
openarray::arrayunit::arrayunit(const arrayunit *in)
{ //拷贝该节点以及后面的所有节点
num=new numunit[unitsize];
for(unsigned i=0;i<unitsize;i++)
{
num[i]=in->num[i];
}
if(in->next!=NULL)
{
next=new arrayunit(in->next);
}
else
{
next=NULL;
}
}
openarray::arrayunit::~arrayunit() //析构该节点以及后面所有节点
{
delete [] num;
if(next!=NULL)
{
delete next; //递归析构
}
}
//#########Class openarray##########################
numunit& openarray::operator[](unsigned offset)
{
arrayunit *p=root; //p记录链表中的位置
unsigned adr=offset/unitsize; //adr是p和所寻找节点的偏移
offset%=unitsize; //所寻找元素在节点中的位置
while(adr!=0)
{
if(p->next!=NULL)
{
p=p->next;
}
else
{
p->next=new arrayunit; //如果访问位置溢出,空间将自动增长
p=p->next;
}
adr--;
}
return p->num[offset]; //找到后返回引用
}
void openarray::resize(unsigned length)//允许在外部维护一个length,
{ //释放掉length位置后的内存
arrayunit *p=root;
unsigned adr=length/unitsize;
while(adr!=0)
{
if(p->next!=NULL)
{
p=p->next;
}
else
{
p->next=new arrayunit;
p=p->next;
}
adr--;
}
if(p->next!=NULL)
{
delete p->next;
p->next=NULL;
}
} |
除法是这个大数类中最令我兴奋的部分,这个算法是我在某个晚上睡着前在床上想到的(记得当时我完全想错了,后来我纠正的算法,比我预想的稍微复杂一点)
先来看第一部,移位。我假设除数为b,被除数为a,我需要对除数移位(这不是必需的,但我觉得这样做比较高效,如果你试图映射位,会给乘法带来一定的压力,而乘法运算比我实现的移位效率低一些),这也是我把传入的除数参数const去掉的原因(放心,它最后会复原的)。移位很简单,把除数的length移得和被除数一样,直接length减一下,令其等于x,则b<<=x;就搞定了。当然,x要保存的,待会儿回复过程依赖于他。
然后是算法的核心,我们来看最高位,令a的最高位为at,b的位bt,at乘以位的权应该小于a,bt+1乘以b的权大于b。所以我用k=at/(bt+1)对商的大小进行预估,然后先把k传递给输出的对象(注意传递的权,这依赖于之前移位过程的x,k应该被传递给第x位(十进制的9位才算1位)),同时在a-=x*b对新的a进行迭代,重复该过程,直到at<bt+1成立。这个时候x自减1,同时b右移一位。接着问题又来了,这次a比b多了一位!跟刚才不一样了。但是这很简单,现在k=at*(1000000000/(bt+1));进行和上面类似的迭代,直到a降位为止(就是a,b的length相等的时候),然后就重复先前的迭代。这这两个迭代交替,直到结束!注意625行我添加的代码,这个是用来避免溢出的(溢出可能来自第二个迭代,注意到确实有可能很大),事实上我没办法证明它不会溢出,所以额外添加了这个代码(当然这没什么损失,顶多就是这代码从来不运行而已)。
当迭代结束,我得保证b复原!(632行的if条件就是基于此的)
迭代结束后就已经得到结果了吗?可能没有,因为我一直是预估的,我把b算大了,a算小了!所以我需要额外的最后判断,判断他们现在的真实值,也就是644行(如果除法仅仅只有这个东西,也是可以的,但是效率可想而知)。然后就是a<b了!注意到没,a就是余数!所以余数的算法也就有了,而且求余数的时候不用保存商,所以简单很多。当然注意参数,哪些东西要做返回,哪些东西不能变,这个自己传递的时候要注意。
最后稍微看一下效率,注意迭代的次数,这是关键,好吧一般情况下,第一个迭代有50%的概率在1次内结束,平均能在5次内就结束。乘法的效率还可以。第二个迭代衰减更快。制约的因素主要是乘法的效率,注意到用到的乘法都是一般的整数乘以大数,为此我专门写了这个函数(否则用构造函数转化的话,效率就不能保证了),这个函数的效率不会比加法低多少!
***注:权,这个其实很好理解,就拿整数135,最高位1的权是100,然后3的权是10。在这里,9位整数是一个位,位的权就是位的(length-1)*10^9
其实多数的算法并不复杂,只要你明白以下几点就基本没问题了:
- 小学学过的竖式运算是我设计加法减法和乘法的主要算法,这一思想很简单
- 注意竖式运算有位的结合性,你可以把好几位看成一个整体!因此你一下子可以处理9位!
- 注意不要溢出,这点很重要!尤其是乘法中,因为两个单位元素的乘积很有可能溢出,因此你需要long long这样的类型来临时保存它
- 注意,有一种情况是虽然溢出了,但是你却没有损失任何数据!因为unsigned int是32位的,9位整数只用了其30空间,你有额外的两位可以让你“溢出”,这可以使你不需要用一个额外的变量存储,比如在加法中,两个单位相加,可能超过10^9,但是它绝对不会超过4*10^9,所以你可以放心地把和赋给保存和的那个对象的位,然后你把该位/10^9的结果传递给更高的位(就是进位,竖式中的进位还知道吧),把%10^9的结果返回给自身。减法也是,需要注意的是,我的单位是unsigned,没有负数!所以不要指望用负数投机取巧!要是那一位不够, 你应该乖乖的从高位上拿下个10^9来给他用,不要以为这样就够了,万一高位上是0怎么办,那就从更高为上拿!具体看我的代码
- 注意,乘法的特殊性,它很容易溢出,所以随时把那个单位上多余的部分转移到高位!在422行那部分我给出了代码,这个相当简单的。当然这里要注意下,因为我把原来那位上的传递给了s,不要以为我把原来那部分给吞了(刚才我自己都差点这么想,仔细一看原来是这样)。
当然,我的除法用了完全不同的算法,刚开始我试图使用移位配合减法搞定,最后我发现了一个更高效的算法,我会在另一篇文章给出完整分析