要約：素数判定に使われるミラーラビン法を解説しながら、Haskell で実装してみる。

フェルマーテスト 大きな数を確実に素数だと判定するには、大変時間がかかるので、実用的には「ほぼ素数だ」と確率的に判定する。確率的な素数判定の代表格がフェルマーテストである。 フェルマーテストには、以下に示すフェルマーの小定理を利用する。 a^p ≡ a (mod p) a は任意の整数。p は素数である。法 p の下で a を p 乗したものは、a と合同であると言う意味だ。a には制限はないが、特に a を p より小さい整数、0 ≦ a ≦ p - 1 とすれば、a を p 乗して、p で割ると a に戻るとも解釈できる。 最初に見たときは、だからどうしたと思われるかもしれない。しかし、有名なフェルマーの大定理が実用上何の役にも立たないのに対し、フェルマーの小定理はいろんな場面で活躍する。 実際に計算してみよう。準備として、base を ex 乗して m で割った余りをとる関数 expmod を以下のように定義する。この関数は、高速で、しかも使用するメモリーが少ない。 expmod :: Integer -> Integer -> Integer -> Integer expmod base ex m | ex == 0 = 1 | even ex = square (expmod base (ex `div` 2 ) m) `mod` m | otherwise = (base * (expmod base (ex - 1 ) m)) `mod` m where square x = x * x 素数の例として 563 を取り上げる。 expmod 2 563 563 → 2 map ( \ x -> expmod x 563 563 ) [ 0 .. 562 ] == [ 0 .. 562 ] → True すべての素数はフェルマーの小定理を満たす。よって大きな整数 n が素数であるかを判定するには、任意の整数 a に対して n がフェルマーの小定理を満たすか調べればよい。 フェルマーの小定理を使って、n が素数か確率的に判定する fermat という関数は以下のように書ける。 fermat :: Integer -> Integer -> Bool fermat n a = expmod a n n == a 素数に対して実験してみよう。 fermat 563 2 → True fermat 563 3 → True fermat 563 4 → True 素数であると判定された。 次に、合成数 560 について実験してみる。 fermat 560 2 → False fermat 560 3 → False fermat 560 4 → False 合成であると判定された。 乱数を自動的に生成して、テストするようにしてみよう。 import Random import Control.Applicative fermatRandom :: Integer -> IO Bool fermatRandom n = fermat n <$> getRandom 2 (n - 2 ) 実際に使ってみる。 fermatRandom 563 → True fermatRandom 560 → False Bool が返ってきているように見えるが、乱数を使っているので IO Bool が返っていることに注意。

疑素数 フェルマーテストを騙す数を疑素数という。たとえば、合成数 341 (11 × 31) は、乱数がたまたま 2 となれば、フェルマーテストを騙すことができる。 >||haskell|

fermat 341 2

→ True

|

つまり、こういうことだ。

整数 n が、フェルマーテストに落ちたら、n は確実に合成数

整数 n が、フェルマーテストに受かったら、n は素数か合成数

1回のテストでは心もとないので、通常は k 回ほどこのテストを繰り返して確率を上げる。

import Control.Monad fermatTest :: Int -> Integer -> IO Bool fermatTest times n = and <$> replicateM times (fermatRandom n)

341 に対して、10 回テストを繰り返してみよう。

fermatTest 10 341 → False

カーマイケル数 最強の疑素数をカーマイケル数という。 疑素数 341 が 2 のときにしかフェイルマーテストを騙せなかったのに対し、疑素数はフェルマーテストを騙す。 最小のカーマイケル数は、561 (3 × 11 × 17) である。実際に試してみよう。 fermat 561 2 → True fermat 561 3 → True fermat 561 4 → True length $ filter (fermat 561 ) [ 0 .. 560 ] → 561 恐ろしいことにすべての場合で、騙していると分かる。

フェルマーの小定理の変形 カーマイケル数に騙されにくくする方法の1つは、変形したフェルマーの小定理を用いることである。 フェルマーの小定理は、以下の通りであった。 a^p ≡ a (mod p) a が p の倍数でなければ、両辺を a で割ってよいので、以下のように変形できる。(法 p の下では、p の倍数は 0 と同じ意味になるので、これで割ってはいけない。) a^(p - 1)≡ 1 (mod p) 我々にとっては、0 を取り除いた、1 ≦ a ≦ p - 1 の範囲で考えれば十分である。変形フェルマーの小定理を使ったテストは以下のように実装できる。 fermat' :: Integer -> Integer -> Bool fermat' n a = expmod a (n - 1 ) n == 1 実際に使ってみよう。 fermat' 561 2 → True fermat' 561 3 → False fermat' 561 4 → True カーマイケル数 n が変形フェルマーテストを騙す a の条件は、gcd(n,a) = 1 である。 たとえば、561 の場合、3、11、17 を因数に持つ a は、変形フェルマーテストを騙せない。 length $ filter (fermat' 561 ) [ 1 .. 561 ] → 320 フェルマーテストを変形すると、騙されにくくなったことが分かる。

ミラーラビン法のアイディア 変形フェルマーテストの確率をさらに上げるために、自明でない 1 の平方根を探す方法がある。これは、ミラーラビン法と呼ばれている。 まず、1 の自明な平方根とはなんだろうか？ その数 x は二乗したら 1 になるはずである。 x^2 = 1 (mod n) 変形すると (x - 1)(x + 1) = 0 (mod n) n が素数であれば、n は (x - 1) も (x + 1) も割り切らない。そこで、x は 1 か -1 であると分かる。これを 1 の自明な平方根という。-1 とは、n - 1 のことである。 合成数の場合は、1 や -1 以外に、二乗して 1 となる整数が存在する。その数を 1 の自明でない平方根と呼ぶ。 意味が分からないと思うので、例を挙げよう。 素数 563 に対して、自明な平方根は、1 と 562 である。 expmod 1 2 563 → 1 expmod 562 2 563 → 1 たとえば、67 は二乗しても 1 にはならない。 expmod 67 2 563 → 548 ここで、カーマイケル数 561 に対しても同じ計算をしてみる。 expmod 1 2 561 → 1 expmod 560 2 561 → 1 expmod 67 2 561 → 1 67 を 2 乗すると 1 になった。これが1の自明でない平方根である。つまり、66 と 68 を素因数分解してみれば、3、11、17 という因数が現れる！ 僕が計算したところ、561 に対して1の自明でない平方根は 67,188,254,307,373,494 であった。

ミラーラビン法のアルゴリズム ミラーラビン法は、変形フェルマーテストを計算する過程で、1の自明でない平方根を見つけるアルゴリズムである。 素数か合成数か判定したい整数 n は奇数である。(n = 2 は自明すぎて対象外。それ以外の偶数は、もちろん合成数である。) よって、n - 1 は偶数であり、因数 2 を持つ。 n - 1 を d × 2^s と表すことにする。d は奇数である。変形フェルマーの小定理より。 a^(d × 2^s) ≡ 1 (mod n) つまり、n が素数であれば、a を d × 2^s 乗したときは 1。素数は、1 の自明な平方根しか持たないので、a を d × 2^(s-1) 乗したときは、1 か -1 である。 a を d × 2^(s-1) 乗したとき 1 となった場合、a を d × 2^(s-2) 乗したときは、1 か -1 となる。 同様の考察を続けると、n が素数であるための条件は以下の通り。 a^d ≡ 1 (mod n)、または

ある r (0 ≦ r ＜ s - 1) に対して a^(d × 2^r) ≡ -1 (mod n) これを実装すると、以下のようになる。 factor2 :: Integer -> (Int, Integer) factor2 n = factor2' ( 0 ,n) where factor2' (s,d) | even d = factor2' (s + 1 , d `div` 2 ) | otherwise = (s,d) millerRabin :: Integer -> Integer -> Bool millerRabin n a = let (s,d) = factor2 (n - 1 ) aExpD = expmod a d n double x = expmod x 2 n evens = take s $ iterate double aExpD in aExpD == 1 || any ( == n - 1 ) evens millerRabinRandom :: Integer -> IO Bool millerRabinRandom n = millerRabin n <$> getRandom 2 (n - 2 ) millerRabinTest :: Int -> Integer -> IO Bool millerRabinTest times n = and <$> replicateM times (millerRabinRandom n) 実際に実験してみよう。 length $ filter (millerRabin 561 ) [ 1 .. 560 ] → 10 fermat' では 320 通り騙されていたのに、millerRabin では 10 通りしか騙されない。つまり、310個の場合の計算過程において、1の自明でない平方根(67,188,254,307,373,494)を見つけたことになる。 もちろん素数は、すべての場合に millerRabin テストにパスする。 length $ filter (millerRabin 563 ) [ 1 .. 562 ] → 562

最終実験 以下のようにカーマイケル数と素数を定義する。 carmichaels :: [Integer] carmichaels = [ 561 , 1105 , 1729 , 2465 , 2821 , 6601 , 8911 , 10585 , 15841 ] primes :: [Integer] primes = [ 563 , 1009 , 1013 , 1019 , 10007 , 10009 , 10037 , 100003 , 100019 ] カーマイケル数は、高い確率で fermatTest を騙す。 mapM (fermatTest 10 ) carmichaels → [True,True,True,True,True,True,True,True,True] mapM (fermatTest 10 ) primes → [True,True,True,True,True,True,True,True,True] 一方、カーマイケル数は、millerRabinTest を騙すこともあるが、その確率は高くない。 mapM (millerRabinTest 10 ) carmichaels → [False,False,False,False,False,False,False,False,False] mapM (millerRabinTest 1 ) carmichaels → [True,False,False,False,False,False,False,False,False] mapM (millerRabinTest 10 ) primes → [True,True,True,True,True,True,True,True,True]