How can I improve my code for the Chudnovsky algorithm to calculate pi more optimised (Python)?

Question

I wrote this code to calculate pi, I'm trying to calculate 2 million digits within at least a minute. This has proven to be really difficult and I can't find any way to further improve this code.

As of now this code can consistently calculate 100,000 digits in around 16 seconds. I've written it based off of the Wikipedia article: Wikipedia's article

My code:

from mpmath import mp, mpf, sqrt
from math import ceil
from time import time

print("Chudnovsky algorithm")

# THE FUNCTION
def chudnovsky(digits):
    mp.dps = digits + 2  # extra precision
    C = mpf(426880) * sqrt(10005)
    L = mpf(13591409)
    M = mpf(1)
    X = mpf(1)
    K = mpf(-6)
    a1 = mpf(0)
    for q in range(ceil(digits / 14.182)):
        a1 += (M * L) / X
        L += 545140134
        X *= -262537412640768000
        K += 12
        M *= (K**3 - 16*K) / (q+1)**3
    result = str(C * a1 ** -1)[:digits + 1]  # cuts off extra potentially inaccurate digits
    mp.dps = 10
    return result




# GET USER INPUT
while True: # loops until an integer is entered
    try:
        user_input = int(input("\nHow many digits?: "))
        break # this only get executed if the user types a valid integer
    except ValueError: # if the user types in something that cant be interpreted as an integer
            print("Enter an integer")

start = time()
print()
print(chudnovsky(user_input))
print()
print(f"Elapsed time: [{time() - start}s]")

Are there any inaccuracies in my code or is this as good as it gets?

Answer 1

I would recommend @Nick Craig-Wood methods to calculate Chudnovsky's Algorithm.

he has explained thoroughly and the timings are present for comparison.

https://www.craig-wood.com/nick/articles/pi-chudnovsky/

A code snippet from the website:

import math
from time import time

def sqrt(n, one):
    """
    Return the square root of n as a fixed point number with the one
    passed in.  It uses a second order Newton-Raphson convgence.  This
    doubles the number of significant figures on each iteration.
    """
    # Use floating point arithmetic to make an initial guess
    floating_point_precision = 10**16
    n_float = float((n * floating_point_precision) // one) / floating_point_precision
    x = (int(floating_point_precision * math.sqrt(n_float)) * one) // floating_point_precision
    n_one = n * one
    while 1:
        x_old = x
        x = (x + n_one // x) // 2
        if x == x_old:
            break
    return x

def pi_chudnovsky(one=1000000):
    """
    Calculate pi using Chudnovsky's series

    This calculates it in fixed point, using the value for one passed in
    """
    k = 1
    a_k = one
    a_sum = one
    b_sum = 0
    C = 640320
    C3_OVER_24 = C**3 // 24
    while 1:
        a_k *= -(6*k-5)*(2*k-1)*(6*k-1)
        a_k //= k*k*k*C3_OVER_24
        a_sum += a_k
        b_sum += k * a_k
        k += 1
        if a_k == 0:
            break
    total = 13591409*a_sum + 545140134*b_sum
    pi = (426880*sqrt(10005*one, one)*one) // total
    return pi

if __name__ == "__main__":
    print(pi_chudnovsky(10**100))
    for log10_digits in range(1,7):
        digits = 10**log10_digits
        one = 10**digits

        start =time()
        pi = pi_chudnovsky(one)
        #print(pi)
        print("chudnovsky: digits",digits,"time",time()-start)

---

output

31415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421172983
chudnovsky: digits 10 time 0.0
chudnovsky: digits 100 time 0.0
chudnovsky: digits 1000 time 0.0
chudnovsky: digits 10000 time 0.027315139770507812
chudnovsky: digits 100000 time 3.12880802154541

---

other methods are much faster that are mentioned on the site.

https://www.craig-wood.com/nick/articles/pi-chudnovsky/