#!/usr/bin/env python from functools import lru_cache import math """ Continued Fraction Representation of Pi This script computes the convergents (rational representation) of some continued fraction representations of Pi. Happy Pi Day! 2019-03-14 """ class GeneralContinuedFraction(object): """Abstract class representing a general continued fraction. This class generates terms in the general continued fraction representation of a number: x = b0 + a1 ------------ b1 + a2 -------- b2 + a3 ---- b3 + ... """ def __init__(self,nterms): self.nterms = nterms def __len__(self): return self.nterms def a(self,k): raise NotImplementedError def b(self,k): raise NotImplementedError @classmethod def str_convergent(cls,n): n,d = cls.convergent(n) return "%d / %d"%(n,d) @classmethod def dec_convergent(cls,n): n,d = cls.convergent(n) return "%0.16f"%(n/d) @classmethod def tex_convergent(cls,n): n,d = cls.convergent(n) return "\dfrac{ %d }{ %d }"%(n,d) @classmethod @lru_cache(maxsize=64) def convergent(cls,n): """Compute and return the nth convergent for x = b0 + a1/( b1 + a2/( b2 + a3/( b3 + a4/( ... by applying the recurrence formula: x0 = A0/B0 = b0 x1 = A1/B1 = (b1 b0 + a1)/(b1) x2 = A2/B2 = (b2 (b1 b0 + a1) + a2 b0)/(b2 b1 + a2) ... where An is the numerator and Bn is the denominator, called continuants, given by the recursion: A_{n} = b_{n} A_{n-1} + a_{n} A_{n-2} B_{n} = b_{n} B_{n-1} + a_{n} B_{n-2} for n > 0 with init values A_{-1} = 1 A_{0} = b0 B_{-1} = 0 B_{0} = 1 """ def A(k): if k==-1: return 1 elif k==0: return cls.b(0) else: return cls.b(k)*A(k-1) + cls.a(k) * A(k-2) def B(k): if k==-1: return 0 elif k==0: return 1 else: return cls.b(k)*B(k-1) + cls.a(k) * B(k-2) return A(n), B(n) class SimpleContinuedFraction(GeneralContinuedFraction): """A simple continued fraction is equivalent to a general continued fraction with all denominators (a_i) set to 1. """ @classmethod def a(cls,k): return 1 class FourOverPi_ContinuedFraction(GeneralContinuedFraction): """Class representing Pi as a general continued fraction """ @classmethod def a(cls,k): if k>0: return (2*k-1)*(2*k-1) @classmethod def b(cls,k): if k==0: return 1 else: return 2 class Pi_ContinuedFraction(FourOverPi_ContinuedFraction): @classmethod def str_convergent(cls,n): # switch numerator and denominator d,n = cls.convergent(n) # multiply numerator by 4 n *= 4 return "%d / %d"%(n,d) @classmethod def dec_convergent(cls,n): # switch numerator and denominator d,n = cls.convergent(n) # multiply numerator by 4 n *= 4 return "%0.16f"%(n/d) @classmethod def tex_convergent(cls,n): # switch numerator and denominator d,n = cls.convergent(n) # multiply numerator by 4 n *= 4 return "\dfrac{ %d }{ %d }"%(n,d) class Pi_ContinuedFraction2(GeneralContinuedFraction): """Class representing Pi as a general continued fraction (redux) """ @classmethod def a(cls,k): if k>0: return (2*k-1)*(2*k-1) @classmethod def b(cls,k): if k==0: return 3 else: return 6 class Pi_Simple(SimpleContinuedFraction): terms = [3, 7, 15, 1, 292, 1, 1, 1, 2, 1, 3, 1, 14, 2] @classmethod def b(cls,k): return cls.terms[k] if __name__=="__main__": print("%40s"%("Pi (pretty bad)")) print("-"*80) for i in range(1,24): print("%-60s\t%-18s"%( Pi_ContinuedFraction.str_convergent(i), Pi_ContinuedFraction.dec_convergent(i), )) print("\n\n") print("%40s"%("Pi (a bit worse)")) print("-"*80) for i in range(1,24): print("%-60s\t%-18s"%( Pi_ContinuedFraction2.str_convergent(i), Pi_ContinuedFraction2.dec_convergent(i), )) print("\n\n") print("%40s"%("Pi (simplest and best)")) print("-"*80) for i in range(1,14): print("%-60s\t%-18s"%( Pi_Simple.str_convergent(i), Pi_Simple.dec_convergent(i), )) print("\n\n")