Skip to content →

ミラー・ラビン素数判定

与えられた数 \(n\) が素数かどうかを 誤り率 \(4^{-k}\)、時間 \(O(\log(n) k)\) で判定する、ミラー・ラビン素数判定 (Miller-Rabin primality test) を紹介します(計算量は \(n\) が word size に収まるとした場合)。

フェルマーテスト

\(n\) が素数のとき、フェルマーの小定理より任意の \(x \neq 0\) について \(x^{n-1} = 1 \bmod n\) が成り立ちます。従って、ランダムな \(x \in [1,n)\) について、この等式が必ず成り立つならば \(n\) が素数であると高い確率で言えそうです。また、この等式を満たさない \(x\) が存在すれば、\(n\) は合成数です。この素数判定法をフェルマーテストと言います。

しかし、残念ながら合成数であるにも関わらず、フェルマーテストが高い確率で失敗する数が存在します。\(n\) と互いに素な任意の \(a\) について \(x^{n-1} = 1 \bmod n\) となるような \(n\) が無限に存在し、カーマイケル数といいます。

最小のカーマイケル数は \(n=561=3 \cdot 11 \cdot 17\) です。カーマイケル数であることは、中国剰余定理より各素因数に対して

  • \(561=2 \cdot 280 + 1 = 1\bmod 2\)
  • \(561=10 \cdot 56 + 1 = 1\bmod 10\)
  • \(561=16 \cdot 35 + 1 = 1\bmod 16\)

が成り立つことから従います。

\(n\) と互いに素でない数 \(a\) をフェルマーテストに掛けると、必ず \(n\) が合成数だと判明します。例えば、\(n=p_1p_2p_3\) がカーマイケル数ならば \((1-1/p_1)(1-1/p_2)(1-1/p_3)\) の確率で合成数だと分かりません。従って、素数と判明するまでには \(p_1\) 回程度試行する必要があり \(n\) の多項式時間になってしまいます。そこで、ミラー・ラビンの素数判定法へと改良します。

ミラー・ラビン素数判定法

素数 \(p\) ならば、\(x^2=1 \bmod p\) は高々 2 個しか解を持たず、\(1\) の平方根は \(\pm 1\) のみです。一方 \(n\) が合成数の場合は複数の平方根を持ち得ます。例えば、カーマイケル数 \(n=561\) を mod とする 1 の平方根は \(x=1,67,188,254,307,373,494,560\) の 8 つあります。

\(n-1=q2^Q\) (\(q\) : 奇数) とします。フェルマーテストで合成数と確定しなかったランダムな数 \(a\) (つまり \(a^{n-1}=1\))に対して、
$$a^{q}, a^{q2^1}, a^{q2^2}\ldots, a^{q2^{Q-1}}$$に \(\pm 1\) 以外の \(1\) の平方根が含まれていないか調べるのがミラー・ラビン素数判定法です。

ミラー・ラビン素数判定法ではランダムな \(a\) に対して、\(1-4^{-1}\) の確率で \(n\) が合成数だと確定させることができます(証明参照)。\(k\) 個のランダムな数に対して行えば、確率 \(1-4^{-k}\) です。計算量は一回当たり \(O(\log(n))\) です。

誤り率の証明

偶数の素数は \(2\) のみですから、\(n=\prod_{i=1}^m p_i^{e_i}\) が奇数と仮定しても一般性を失いません。

\(a \not \in (\mathbb{Z}/n\mathbb{Z})^\times\) は \(\gcd(a,n) \neq 1\) より \(a^{n-1} \neq 1 \bmod n\) です。従って、\(n\) が合成数だと判定できます。

\((\mathbb{Z}/n\mathbb{Z})^\times\) の元を考えましょう。

\(x^{n-1}=1 \pmod n\) を満たす \(x\) からなる部分集合 \(H\) は \((\mathbb{Z}/n\mathbb{Z})^\times\) の部分群をなします。ラグランジュの定理より、\((\mathbb{Z}/n\mathbb{Z})^\times\) は互いに共通部分を持たない \(H,g_1H,\ldots,g_sH\) の和集合$$(\mathbb{Z}/n\mathbb{Z})^\times = H \cup g_1 H \cup \cdots \cup g_s H$$ で表せます。\(a \not \in H\) だと、\(a^{n-1} \neq 1\) より \(n\) が合成数だと判定できます。\(n-1=q2^Q\) (\(q\) : 奇数) とします。\(x^q=1\) を満たす \(x\) からなる部分集合 \(F \subset H\) も部分群をなし、互いに共通部分を持たない \(F,a_1F,\ldots,a_tF\) の和集合 $$ H = F \cup a_1 F \cup \cdots \cup a_t F$$ で表せます。\(q\) が奇数より \(b \in F \Rightarrow -b \not\in F\) にもかかわらず \(-b \in H\) です。中国剰余定理で素冪の mod にバラシて適用することで、\(t \geq 2^m\) が得られます(\(n=\prod_{i=1}^m p_i^{e_i}\))。つまり \(|H| \geq 2^m|F|\) です。

\(H\) を次の3つに分割します。

  • \(F = \{x : x^q=1\}\)
  • \(F_1 = \{x : \) \(x^{q2^s}=-1, x^{q2^{s+1}}=1\) なる \(s\) が存在\(\}\)
  • \(F_2 = \{x : \) \(x^{q2^s}\neq \pm1, x^{q2^{s+1}}=1\) なる \(s\) が存在\(\}\)

\(F,F_1,F_2\) のうち \(a \in F_2\) の場合だけ \(a\) が合成数だと判定できます。\(F,F_1\) の全体に占める割合を知りたいです。中国剰余定理より、素因数分解 \(n=\prod_{i=1}^m p_i^{e_i}\) の各素冪 \(p_i^{e_i}\) の mod に分解できます。つまり、\(x \in F_1\) ならば $$x^{q2^s}=-1 \bmod n \Leftrightarrow \begin{cases}
x^{q2^s}=-1 \bmod p_1^{e_1} \\
x^{q2^s}=-1 \bmod p_2^{e_2} \\
\vdots \\
x^{q2^s}=-1 \bmod p_m^{e_m} \\
\end{cases}$$ です。各素冪に対する mod を独立に 1 にすることができるので、\(|F_2| \geq (2^m-1)|F_1| \) です。\(2^m|F|\leq |H|\) と併せて
$$\begin{align}
|F|+|F_1| \leq& |F| + (|H|-|F|)/(2^m-1) \\
\leq& |H|/(2^m-1)+|F|(2^m-2)/(2^m-1) \\
=& \frac{2-2^{-(m-1)}}{2^m-1} |H|\\
\end{align}$$を得ます。\(m\geq3\) のとき \(|F|+|F_1| \leq 4^{-1}|H|\)、\(m=2\) のとき \(|F|+|F_1| \leq 2^{-1}|H|\) です。

(i) 2つの素数の積 \(n=p^2\) の場合(\(m=1\)):
\(n\) が平方因子 \(p^2\) を持つ場合、\(a_k=1+kp\) (\(k \neq 0\):原始根)とすると $$\begin{align}a_k^{n-1} &= 1 + k(n-1)p \bmod p^2\\ &=1-kp\bmod p^2\end{align}$$ です。\(k=1,2,\ldots,p-1\) それぞれに対して、\(a_k\) が異なる \(g_iH\) に属するので、\(s \geq p-1\) です。\(p=3\) のとき、テストをすり抜けるのは \(x=1\) のみなので失敗確率 \(1/8\) です。\(p \geq 5\) のとき、\(s \geq 4\) より失敗確率 \(5^{-1}\) 以下です。

(ii) 相異なる素数の積 \(n = pq\) の場合(\(m=2\)):
\(n\) が相異なる素数の積 \(n=pq\) であるとします( \(p < q\) )。\(a^{n-1}=a^{pq-1}=a^{q-1}\bmod p\) ですが、\(x^{q-1} = 1 \bmod p\) は解を \(\gcd(q-1,p-1)\) 個しか持たないので \(s \geq 1\) です。よって、失敗確率 \(4^{-1}\) 以下です。

(iii) \(n\) が3つ以上の素因数を持つ場合(\(m\geq3\)):
失敗確率 \(4^{-1}\) 以下です。

実装

C++の実装は次の通りです:

#include <iostream>
#include <random>
using namespace std;

long long pow(long long a, long long n, long long p) {
  return n!=0?pow(a*a%p,n/2,p)*(n%2==1?a:1)%p:1;
}

// n が素数か判定
bool miller_rabin(long long n) {
  if (n % 2 == 0) return n == 2;
  long long q = n, Q = 0;
  while (q % 2 == 0) {
    q/=2;
    ++Q;
  }
  
  random_device seed_gen;
  default_random_engine engine(seed_gen());
  uniform_int_distribution<> dist(1, n - 1);
  for (int i = 0; i < 20; ++i) {
    long long a = dist(engine);
    if (pow(a, n - 1, n) != 1) {
      return false;
    }
    a = pow(a, q, n);
    for (int j = 0; j < Q; ++j) {
      long long na = a * a % n;
      if (na == 1 && a == n - 1) return false;
      a = na;
    }
    
  }
  return true;
}

int main() {
  for (int i=2;i<100;++i)
  cout<<i<<" "<<miller_rabin(i)<<endl;
}

参考文献

  1. はじめての数論 原著第3版 発見と証明の大航海‐ピタゴラスの定理から楕円曲線まで
  2. Bobby Kleinberg氏の授業資料
フェルマーテストだと不十分で、1 の平方根を調べると十分なのは、何か上手く説明できるのでしょうか。今のところ、偶然にしか思えません汗。

Published in データ構造

Comments

コメントを残す

メールアドレスが公開されることはありません。 * が付いている欄は必須項目です

CAPTCHA