10000 decimal places after calculating pi

In the data structure class, the teacher asked us to calculate the 10000 bits of pi and review the object orientation.

First Try: define a high precision class hp(high precision), implement addition, subtraction, multiplication and division, input and output, and then find any formula, such as Wallis formula, all of which are calculated by hp.

The reality is: it's too slow to multiply a large number by a large number or divide a large number by a large number. And Wallis formula converges too slowly.

Second Try first finds a company with fast convergence: BBP formula. According to the observation formula, it is only necessary to define a high-precision class hp with the following characteristics: the integer part is represented by an int, and every 8 digits of the decimal part are saved by an int to realize hp addition and subtraction, hp multiplication or division of 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 for double precision

class hp{
    public:
        int prec,s[maxn+10];   //s[0] keeps the integer part, followed by 8 decimal parts
        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;}
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 will overflow, use long long to transition
    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;
}

Posted by asuperstar103 on Fri, 03 Apr 2020 00:41:11 -0700