SICPの勉強 問題1.29

SICP1.3.1の問題1.29を解いてみる。

とりあえずSimpson Ruleを使わない積分法：

(defn sum [f a next b]
  (if (> a b)
    0
    (+ (f a)
       (sum f (next a) next b))))

(defn integral [f a b dx]
  (letfn [(add-dx [x] (+ x dx))]
    (* (sum f (+ a (/ dx 2.0)) add-dx b)
       dx)))

これにcubeを定義して入れてやる：

(defn cube [x] (* x x x))
(integral cube 0 1 0.001)

結果は0.249999875000001。ぴったりとはいかない。ちなみにdxを0.0001にするとStackOverflowになる。解決法は1.30問で。

それに対してSimpson RuleをClojureの高階関数を利用して書いてみる：

(defn simpson [f a b n]
  (let [h (* 1.0 (/ (- b a) n))
        y (fn [k] (f (+ a (* k h))))]
    (* (/ h 3)
       (+ (f a)
          (* 4 (apply + (map y (range 1 n 2))))
          (* 2 (apply + (map y (range 2 n 2))))
          (f b)))))

これだとn=10ですらぴったり0.25と出る。（hの定義で無理やり浮動小数点にしているのに・・・）

ちなみにさっき定義したsumを使うなら：

(defn simpson-sum [f a b n]
  (let [h (* 1.0 (/ (- b a) n))
        y (fn [k] (f (+ a (* k h))))]
    (* (/ h 3)
       (+ (f a)
          (* 4 (sum y 1 #(+ 2 %) (- n 1)))
          (* 2 (sum y 2 #(+ 2 %) (- n 2)))
          (f b)))))

あまり好みじゃない。sumを数学のΣに寄せてあるのはわかるのだが、clojure（というか関数型言語）の高階関数の組み合わせで表現したほうが何が起きているか、すっきりとわかる気がする。