프로젝트 오일러 64

내 이 세상 도처에서 쉴 곳을 찾아보았으나, 마침내 찾아낸, 컴퓨터가 있는 구석방보다 나은 곳은 없더라.

on Project-Euler Clojure

프로젝트 오일러 64

제곱근을 연분수로 나타낼 때 반복 주기가 홀수인 경우 세기

문제 자세히 보기: [국어] [영어]

위키피디아에 연분수로 제곱근을 구하는 방법이 자세히 설명되어 있다.

%math \begin{aligned} \sqrt{N} = a_0 + \cfrac{1}{a_1 + \cfrac{1}{a_2 + \cfrac{1}{a_3 + \cfrac{1}{\ddots}}}} \end{aligned}

제곱근을 위와 같이 연분수로 나타냈을 때 다음 알고리즘으로 a_n 을 구할 수 있다. N 은 자연수이고 완전제곱수가 아니어야 한다.

%math \begin{aligned} m_0 &= 0,\quad d_0 = 1,\quad a_0 = \lfloor{\sqrt{N}}\rfloor \newline m_{n+1} &= d_n a_n - m_n \newline d_{n+1} &= \frac{N - m^2_{n+1}}{d_n} \newline a_{n+1} &= \left \lfloor \frac{a_0 + m_{n+1}}{d_{n+1}} \right \rfloor \end{aligned}

a_i = 2a_0 가 되면 반복이 끝난다. 이 알고리즘을 다음과 같이 Clojure 함수로 옮길 수 있다.

(defn expand-continued-fraction [n]
  (let [a0 (int (Math/sqrt n))]
    (loop [m 0, d 1, a a0, acc [a0]]
      (if (= a (* 2 a0))
        acc
        (let [m (- (* d a) m), d (/ (- n (* m m)) d), a (int (/ (+ a0 m) d))]
          (recur m d a (conj acc a)))))))

이제 N \le 10000 일 때 반복주기가 홀수인 경우를 세기만 하면 된다. N 은 완전제곱수이면 계산할 때 분모가 0 이 되는 경우가 생겨 ArithmeticException이 발생하므로 완전제곱수는 제외하고 계산해야 한다.

(defn square [x] (* x x))

(defn solve []
  (->> (range 2 (inc 10000))
       (filter #(not= % (square (int (Math/sqrt %)))))
       (map expand-continued-fraction)
       (filter #(even? (count %)))
       count))

실행 결과는 다음과 같다.

p064> (time (solve))
"Elapsed time: 135.54102 msecs"
1?22

참고