SICP 1.21-1.28

寻找因子法

(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 (divides? a b)
  (= (remainder b a) 0))

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

练习1.21

(smallest-divisor 199)
199

(smallest-divisor 1999)
1999

(smallest-divisor 19999)
7

---

练习1.22
检测prime的运行时间

(define (runtime)
        (time->seconds (current-time)))

(define (timed-prime-test n)
  (newline)
  (display n)
  (start-prime-test n (runtime)))

(define (start-prime-test n start-time)
  (if (prime? n)
      (report-prime (- (runtime) start-time))))

(define (report-prime elapsed-time)
  (display " *** ")
  (display elapsed-time))

---

解题

(define (search-for-primes start count)
  (define (search-for-primes-iter n count)
    (if (= 0 count)
        (newline)
        (if (prime? n)
            (begin (timed-prime-test n) (search-for-primes-iter (+ n 2) (- count 1)))
            (search-for-primes-iter (+ n 2) count))))

  (if (odd? start)
      (search-for-primes-iter start count)
      (search-for-primes-iter (+ start 1) count)))

> (search-for-primes 1000 3)
1009 *** 2.8848648071289062e-5
1013 *** 2.6702880859375e-5
1019 *** 2.6464462280273438e-5

> (search-for-primes 10000 3)
10007 *** 8.463859558105469e-5
10009 *** 8.225440979003906e-5
10037 *** 8.988380432128906e-5

> (search-for-primes 100000 3)
100003 *** 2.63214111328125e-4
100019 *** 2.5534629821777344e-4
100043 *** 2.7489662170410156e-4

> (search-for-primes 1000000 3)
1000003 *** 8.673667907714844e-4
1000033 *** 8.487701416015625e-4
1000037 *** 8.304119110107422e-4

> (search-for-primes 10000000 3)
10000019 *** .0018930435180664062
10000079 *** .0021266937255859375
10000103 *** .001911163330078125

---

练习1.23

(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 (next test-divisor)))))

(define (next test-divisor)
  (if (= test-divisor 2)
      3
      (+ test-divisor 2)))

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

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

> (search-for-primes 1000 3)
1009 *** 2.0265579223632812e-5
1013 *** 2.1696090698242188e-5
1019 *** 2.1696090698242188e-5

> (search-for-primes 10000 3)
10007 *** 6.079673767089844e-5
10009 *** 5.698204040527344e-5
10037 *** 5.9604644775390625e-5

> (search-for-primes 100000 3)
100003 *** 1.6808509826660156e-4
100019 *** 1.633167266845703e-4
100043 *** 1.6307830810546875e-4

> (search-for-primes 1000000 3)
1000003 *** 5.249977111816406e-4
1000033 *** 5.125999450683594e-4
1000037 *** 5.097389221191406e-4

> (search-for-primes 10000000 3)
10000019 *** .0016429424285888672
10000079 *** .0016353130340576172
10000103 *** .0010938644409179688

---

练习1.24

费马定理: 如果n是一个素数，a是小于n的任意正整数，那么a的n次方与a模n同余。
(如果两个数除以n的余数相同，称为模n同余)
0 < a < n => a$^n$ mod n = a mod n

;计算一个数的幂对另一个数取模的结果
(define (expmod base exp m)
  (cond ((= exp 0) 1)
        ((even? exp)
         (remainder (square (expmod base (/ exp 2) m))
                    m))
        (else
         (remainder (* base (expmod base (- exp 1) m))
                    m))))

(define (fermat-test n)
  (define (try-it a)
    (= (expmod a n n) a))
  (try-it (random-integer n)));random-integer是Gambit-c提供的函数

(define (fast-prime? n times)
  (cond ((= times 0) #t)
        ((fermat-test n) (fast-prime? n (- times 1)))
        (else #f)))

(define (runtime)
        (time->seconds (current-time)))

(define (timed-prime-test n)
  (newline)
  (display n)
  (start-prime-test n (runtime)))

(define (start-prime-test n start-time)
  (if (fast-prime? n 5)
      (report-prime (- (runtime) start-time))))

(define (report-prime elapsed-time)
  (display " *** ")
  (display elapsed-time))

(define (search-for-primes start count)
  (define (search-for-primes-iter n count)
    (if (= 0 count)
        (newline)
        (if (fast-prime? n 3)
            (begin (timed-prime-test n) (search-for-primes-iter (+ n 2) (- count 1)))
            (search-for-primes-iter (+ n 2) count))))

  (if (odd? start)
      (search-for-primes-iter start count)
      (search-for-primes-iter (+ start 1) count)))

> (search-for-primes 1000 3)
1009 *** 7.772445678710938e-5
1013 *** 8.20159912109375e-5
1019 *** 8.988380432128906e-5

> (search-for-primes 10000 3)
10007 *** 1.010894775390625e-4
10009 *** 9.560585021972656e-5
10037 *** 9.894371032714844e-5

> (search-for-primes 100000 3)
100003 *** 1.1563301086425781e-4
100019 *** 1.1944770812988281e-4
100043 *** 1.1801719665527344e-4

> (search-for-primes 1000000 3)
1000003 *** 1.354217529296875e-4
1000033 *** 1.3875961303710938e-4
1000037 *** 1.361370086669922e-4

> (search-for-primes 10000000 3)
10000019 *** 1.5854835510253906e-4
10000079 *** 1.671314239501953e-4
10000103 *** 1.6617774963378906e-4

---

练习1.25

(define (expmod base exp m)
  (remainder (fast-expt base exp) m))

这样会比1.24运行的慢

---

练习1.26
(expmod base (/ exp 2) m)被递归两次

---

练习1.27
能够骗过费马检查的数，称为Carmichael数，在100000000之内有255个Carmichael数，其中最小的几个是561、1105、1729、2465、2821和6601。

> (fast-prime? 561 10)
#t

> (fast-prime? 1105 5)
#t

> (fast-prime? 1729 5)
#t

> (fast-prime? 2465 5)
#t

> (fast-prime? 2821 5)
#t

> (fast-prime? 6601 5)
#t

---

练习1.28

(define (mr n)
  (define (test-all-a count result)
    (display count) (display result) (newline)
    (if result
      (if (> count 1)
        (test-all-a (- count 1) (try-it (- count 1))) ;test all a, 1 < a < n
        #t)
      #f))
  (define (try-it a)
    (= (expmod a (- n 1) n) 1))
  (test-all-a n #t)) 

(define (expmod base expo m)
  (cond ((= expo 0) 1)
        ((even? expo) (let* ((x (expmod base (/ expo 2) m)) (r (remainder (square x) m)))
                          (if (and (= r 1) (not (= x 1)) (not (= x (- m 1)))) 0 r)))
        (else (remainder (* base (expmod base (- expo 1) m))
                         m))))