Skip to content

Instantly share code, notes, and snippets.

@rubik
Created December 10, 2011 10:32
Show Gist options
  • Star 1 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save rubik/1454917 to your computer and use it in GitHub Desktop.
Save rubik/1454917 to your computer and use it in GitHub Desktop.
Creating a simple continued fraction from a square root
import math
def from_root(n):
'''
Construct a continued fraction from a square root. The argument
`n` should be an integer representing the radicand of the root:
>>> from_root(2)
(1, [2])
>>> from_root(4)
(2,)
>>> from_root(97)
(9, [1, 5, 1, 1, 1, 1, 1, 1, 5, 1, 18])
'''
a0 = int(math.floor(math.sqrt(n)))
r = (a0, [])
a, b, c = 1, 2 * a0, a0 ** 2 - n
delta = math.sqrt(4*n)
while True:
try:
d = int((b + delta) / (2 * c))
except ZeroDivisionError: # a perfect square
return (r[0],)
a, b, c = c, -b + 2*c*d, a - b*d + c*d ** 2
r[1].append(abs(d))
if abs(a) == 1:
break
return r
import math
import fractions
def from_root(n):
'''
Construct a continued fraction from a square root. The argument
`n` should be an integer representing the radicand of the root:
>>> from_root(2)
(1, [2])
>>> from_root(4)
(2,)
>>> from_root(97)
(9, [1, 5, 1, 1, 1, 1, 1, 1, 5, 1, 18])
'''
coeff = 1
floor_part = floor_ = math.floor(math.sqrt(n))
denom = n - floor_part ** 2
exp = (int(floor_), [])
while True:
try:
floor_ = math.floor((math.sqrt(n) + floor_part) / float(denom))
except ZeroDivisionError: # perfect square
return (exp[0],)
if denom != 1:
exp[1].append(int(floor_))
floor_part = denom * floor_ - floor_part
coeff = denom
denom = n - floor_part ** 2
common_div = fractions.gcd(coeff, denom)
coeff /= common_div
denom /= common_div
if denom == 1:
exp[1].append(int(floor_part + exp[0]))
break
return exp
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment