Calculating pi with Catalan Numbers: A Novel Approach
Introduction
Calculating the digits of pi has been a topic of interest for mathematicians and computer scientists alike. In this post, I introduce an algorithm for pi calculation that leverages Catalan numbers and fixed-point arithmetic.
Background
History of pi Calculation Methods
Historically, various methods have been used to approximate the digits of pi, from Archimedes' method to the more recent Chudnovsky algorithm.
Catalan Numbers
Catalan numbers are a sequence of natural numbers that have many applications in combinatorial mathematics. They can be defined using the recurrence relation:
Cn = [2(2n - 1) / (n + 1)] Cn-1, with C0 = 1.
The Algorithm
Series Representation
The series representation for pi used in this algorithm is derived from the inverse sine function:
pi = 6 arcsin(0.5)
By expanding arcsin(x) as a series, we get:
arcsin(x) = sum from n = 0 to infinity of [(2n)! / (4n(n!)2(2n + 1))] x2n + 1
Substituting x = 0.5 and multiplying by 6, we derive the series:
pi = 3 + 6 sum from n = 1 to infinity of [(2n - 1) Cn-1 / ((2n + 1)16n)]
Where Cn is the nth Catalan number.
Try this in the Lab to calculate a Catalan value live from the same public endpoint.
Optimization
To enhance computational efficiency, the implementation accumulates scaled integer terms and performs one long division at the end. This keeps the calculation in fixed-point arithmetic and reduces repeated high-precision division inside the loop.
Implementation
The implementation in Python leverages the gmpy2 library for high-precision integer arithmetic:
from gmpy2 import mpz
import sys
# Input number of digits
digits = mpz(eval(input('How many digits of pi? :\t')))
iters = digits * 5 // 3
# Initialize variables
pi = mpz()
c = d1 = d2 = m = mpz(1)
# Iterative calculation
for i in range(2, iters + 2):
pi *= 16 * (m + 2)
pi += c * m * d2
c = 2 * c * m // i
d1 *= 16
m += 2
d2 *= m
if i % 1000 == 0:
print("\r%d terms..." % i, end='')
# Final calculation
one = mpz(10) ** digits
pi = 6 * pi * one // (d1 * d2) + 3 * one
pi_digits = f'{pi // one}.' + str(pi % one).zfill(digits)
# Output result
if sys.argv[1:2] in (['-l'], ['--log']):
with open((sys.argv[2:3] or [f'pi{digits}.txt'])[0], 'w') as f:
f.write(pi_digits)
else:
print(f'\rpi = {pi_digits}...infinity')Explanation of the Code
This code performs the calculation with integers, carrying the denominator growth through d1 and d2, then scaling by 10 ** digits before the final division. That matches the fixed-point approach used in the implementation without claiming that the code explicitly computes an LCM.
Conclusion
This approach to calculating the digits of pi is both efficient and intriguing due to its use of Catalan numbers. Further optimization and analysis could lead to even more efficient algorithms.