第5回社内SICP勉強会メモ

SICP

第４回の話題

こちらに書きました

範囲

P26の1.2.5最大公約数からP30の問題1.24まで

内容

aをbで割った剰余をrとすると、aとbの公約数はbとrの公約数と全く同じという考え方

(define (gcd a b)
  (if (= b 0)
      a
      (gcd b (remainder a b))))
  • Lameの定理
  • 除数の探索
    • 2から初めて√nまで除算をしていって、途中で割り切れるものが出て来たら素数ではない、という方法
    • 1から√nまでの間の数を１ずつずらしながら探索するのでステップ数はΘ(√n)となる。
  • Fermatテスト
    • nが素数のとき、a
    • n=5, a=3で確かめると、確かに合同になる。
    • これを利用して、与えられたnより小さい様々な値において、a^nの剰余を計算してaと比較するという作業を繰り返す
      • 繰り返し回数が多いほど信頼性が増す

問題

問題1.20

GCD手続きを正規順序と作用的順序の２パターンで置き換えモデルを使って展開する。
正規順序の場合は、remainderも展開せずに次を呼び出す。
ifは特殊形式なので、そこに与えられる(remainder a b)のみ計算されて、そこが0になったときに初めて全ての式が評価される

作用的順序は適宜remainderを計算しながら展開していく。展開済みの値がifに与えられて0になったら終了で値を返していく

あとで書こうと思ったけど、紙から転載するのがめんどくさすぎるので諦める。
上で書いた要点に従って展開していくと

  • 正規順序だとremainder手続きが18回
  • 作用的順序だとremainder手続きが4回実行される

問題1.21

1.2.6で書いたsmallest-divisorをそのまま利用する

;; 問題1.21
(smallest-divisor 199)
;; gosh> 199

(smallest-divisor 1999)
;; gosh> 1999

(smallest-divisor 19999)
;; gosh> 7

問題1.22

;; 問題1.22
;; srsi-19を利用する
;; あと順次実行のためにbeginも利用している
;; ので、なんか負けた気分
;; 他にいい方法はないものか
(use srfi-19)

(define (smallest-divisor n)
  (find-divisor n 2))

(define (find-divisor n test-divisor)
  (cond [(> (square test-divisor) n) n]
	[(divides? test-divisor n) test-divisor]
	[else (find-divisor n (+ test-divisor 1))]))

(define (square x)
  (* x x))

(define (divides? a b)
  (= (remainder b a) 0))

(define (prime? n)
  (= n (smallest-divisor n)))

(define (report-time elapse-time)
  (display " *** ")
  (display (time-nanosecond elapse-time))
  (newline))

;; 方針
;; 偶数なら飛ばして奇数にする
;; 奇数で素数かどうかを判定
(define (search-for-primes start end)
  (define (timed-search-for-primes n start-time cnt)
    (cond [(or (> n end) (= cnt 3)) #t]
	  [(even? n) (timed-search-for-primes (+ n 1) (current-time) cnt)]
	  [(prime? n)
	   (begin (display n)
		  (report-time (time-difference (current-time) start-time))
		  (timed-search-for-primes (+ n 2) (current-time) (+ cnt 1)))]
	  [else (begin (display n)
		       (newline)
		       (timed-search-for-primes (+ n 2) (current-time) cnt))]))
  (timed-search-for-primes start (current-time) 0))

(define (start-test start)
  (newline)
  (search-for-primes start (* start 10)))

(start-test 1000)
;;gosh> 
;;1001
;;1003
;;1005
;;1007
;;1009 *** 27000
;;1011
;;1013 *** 26000
;;1015
;;1017
;;1019 *** 26000

(start-test 10000)
;;gosh> 
;;10001
;;10003
;;10005
;;10007 *** 55000
;;10009 *** 55000
;;10011
;;10013
;;10015
;;10017
;;10019
;;10021
;;10023
;;10025
;;10027
;;10029
;;10031
;;10033
;;10035
;;10037 *** 54000

(start-test 100000)
;;gosh> 
;;100001
;;100003 *** 132000
;;100005
;;100007
;;100009
;;100011
;;100013
;;100015
;;100017
;;100019 *** 132000
;;100021
;;100023
;;100025
;;100027
;;100029
;;100031
;;100033
;;100035
;;100037
;;100039
;;100041
;;100043 *** 132000

(start-test 1000000)
;;gosh> 
;;1000001
;;1000003 *** 416000
;;1000005
;;1000007
;;1000009
;;1000011
;;1000013
;;1000015
;;1000017
;;1000019
;;1000021
;;1000023
;;1000025
;;1000027
;;1000029
;;1000031
;;1000033 *** 430000
;;1000035
;;1000037 *** 423000

;;gosh> (/ 55000.0 27000.0)
;;2.037037037037037
;;gosh> (/ 132000.0 55000.0)
;;2.4
;;gosh> (/ 423000.0 132000.0)
;;3.20454545454546

;; (sqrt 10) = 3.162227
;; なので、数が小さい時は予想の増加の程度を下回っている
;; 一方、数が大きくなると予想の精度に近似的になる
;; これは全体の計算量に占める素数判定以外の部分が原因と思われる
;; 数が少ない時には、素数判定以外の割合が大きく、その部分の影響のため、
;; 単純に素数判定部分の計算量の見積もり以上に、増加量の割合が小さくなって
;; しまっているように見えてしまう。
;; 値が大きくなるほど、計算量全体に占める素数判定の計算量が大きくなり
;; 逆に素数判定以外の計算量が全体から見て無視出来るほどの割合になり、
;; 理論的な見積もりに近い値になっていると考えられる

問題1.23

;; 問題2.13

(define (next n)
  (if (= n 2) 3
      (+ n 2)))

(define (find-divisor n test-divisor)
  (cond [(> (square test-divisor) n) n]
	[(divides? test-divisor n) test-divisor]
	[else (find-divisor n (next test-divisor))]))


(start-test 1000)
;; 1009 *** 13000
;; 1013 *** 13000
;; 1019 *** 14000
(start-test 10000)
;; 10007 *** 39000
;; 10009 *** 35000
;; 10037 *** 34000
(start-test 100000)
;; 100003 *** 103000
;; 100019 *** 103000
;; 100043 *** 103000
(start-test 1000000)
;; 1000003 *** 367000
;; 1000033 *** 324000
;; 1000037 *** 342000

;; 問題1.22の値
;;1009 *** 27000
;;1013 *** 26000
;;1019 *** 26000
;;10007 *** 55000
;;10009 *** 55000
;;10037 *** 54000
;;100003 *** 132000
;;100019 *** 132000
;;100043 *** 132000
;;1000003 *** 416000
;;1000033 *** 430000
;;1000037 *** 423000

;;gosh> (/ 26000 13000.0)
;;2.0
;;gosh> (/ 55000 35000.0)
;;1.5714285714285714
;;gosh> (/ 132000 103000.0)
;;1.2815533980582525
;;gosh> (/ 423000 342000.0)
;;1.236842105263158

;; 小さな値の時は期待通り1/2になっているが、
;; 値が大きくなるにつれて1/2ではなくなっている
;; 値が大きくなるに連れて、nextのオーバーヘッドが
;; 1つずつ探索することと比較する効果が強くなってしまっている
;; のではないかと考えられる

問題1.24

;; 問題1.24
(define (search-for-primes start end)
  (define (timed-search-for-primes n start-time cnt)
    (cond [(or (> n end) (= cnt 3)) #t]
	  [(even? n) (timed-search-for-primes (+ n 1) (current-time) cnt)]
	  [(fast-prime? n 1)
	   (begin (display n)
		  (report-time (time-difference (current-time) start-time))
		  (timed-search-for-primes (+ n 2) (current-time) (+ cnt 1)))]
	  [else (begin (display n)
		       (newline)
		       (timed-search-for-primes (+ n 2) (current-time) cnt))]))
  (timed-search-for-primes start (current-time) 0))


(start-test 1000)
;; 1009 *** 19000
;; 1013 *** 12000
;; 1019 *** 18000
(start-test 10000)
;; 10007 *** 21000
;; 10009 *** 15000
;; 10037 *** 16000
(start-test 100000)
;; 100003 *** 23000
;; 100019 *** 17000
;; 100043 *** 16000
(start-test 1000000)
;; 1000003 *** 28000
;; 1000033 *** 18000
;; 1000037 *** 19000

;; gosh> (/ 26000 18000.0)
;; 1.4444444444444444
;; gosh> (/ 21000 18000.0)
;; 1.1666666666666667
;; gosh> (/ 17000 16000.0)
;; 1.0625
;; gosh> (/ 19000 17000.0)
;; 1.1176470588235294

;; log_2(10)は約10であるが、
;; 計算時間は10倍程度の増加となっていない
;; これは、fast-prime?では、prime?に比べて
;; 複雑な演算が多く含まれているため、
;; 数が小さい場合はそれらの影響が大きくなり
;; 数が小さい時の演算は処理に時間がかかるからと思われる
;; 数が大きくなるとこれらの処理の影響がなくなっていく
;; と考えられる

所感

ギリギリすぎて論理が破綻してる気がしなくもない