Python code: data structures, algorithms, and OOP.
You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

234 lines
4.8KB

  1. #!/usr/bin/env python
  2. from functools import lru_cache
  3. import math
  4. """
  5. Continued Fraction Representation of Pi
  6. This script computes the convergents
  7. (rational representation) of some continued
  8. fraction representations of Pi.
  9. Happy Pi Day!
  10. 2019-03-14
  11. """
  12. class GeneralContinuedFraction(object):
  13. """Abstract class representing a
  14. general continued fraction.
  15. This class generates terms in the
  16. general continued fraction representation
  17. of a number:
  18. x = b0 + a1
  19. ------------
  20. b1 + a2
  21. --------
  22. b2 + a3
  23. ----
  24. b3 + ...
  25. """
  26. def __init__(self,nterms):
  27. self.nterms = nterms
  28. def __len__(self):
  29. return self.nterms
  30. def a(self,k):
  31. raise NotImplementedError
  32. def b(self,k):
  33. raise NotImplementedError
  34. @classmethod
  35. def str_convergent(cls,n):
  36. n,d = cls.convergent(n)
  37. return "%d / %d"%(n,d)
  38. @classmethod
  39. def dec_convergent(cls,n):
  40. n,d = cls.convergent(n)
  41. return "%0.16f"%(n/d)
  42. @classmethod
  43. def tex_convergent(cls,n):
  44. n,d = cls.convergent(n)
  45. return "\dfrac{ %d }{ %d }"%(n,d)
  46. @classmethod
  47. @lru_cache(maxsize=64)
  48. def convergent(cls,n):
  49. """Compute and return the nth convergent for
  50. x = b0 + a1/( b1 + a2/( b2 + a3/( b3 + a4/( ...
  51. by applying the recurrence formula:
  52. x0 = A0/B0 = b0
  53. x1 = A1/B1 = (b1 b0 + a1)/(b1)
  54. x2 = A2/B2 = (b2 (b1 b0 + a1) + a2 b0)/(b2 b1 + a2)
  55. ...
  56. where An is the numerator and Bn is the denominator,
  57. called continuants, given by the recursion:
  58. A_{n} = b_{n} A_{n-1} + a_{n} A_{n-2}
  59. B_{n} = b_{n} B_{n-1} + a_{n} B_{n-2}
  60. for n > 0 with init values
  61. A_{-1} = 1
  62. A_{0} = b0
  63. B_{-1} = 0
  64. B_{0} = 1
  65. """
  66. def A(k):
  67. if k==-1:
  68. return 1
  69. elif k==0:
  70. return cls.b(0)
  71. else:
  72. return cls.b(k)*A(k-1) + cls.a(k) * A(k-2)
  73. def B(k):
  74. if k==-1:
  75. return 0
  76. elif k==0:
  77. return 1
  78. else:
  79. return cls.b(k)*B(k-1) + cls.a(k) * B(k-2)
  80. return A(n), B(n)
  81. class SimpleContinuedFraction(GeneralContinuedFraction):
  82. """A simple continued fraction is equivalent
  83. to a general continued fraction with all
  84. denominators (a_i) set to 1.
  85. """
  86. @classmethod
  87. def a(cls,k):
  88. return 1
  89. class FourOverPi_ContinuedFraction(GeneralContinuedFraction):
  90. """Class representing Pi
  91. as a general continued fraction
  92. """
  93. @classmethod
  94. def a(cls,k):
  95. if k>0:
  96. return (2*k-1)*(2*k-1)
  97. @classmethod
  98. def b(cls,k):
  99. if k==0:
  100. return 1
  101. else:
  102. return 2
  103. class Pi_ContinuedFraction(FourOverPi_ContinuedFraction):
  104. @classmethod
  105. def str_convergent(cls,n):
  106. # switch numerator and denominator
  107. d,n = cls.convergent(n)
  108. # multiply numerator by 4
  109. n *= 4
  110. return "%d / %d"%(n,d)
  111. @classmethod
  112. def dec_convergent(cls,n):
  113. # switch numerator and denominator
  114. d,n = cls.convergent(n)
  115. # multiply numerator by 4
  116. n *= 4
  117. return "%0.16f"%(n/d)
  118. @classmethod
  119. def tex_convergent(cls,n):
  120. # switch numerator and denominator
  121. d,n = cls.convergent(n)
  122. # multiply numerator by 4
  123. n *= 4
  124. return "\dfrac{ %d }{ %d }"%(n,d)
  125. class Pi_ContinuedFraction2(GeneralContinuedFraction):
  126. """Class representing Pi
  127. as a general continued fraction
  128. (redux)
  129. """
  130. @classmethod
  131. def a(cls,k):
  132. if k>0:
  133. return (2*k-1)*(2*k-1)
  134. @classmethod
  135. def b(cls,k):
  136. if k==0:
  137. return 3
  138. else:
  139. return 6
  140. class Pi_Simple(SimpleContinuedFraction):
  141. terms = [3, 7, 15, 1, 292, 1, 1, 1, 2, 1, 3, 1, 14, 2]
  142. @classmethod
  143. def b(cls,k):
  144. return cls.terms[k]
  145. if __name__=="__main__":
  146. print("%40s"%("Pi (pretty bad)"))
  147. print("-"*80)
  148. for i in range(1,24):
  149. print("%-60s\t%-18s"%(
  150. Pi_ContinuedFraction.str_convergent(i),
  151. Pi_ContinuedFraction.dec_convergent(i),
  152. ))
  153. print("\n\n")
  154. print("%40s"%("Pi (a bit worse)"))
  155. print("-"*80)
  156. for i in range(1,24):
  157. print("%-60s\t%-18s"%(
  158. Pi_ContinuedFraction2.str_convergent(i),
  159. Pi_ContinuedFraction2.dec_convergent(i),
  160. ))
  161. print("\n\n")
  162. print("%40s"%("Pi (simplest and best)"))
  163. print("-"*80)
  164. for i in range(1,14):
  165. print("%-60s\t%-18s"%(
  166. Pi_Simple.str_convergent(i),
  167. Pi_Simple.dec_convergent(i),
  168. ))
  169. print("\n\n")