자연수 n의 제곱근은 연분수 형태로 나타낼 수 있다. 정수부와 실수부를 각각 나눈다고 생각하고 n의 제곱근을 넘지 않는 정수와 실수부의 합으로 나눈 후, 실수부를 분자가 1이 되는 분수의 형태로 바꾼다. 즉 실수부의 역수를 다시 정수부와 실수부로 나누는 작업을 반복하면서 연분수로 전개한다. 그 과정은 아래와 같이 단계를 일반화 할 수 있다.
\begin{array}{llll} \sqrt{n} & = & A_0 + (\sqrt{n} - A_0) \\ & = & A_0 + \frac{1}{A_1 + \frac{1}{A_2 + \cdots }} \\ \\ \\ A_k + \frac{\sqrt{n} + B_k}{C_k} & = & A_k + \frac{1}{A_{k+1} + \frac{\sqrt{n} + B_{k+1}}{C_{k+1}}}\\ \\ \frac{C_k}{\sqrt{n} + B_k} & = & A_{k+1} + \frac{\sqrt{n} + B_{k+1}}{C_{k+1}} \\ \frac{C_k}{\sqrt{n} + B_k} & = & A_{k+1} + \frac{C_k(\sqrt{n} - B_k}{n - B^2_k}- A_{k+1} \\ & = & \frac{\sqrt{n} - B_k - (A_{k+1}(n - B^2_k) / C_k)}{(n - B^2_k) / C_k} \end{array}
위 계산 과정으로부터 초기값 및 점화식을 아래와 같이 얻을 수 있다.
\begin{array}{lll} & M & = & \lfloor \sqrt{n} \rfloor \\ & A_0 & = & M \\ & B_0 & = & -M \\ & C_0 & = & 1\\ & A_{k+1} & = & \lfloor C_{k}(M - B_k) / (n - B^2_k) \rfloor \\ & B_{k+1} & = & -B_k - (A_{k + 1}(n - B^2_k) / C_k) \\ & C_{k+1} & = & (n - B^2_k) / C_k \end{array}
위 점화식을 사용하여, A, B, C의 다음 항을 각각 구해나간다. 이 때 완전 제곱인 형태이거나, 이전의 a, b, c 값의 조합이 다시 반복되면 순환마디를 발견한 것이라고 보면 된다.
def sqrt(n: int) -> list[int]:
M = int(n ** .5)
if M * M == n:
return [M]
cache = set()
res, a, b, c = [], M, -M, 1
while not ((a, b, c) in cache or b * b == n):
res.append(a)
cache.add((a, b, c))
A = (c * (M - b)) // (n - b**2)
B = -b - (A * (n - b**2) // c)
C = (n - b**2) // c
a, b, c = A, B, C
return res
sqrt(23)
# [4, 1, 8, 1, 3] --> [4; (1 8 1 3)]
# e064 : 연분수 전개의 주기가 홀수인 수의 개수
%%time
xs = map(lambda x: len(sqrt(x + 1)), range(10000))
sum(1 for x in xs if x % 2 == 0)
# 1322
# 266ms