您的位置:首页 > 其它

计算pi小数点后10000位

2018-03-11 21:53 288 查看
数据结构课上,老师让我们计算一下pi的10000位,复习一下面向对象。

First Try:定义一个高精度类hp(high precision),实现加减乘除,输入输出,然后随便找一个公式,比如Wallis公式,全都用hp进行计算。

现实是:大数乘大数或者大数除大数,实在是太慢了。而且Wallis公式收敛太慢。

Second Try先找到一个收敛较快的公司:BBP公式,观察公式发现,其实只需要定义一个具有如下特性的高精度类hp即可:整数部分用一个int表示,小数部分每8位用一个int保存,实现hp相加减,hp乘或除int即可。

HP.h

#ifndef HP_H
#define HP_H
#include<iostream>
#include<cstring>
#include<algorithm>
#include<iomanip>
using namespace std;
const int maxp=1250,maxn=2*maxp,base=1e8,width=8; //maxn表示双精度

class hp{
public:
int prec,s[maxn+10];   //s[0]保留整数部分,后面依次保留8位小数部分
hp(){ *this = 0;}
hp(int num){ *this = num;}
hp operator=(int num);
hp operator+( const hp&);
hp operator-( const hp&);
hp operator*(const int&);
hp operator/(const int&);
friend ostream& operator<<(ostream &out,const hp&);
};

ostream& operator<<(ostream &out,const hp& x){
out<<x.s[0]<<"."<<endl;;
for(int i=1;i<=x.prec;++i) {
out<<setfill('0')<<setw(width)<<x.s[i]<<" ";
if(i%10==0) out<<endl;
}
return out;
}

hp hp::operator=(int num){prec=0;s[0]=num;return this;}
hp hp::operator+( const hp& b){
hp c;
c.prec=max(prec,b.prec);
int  x=0;
for(int i=c.prec;i>=0;--i){
if(i<=prec) x+=s[i];
if(i<=b.prec) x+=b.s[i];
c.s[i]=x%base;
x/=base;
}
return c;
}

hp hp::operator -(const hp&b){
hp c;
int x=0;
c.prec=max(prec,b.prec);
for(int i=c.prec;i>=0;--i){
if(i<=prec) x+=s[i];
if(i<=b.prec) x-=b.s[i];
if(x<0){c.s[i]=x+base;x=-1;}
else{c.s[i]=x;x=0;}
}
return c;
}
hp hp::operator *(const int &b){
hp c;
long long x=0;  //int * int会溢出,用long long过渡
c.prec=prec;
for(int i=prec;i>=0;--i){
x=x+s[i]*(long long )b;
c.s[i]=x%base;
x/=base;
}
return c;
}
hp hp::operator/(const int &b){
hp c;
long long x=0;
int i;
for( i=0;i<=maxp;++i){
if(x==0&&i>prec) break;
x=x*base;
if(i<=prec) x+=s[i];
c.s[i]=x/b; x=x%b;
}
c.prec=i-1;
return c;
}

#endif


main()

#include "HP.h"
#include<iostream>
#include<cstdio>
using namespace std;
int main(){
hp pi=0,e=1;
int tmp=0;
for(int i=0;i<10000;++i){
pi=pi+e*4/(tmp+1);
pi=pi-e*2/(tmp+4);
pi=pi-e/(tmp+5);
pi=pi-e/(tmp+6);
tmp+=8;
e=e/16;
}
cout<<"Pi's first "<<pi.prec*8<<" "<<"digits\n";
cout<<pi<<endl;
cout << "Totle Time : " << (double)clock() /CLOCKS_PER_SEC<< "s" << endl;
return 0;
}
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: