計算pi小數點後10000位
阿新 • • 發佈:2019-01-25
資料結構課上,老師讓我們計算一下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; }