# OK to post homework
# James Betti, 12 Mar 2025, Assignment 14 (Pi Day)

JesusG := proc(N) local total,i,x;
    total := 0;
    for i from 0 to N do
        total := total+(-1)^i*((2*i)!/(2^i*i!)^2)^5*(20*i*(9+41*i)+13)/2^(10*i);
    od;
    evalf(sqrt(128/total));
end:

# https://en.wikipedia.org/wiki/Chudnovsky_algorithm
# Chudnovsky algorithm using binary splitting as in the above article.
Chudnovsky := proc(N) local split,P,Q,R;
    split := proc(a,b) local P1,Q1,R1,P2,Q2,R2;
        P,Q,R;
        if b=a+1 then
            P := -(6*a-1)*(2*a-1)*(6*a-5);
            Q := 10939058860032000*a^3;
            R := P*(545140134*a+13591409);
        else
            P1,Q1,R1 := split(a,floor((a+b)/2));
            P2,Q2,R2 := split(floor((a+b)/2),b);
            P := P1*P2;
            Q := Q1*Q2;
            R := Q2*R1+P1*R2;
        fi;
        evalf(P),evalf(Q),evalf(R);
    end;
    P,Q,R := split(1,N);
    evalf(426880*sqrt(10005)*Q/(13591409*Q+R));
end:

Digits := 100:
JesusG(38);
Chudnovsky(8);

Chudnovsky(1000000);             # 6.92 seconds