如何算大数除法？

[rickone]

我在求e的时候得到了e=2+a/b
a,b是两个大数，为了精确度的需要成了大数，可是我现在不知道怎么算除法了。。。。

[NowCan]

有没有小数部分？这很关键。
如果都是整数，可以看看这个，效率虽然不高，不过算法比较好懂。
http://www.pediy.com/tools/Cryptography/RsaKit/RsaKit_V1.0.rar

[rickone]

怎么只有取整除法和取余啊
我要的就是要小数部分的，看我上面的,e=2+a/b，显然是要一个小数啊……

[kerbcurb]

这是以前写的，其中被除数以c[]表示，除数以a[]表示，b[]表示商，BASE = 10表示10进制，LENGTH = 500表示位数。你还需要稍微改一下，仅供参考

const int LENGTH = 500;
const int BASE = 10;

int main()
{
int a[LENGTH] = {1,2,3,4,5,6,7,8,9};
int c[LENGTH] = {9,8,7,6,5,4,3,2,1};
int b[LENGTH] ;
for(int j = 0;j < LENGTH;j++)b[j] = 0;
int temp;
b[0] = c[0] / a[0];
for(int k = 0;k < LENGTH;k++)
{
temp = c[k];
for(int i = 0;i < k;i++)
temp -= a[k - i] * b[i];

b[k] = temp / a[0];
}
int s = 0;

for(int k = LENGTH - 1;k >0;k--)
{
s = b[k];
b[k] %= BASE;
b[k - 1] += s / BASE;
}

for(int k = LENGTH - 1;k >0;k--)
{
if(b[k] < 0)
{
s = b[k];
b[k] = (BASE + b[k] % BASE);
b[k - 1] += s / BASE - 1 ;
}
else
{
s = b[k];
b[k] %= BASE;
b[k - 1] += s / BASE;
}

}
cout<<endl;
for(int k = 0;k < LENGTH;k++)printf("%d",b[k]);

return 0;
}

[mathe]

计算e*x^k (如果你要计算到小数点后面k位），x是你使用的进制（比如x=10或16)
e*b^k=2*x^k + (a*x^k)/b
所以计算(a*x^k)/b就可以了，现在使用整除法就可以了。

[ygd]

看看gmp库

[yaos]

求e不用大数运算的

[rickone]

不用大数怎么算？

[rickone]

我是用的公式：
e=2+1/2!+1/3!+1/4!+...
但是这么算能精确到小数点后多少位呢？我想要很精确的，怎么办？

[mathe]

想多少位都可以，只要时间够长

[mathe]

计算到一万位
#define M (10000/4)
#define N 3250 // 3250*(log(3250)-log(e))>10000, so 3250!<10^10000
int x[M],y[M];
void init(){//initialize y=x to 0.5
int i;
x[ 0 ]=5000;
for(i=1;i<M;i++)x[ i ]=0;
for(i=0;i<M;i++)y[ i ]=x[ i ];
}

void step(int k){
//divide y by k
int i;
unsigned c=0;
for(i=0;i<M;i++){
int u=c*10000+y[ i ];
y[ i ]=u/k;
c=u%k;
}
//Add x by y
c=0;
for(i=N-1;i>=0;i--){
int u=x[ i ]+y[ i ]+c;
if(u>=10000){
c=1;u-=10000;
}else c=0;
x[ i ]=u;
}
}

int main(){
int i;
init();
for(i=3;i<=N;i++){
step(i);
}
printf("2.");
for(i=0;i<M;i++)printf("%04d",x[ i ]);
printf("\n");
}

[mathe]

呵呵，上面程序有一处Bug,你自己修改吧:)

[mathe]

直接用数组表示一个很大的数，自己用计算机模拟手工计算乘除法的方法。使用泰勒级数来算e就可以了。

[rickone]

能描述一下算法吗？

我刚在网上下了个牛逼的算pi的程序：
long a=10000,b,c=2800,d,e,f[2801],g;
　　
　　main()
　　{
　　for(;b-c;)
　　f[b++]=a/5;
　　for(;d=0,g=c*2;c-=14,printf("%.4d",e+d/a),e=d%a)
　　for(b=c;d+=f[b]*a,f[b]=d%--g,d/=g--,--b;d*=b);
　　}
还有算e的，类似：
long a=10000, b,c=2800,d,e=10000,f[2801],g;
　　
　　int main()
　　{
　　for(;b-c;)
　　f[b++]=a;
　　for(;d=0,g=c;c-=14,printf("%.4d ",e+d/a),e=d%a)
　　for(b=c;d+=f[b]*a,f[b]=d%b,d/=b,--b;);
　　return 0;
　　
　　}

我靠，这是怎么算的啊，我没个思路，怎么把小数位一位位算出来的。请大侠帮我分析一下！

[rickone]

但是程序是怎么做到的，它似乎好像用的是迭代，怎么把和式转化成迭代式，还有误差怎么处理的？

[NowCan]

你先看mathe的程序吧，小数都放在了一个数组里，每个元素存放三位小数。好像是这样，我没仔细看，错了请mathe指出。
至于你贴的两个程序，你现在功力还不够。不可操之过急。

[mathe]

就是用公式e=2+1/2!+1/3!+1/4!+...
呀

[mathe]

先把我的代码变成：
for(k=3;k<N;k++){
c[0]=0;
for(i=0;i<M;i++){
int u=c[ i ]*10000+y[ i ];
y[ i ]=u/k;
c[ i+1 ]=u%k;
}
//Add x by y
d[ i ]=0;
for(i=M-1;i>=0;i--){
int u=x[ i ]+y[ i ]+d[ i ];
d[ i+1 ]=u/10000;
x[ i ]=u%10000;
}
}
然后改变代码的执行顺序，要让k变成内循环先执行。为此，需要先将两个for(i的循环合并，
为此，我们需要将两个for(i=...)的循环都要变成i递减的循环。
但是由于第一个循环中有跨循环的数据依赖关系(Loop-carried data dependence),我们无法直接交换
所以对循环还要采取线性变换，就是那种可以将矩形变成平行四边行方式的变换。
在这里，它的代码还使用了一个特性，也就是1/k!的计算中，对于比较大的k,算出来的数字很小，所以前面会有很多0,所以计算前面的位数时，我们不需要累加的1/N!,而只要比较小的就可以了，而它的程序是表明，对于计算道第4k位(10进制位），只要计算道1/(14k)!就可以了。

[mathe]

贴一个比较完整的：
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
int M,N;
int starty;
int *x,*y;
void init_N()
{
double n=M*8;
double nn;
do{
nn=n;
n=(log(10)*M*8+nn)/log(nn);
}while(fabs(nn-n)>=0.1);
N=(int)n+1;
}

void init(){//initialize y=x to 0.5
int i;
x=(int *)malloc(sizeof(int)*M);
y=(int *)malloc(sizeof(int)*M);
if(x==NULL||y==NULL){
printf("Out of memory\n");
exit(-2);
}
x[ 0 ]=50000000;
for(i=1;i<M;i++)x[ i ]=0;
for(i=0;i<M;i++)y[ i ]=x[ i ];
init_N();
}
#ifdef _WIN32
typedef __int64 longlong;
#else
typedef long long longlong;
#endif
void step(int k){
//divide y by k
int i;
longlong c=0;
for(i=starty;i<M;i++){
longlong u=c*100000000+y[ i ];
y[ i ]=(int)(u/k);
c=(int)(u%k);
}
while(y[starty]==0&&starty<M)starty++;
//Add x by y
c=0;
for(i=M-1;i>=starty;i--){
int u=x[ i ]+y[ i ]+c;
if(u>=100000000){
c=1;u-=100000000;
}else c=0;
x[ i ]=u;
}
while(c&&x[i]==99999999){x[i]==0;i--;}
if(c)x[i]++;
}

int main(int argc, char *argv[]){
int i,realM;
if(argc>=2){
M=atoi(argv[1]);
if(M<=1)return -1;
}else{
printf("Please input number of digits of e to be calculated:\n");
if(scanf("%d\n",&M)!=1||M<=1)
return -1;
}
realM=M;
M=(M+7)/8+1;
init();
for(i=3;i<=N;i++){
step(i);
}
printf("2.");
if(x[M-1]>=50000000){
for(i=M-2;x[i]==99999999&&i>=0;i--){
x[i]=0;
}
x[i]++;
}
for(i=0;i<realM/8;i++)printf("%08d",x[ i ]);
if(realM%8>0){
int left=realM%8;
int m=1,j;
for(j=0;j<8-left;j++)m*=10;
printf("%0*d",left,(x[i]+m/2)/m);
}
printf("\n");
}