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

ochalog

Ruby と MediaWiki が好きな電子・情報系の学生のブログ。

SICP: Exercise 1.28

問題 1.28

だまされない Fermat テストの変形の一つは Miller-Rabin テスト(Miller-Rabin test)(Miller 1976; Rabin 1980)という。これは、\( n \) が素数であり、\( a \) を \( n \) より小さい任意の正の整数とすると、\( a \) の \( (n - 1) \) 乗は \( n \) を法として 1 と合同であるという、Fermat の小定理のもう一つの方から始める。Miller-Rabin テストで数 \( n \) の素数性をテストするには、\( a < n \) をランダムにとり、expmod 手続きを使って \( n \) を法とした \( a \) の \( (n - 1) \) 乗を作る。しかし、expmod で二乗を計算しながら「\( n \) を法とした 1 の自明でない平方根」、つまり 1 でも \( n - 1 \) でもないがその二乗が \( n \) を法として 1 になる数がなかったか調べる。そういう、\( n \) を法とした 1 の自明でない平方根があれば \( n \) が素数でないことは証明出来る。また \( n \) が素数でない奇数なら、\( a < n \) の少なくても半分は \( a^n-1 \) を計算する時、\( n \) を法とした 1 の自明でない平方根が出てくることも証明出来る。(これが Miller-Rabin テストがだまされない理由である。)\( n \) を法とした 1 の自明でない平方根を見つけたらシグナルを出すように expmod 手続きを修正せよ。これを使い fermat-test に似た Miller-Rabin テストの手続きを実装せよ。分っている素数、非素数をテストし、手続きを検査せよ。ヒント:expmod にシグナルを出させる方法の一つは 0 を返させることである。

計算機プログラムの構造と解釈 第二版 1.2.6 例:素数性のテスト」より

素数性のテスト」の最後の問題は、改良された素数性テストを実装するもの。問題文は長いが、丁寧に説明されているので、それに従えば素直に実装することができる。

「1 の自明でない平方根」の話は、Wikipedia に載っていた証明が分かりやすかった。

まず、有限体 \( \mathbb{Z}_p \) の単位元平方根についての補題を考える。ここで、\( p \) は奇素数である。\( p \) を法とした剰余\( \pmod{p} \) において、1 と -1 の二乗は常に 1 となる。これらを 1 の「自明な」平方根と呼ぶ。\( p \) は素数なので 1 の「自明でない」平方根は存在しない。これを示すため、1 mod \( p \) の非自明な平方根を \( x \) とする。このとき、以下が成り立つ: $$ x^2 \equiv 1\pmod{p} $$ $$ x^{2} - 1 \equiv 0\pmod{p} $$ $$ (x - 1)(x + 1) \equiv 0\pmod{p} $$ \( p \) は素数なので、これは \( x - 1 \) または \( x + 1 \) が \( p \) で割り切れることを示す。一方、\( x \) は 1 でも -1 でもないので\( \pmod{p} \)、\( x - 1 \) も \( x + 1 \) も \( p \) で割り切れない。すなわち、二乗して \( 1 \pmod{p}\) となるのは、1 および -1 だけということになる。

ミラー-ラビン素数判定法 - Wikipedia」より

解答

(define (miller-rabin-test n)
  (define (expmod base exp m)
    ;; 1 の非自明な平方根かどうかを調べる
    (define (check-nontrivial-sqrt-1 x)
      (cond ((= x 1) 1)       ; 1 の自明な平方根
            ((= x (- m 1)) 1) ; 1 の自明な平方根
            (else
              (let ((r (remainder (square x) m)))
                (if (= r 1)
                  0 ; 1 の非自明な平方根であることを示すシグナル
                  r)))))
    (cond ((= exp 0) 1)
          ((even? exp) (check-nontrivial-sqrt-1
                         (expmod base (/ exp 2) m)))
          (else (remainder (* base (expmod base (- exp 1) m))
                           m))))
  (define (try-it a)
    (= (expmod a (- n 1) n) 1))
  (try-it (+ 1 (random (- n 1)))))

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

(define (prime? n)
  ;; Miller-Rabin テストを 20 回行うので誤りの確率は
  ;; (1/4)^20 ≈ 10^(-12)
  (fast-prime? n 20))

テスト

素数、非素数でテストせよとのことなので、とりあえず 100 未満の数から素数を選んでみる。

(use srfi-1)
(define numbers (iota 99 2)) ; 2〜100 までの整数のリスト
(filter prime? numbers)
; (2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97)

何回か試したが、結果は 100 未満の素数のリストになり、大丈夫そう。

問題 1.27 で扱った Carmichael 数はどうだろうか。

ochaochaocha3.hateblo.jp

(define carmichael-numbers '(561 1105 1729 2465 2821 6601))
(map prime? carmichael-numbers)
; (#f #f #f #f #f #f)

Carmichael 数にもだまされていない! したがって、Fermat テストを改良することができた。