Chudnovsky algorithm produces -nan

Question

I'm trying to implement the Chudnovsky algorithm for calculating pi.

Here is my implementation:

int fact(int n)
{
    if(n<=1)
        return 1;
    else
        return fact(n-1)*n;
}

double calcPi(long n)
{
    double z=0;
    for(int k=0; k<n; k++)
    {
        z+=(pow(-1, k)*fact(6*k)*(13591409 + 545140134.0*k))/(fact(3*k)*pow(fact(k), 3)*pow(640320.0, 3.0*k+3.0/2));
    }
    z*=12;
    return 1/z;
}

I'm running into a tiny error though. When I plug in values of N that are greater than 12, I get -nan. I'm guessing this has to do with limited precision, some sort of integer overflow, or my absolutely terrible factorial implementation (yes, I was lazy and used recursion. It's 2am).

Anyways, if you've been through this before and can suggest a quick fix, that would be nice.

Maybe I should just use Python, and stop worrying about the overflows.

Happy (almost) New Years!

Answer 1

floating point arithmetics it's not trivial and considering your problem i prefer to answer your question with some tip.

You can solve this with a library such as GMP or MPFR , and this is a good FAQ for both.

If you really want to master this, on almost every major programming language, you should absolutely start from reading the IEEE 754.