If its true that irrational numbers goes on forever then building one of the Pi computation engine would be a tough thing. Then a bright light of hope from far googleeast revealed the almighty Spigot Algorithm. Spigot Algorithm is simple and only needs integer arithmetics to compute the value of Pi. You can check wiki about the algorithm here . Just follow the links on references. I read one from Rabinowitz, Stanley; Wagon, Stan (1995). The algorithm presented on the paper uses a mixed base representation of PI to calculate its value.
Firstly I used python code as protype to compute first 1000 digits. As people never use python for its speed, I shifted to C++. Using C++ code, it took near about 0.4s in my laptop(i3 dual core) to compute first 1000 digits which was happily far better than python code. Then I tried 10000 digits, the time was near about 18 seconds. Then 100000 digits. This one made me wait quite long time. Waiting was tough. While my laptop was busy calculating digits, I grabbed a sheet of paper and roughly calculated the time. Considering 2Ghz processor and total loop counts, 27 minutes was my conclusion. And I was pretty much sure that it would take more than that. To my surprise, the program stopped with time 27.906 minutes :D. However the program got final 3 digits wrong. A quick google showed it was a perfect valid condition for a spigot algorithm. They inherit such problems in general and are less likely to be used for calculating higher digits.
And here I present my code snippet I had built. First a quick python prototype.
def reform():
    i = -1;
    while(pirepo[i] >= 10):
        q = (int) (pirepo[i]/10)
        pirepo[i] %= 10;
        i -= 1
        pirepo[i] += q
    return
 
 
digits = 1000
bucketsize = (int)(10 * digits / 3) + 1
bucket = [2 for i in xrange(bucketsize)]
 
pirepo = []
 
for i in xrange(digits):
    bucket = [(bucket[j] * 10) for j in xrange(bucketsize)]
 
    quotent = 0
    for j in xrange(bucketsize):
        index = bucketsize - j - 1
        bucket[index] = ( bucket[index] + quotent * (index + 1))
        quotent = (int)(bucket[index] / (2 * index + 1))
        bucket[index] = bucket[index] % (2 * index + 1)
 
    bucket[0] = quotent % 10
    output = (int)(quotent / 10)
 
    pirepo.append(output)
 
    reform()
 
print pirepo
#include 
class PI{
private:
    const static int digits = 100000;
    const static int bucketSize = (10*digits/3) + 1;
 
    int piStore[digits];
    int bucket[bucketSize];
 
public:
    void calculate(){
        for(int i = 0; i < bucketSize ; i++)
            bucket[i] = 2;
 
        int q = 0;
        for (int i = 0; i < digits ; i++){
            q = 0;
            for (int j = bucketSize-1; j >= 0 ; j--){
                bucket[j] *= 10;
                bucket[j] = bucket[j] + (q * (j + 1));
                q = bucket[j] / (2 * j + 1);
                bucket[j] = bucket[j] % (2 * j + 1);
            }
            piStore[i] = q/10;
            bucket[0] = q%10;
 
        }
 
        revalidate();
 
        for ( int i = 0 ; i < digits ; i++)
            std::cout << piStore[i];
    }
    void revalidate(){
        int q;
        int index;
        for (int i = 0 ; i < digits ; i++){
            q = 0;
            index = 0;
            while (piStore[i-index]>9){
                q = piStore[i-index]/10;
                piStore[i-index] %= 10;
                index++;
                piStore[i-index] += q;
            }
        }
    }
};
int main(){
    PI piObj;
    piObj.calculate();
    return 0;
}
  
Sir/mam Very nice and educational article. I am your daily reader and follower in every social media. Thank You So much and please keep it up. Keep Providing us educational articles like this. Thank You . I just got motivated from this.
ReplyDeleteCheers : Raj Mahat
All It Feed
Nice. https://www.neptechofficial.com/2019/01/how-to-read-and-write-files-in-java.html
ReplyDelete