Borwein quadratic pi

borwein_pi.py

""" Borwein + Borwein quadratic pi
    Algorithm 2.1 from "Pi and the AGM" (Wiley, 1987)

    Written by PM 2Ring 2023.1.2
"""
from math import floor, log, pi
from decimal import Decimal as D, getcontext

def borwein_pi(loops):
    x = D(2).sqrt()
    y = rx = x.sqrt()
    hipi = 2 + x
    
    x = (rx + 1/rx) / 2
    for i in range(1, loops+1):
        hipi *= (x + 1) / (y + 1)
        print(i, hipi)
        if i == loops:
            break
        rx = x.sqrt()
        irx = 1 / rx
        x, y = (rx + irx) / 2, (y*rx + irx) / (y + 1)

def num_digits(k):
    b = pi * (1 << (k+1)) - (k + 4) * log(2) - 2 * log(pi)
    return floor(b / log(10))

@interact
def _(loops=4):
    digits = num_digits(loops)
    print(digits, "digits")
    ctx = getcontext()
    ctx.prec = digits
    borwein_pi(loops)

Live demo running on SageMathCell