#include <iostream>
#include <mpir.h>

using namespace std;

int main(){
    int digits;
    int iterations;
    int prec;
    double lg;

    cout << "Digits: ";
    cin >> digits;

    iterations = digits / 14;
    iterations++;
    lg = log(10) / log(2);
    prec = (lg * digits) + 100;

    mpf_set_default_prec(prec);
    mpf_t sigma;
    mpf_t pi;
    mpz_t topfact;
    mpf_t minusone;
    mpz_t bottomfact1;
    mpz_t bottomfact2;
    mpf_t top;
    mpf_t bottom;
    mpf_t temp;
    mpf_t ftopfact;
    mpf_t fbottomfact1;
    mpf_t fbottomfact2;
    mpf_t k2;
    mpf_init(sigma);
    mpf_init(pi);
    mpz_init(topfact);
    mpf_init(minusone);
    mpz_init(bottomfact1);
    mpz_init(bottomfact2);
    mpf_init(top);
    mpf_init(bottom);
    mpf_init(temp);
    mpf_init(ftopfact);
    mpf_init(fbottomfact1);
    mpf_init(fbottomfact2);
    mpf_init(k2);
    mpf_set_ui(k2, 545140134);
    int k1 = 13591409;
    int k3 = 640320;
    int exp;
    int k = 0;
    int i;
    int target;

    while(k <= iterations){
        cerr << "[" << k << "|" << iterations << "]\r";
        //-1^k
        mpf_set_si(minusone, -1);
        mpf_pow_ui(minusone, minusone, k);
        if(k == 0){
            mpf_set_ui(ftopfact, 1);
            mpf_set_ui(fbottomfact1, 1);
            mpf_set_ui(fbottomfact2, 1);
        }
        else if(k == 1){
            mpf_set_ui(ftopfact, 720);
            mpf_set_ui(fbottomfact1, 6);
        }
        else
        {
            target = 6 * k;
            i = 6 * (k - 1);
            i++;
            for(;i <= target; i++){
                mpf_mul_ui(ftopfact, ftopfact, i);
            }
            target = 3 * k;
            i = 3 * (k - 1);
            i++;
            for(;i <= target; i++){
                mpf_mul_ui(fbottomfact1, fbottomfact1, i);
            }
            mpf_mul_ui(fbottomfact2, fbottomfact2, k);
        }
        //(k1 + k2*k)
        mpf_mul_ui(top, k2, k);
        mpf_add_ui(top, top, k1);
        //top = minusone * topfact * top
        mpf_mul(top, top, minusone);
        mpf_mul(top, top, ftopfact);

        //calculate exponent
        exp = 3*k;
        //k3 ^ exp
        mpf_set_ui(bottom, k3);
        mpf_pow_ui(bottom, bottom, exp);

        //bottom = bottomfact1 * bottomfact2 * bottom
        mpf_mul(bottom, bottom, fbottomfact1);
        mpf_mul(bottom, bottom, fbottomfact2);

        //sigma += top/bottom
        mpf_div(temp, top, bottom);
        mpf_add(sigma, sigma, temp);
        k++;
    }

    mpf_sqrt_ui(temp, 10005);
    mpf_mul_ui(temp, temp, 426880);
    mpf_ui_div(temp, 1, temp);
    mpf_mul(pi, sigma, temp);
    mpf_ui_div(pi, 1, pi);

    gmp_printf("%.Ff", pi);
    return 0;
}