読者です 読者をやめる 読者になる 読者になる

Arantium Maestum

プログラミング、囲碁、読書の話題

SICPの勉強 問題1.36

SICP1.3.3の問題1.36を解いてみる。

と、その前にまずfixed-pointをClojureっぽく書いてみる。

シーケンス[a0 a1 ... an-1 an]から[a0 a1 ... ai-2 ai-1] [a1 a2 ... ai-1 ai] [a2 a3 ... ai ai+1]と続くベクトルのシーケンスを作成する関数:

(defn conseq [i seq]
  (->> seq
       (iterate rest)
       (take i)
       (apply map vector)))

これってcoreで定義されてても良くない?といつも思う、シーケンスから条件に合致する最初の要素を取る関数:

(defn first-match [pred seq]
  (first (filter pred seq)))

これらを合わせてfixed-point関数を作る:

(defn fixed-point [f guess]
  (let [tolerance     0.00001
        close-enough? (fn [[a b]] (< (Math/abs (- a b)) tolerance))]
    (->> (iterate f guess)
         (conseq 2)
         (first-match close-enough?)
         first)))

これで平方根を定義してみる:

(defn average [x y]
  (/ (+ x y) 2))

(defn sqrt [x]
  (fixed-point (fn [y] (average y (/ x y)))
               1.0))

と、ここまでが問題を解決する前のコード。

次にfixed-pointを修正して、最終的な答えが出るまでの中間結果もリストにして返すようにしていく。

シーケンスからpredが合致する前の要素全てと、合致した最初の要素一つを返す関数:

(defn take-upto [pred seq]
  (let [[head [t & _]] 
        (split-with (complement pred) seq)]
    (concat head [t])))

これをfixed-point関数に組み入れて中間結果も返す関数を定義:

(defn intermediate-fixed-point [f guess]
  (let [tolerance     0.00001
        close-enough? (fn [[a b]] (< (Math/abs (- a b)) tolerance))]
    (->> (iterate f guess)
         (conseq 2)
         (take-upto close-enough?)
         (map first))))

x^x = 1000の答えを求めてみる:

(intermediate-fixed-point 
  #(/ (Math/log 1000) (Math/log %))
  2)

(2 9.965784284662087 3.004472209841214 6.279195757507157 3.759850702401539 5.215843784925895 4.182207192401397 4.8277650983445906 4.387593384662677 4.671250085763899 4.481403616895052 4.6053657460929 4.5230849678718865 4.577114682047341 4.541382480151454 4.564903245230833 4.549372679303342 4.559606491913287 4.552853875788271 4.557305529748263 4.554369064436181 4.556305311532999 4.555028263573554 4.555870396702851 4.555315001192079 4.5556812635433275 4.555439715736846 4.555599009998291 4.555493957531389 4.555563237292884 4.555517548417651 4.555547679306398 4.555527808516254 4.555540912917957)