Miller-Rabin_質數測試

**米勒 - 拉賓質數判定法**（英語：Miller–Rabin primality test）是一種質數判定法則，利用隨機化算法判斷一個數是合數還是_可能是_質數。1976 年，卡內基梅隆大學的計算機系教授蓋瑞・米勒首先提出了基於廣義黎曼猜想的確定性算法，由於廣義黎曼猜想並沒有被證明，於 1980 年，由以色列耶路撒冷希伯來大學的麥可・拉賓教授作出修改，提出了不依賴於該假設的隨機化算法。

## 概念

首先介紹一個相關的引理。我們發現 $1^{2}{\bmod {p}}$ 和 $(-1)^{2}{\bmod {p}}$ 總是得到 $1$ ，我們稱 $-1$ 和 $1$ 是 $1{\bmod {p}}$ 的「平凡平方根」，當 $p$ 是質數且 $p>2$ 時，不存在 $1{\bmod {p}}$ 的「非平凡平方根」。為了證明該引理，首先假設 $x$ 是 $1{\bmod {p}}$ 的平方根之一，於是有

也就是說，質數 $p$ 能夠整除 $(x-1)(x+1)$ 或者 $x=p-1$ ，根據歐幾里得引理， $x-1$ 或者 $x+1$ 能夠被 $p$ 整除，即 $x\equiv 1{\pmod {p}}$ 或 $x\equiv -1{\pmod {p}}$ ,

即 $x$ 是 $1{\bmod {p}}$ 的平凡平方根。

現在假設 $n$ 是一個質數，且 $n>2$ 。於是 $n-1$ 是一個偶數，可以被表示為 $2^{s}*d$ 的形式， $s$ 和 $d$ 都是正整數且 $d$ 是奇數。對任意在 $(Z/nZ)^{*}$ 範圍內的 $a$ ，必須滿足以下兩種形式的一種：

其中 $r$ 是某個滿足 $0\leq r\leq s-1$ 的整數。

因為由於 費馬小定理 ，對於一個質數 $n$ ，有

又由於上面的引理，如果我們不斷對 $a^{n-1}$ 取平方根後，我們總會得到 $1$ 和 $-1$ 。如果我們得到了 $-1$ ，意味著②式成立， $n$ 是一個質數。如果我們從未得到 $-1$ ，那麼通過這個過程我們已經取遍了所有 $2$ 的冪次，即①式成立。

米勒 - 拉賓質數測試就是基於上述原理的逆否，也就是說，如果我們能找到這樣一個 $a$ ，使得對任意 $0\leq r\leq s-1$ 以下兩個式子均滿足：

那麼 $n$ 就不是一個質數。這樣的 $a$ 稱為 $n$ 是合數的一個憑證（witness）。否則 $a$ 可能是一個證明 $n$ 是質數的「強偽證」（strong liar），即當 $n$ 確實是一個合數，但是對當前選取的 $a$ 來說上述兩個式子均不滿足，這時我們認為 $n$ 是基於 $a$ 的大機率質數。每個奇合數 $n$ 都有很多個對應的憑證 $a$ ，但是要生成這些 $a$ 並不容易。當前解決的辦法是使用機率性的測試。我們從集合 $(Z/nZ)^{*}$ 中隨機選擇非零數 $a$ ，然後檢測當前的 $a$ 是否是 $n$ 為合數的一個憑證。如果 $n$ 本身確實是一個合數，那麼大部分被選擇的 $a$ 都應該是 $n$ 的憑證，也即通過這個測試能檢測出 $n$ 是合數的可能性很大。然而，仍有極小機率的情況下我們找到的 $a$ 是一個「強偽證」（反而表明了 $n$ 可能是一個質數）。通過多次獨立測試不同的 $a$ ，我們能減少這種出錯的機率。

對於測試一個大數是否是質數，常見的做法隨機選取基數 $a$ ，畢竟我們並不知道憑證和偽證在 $[1,n-1]$ 這個區間如何分布。典型的例子是 Arnault 曾經給出了一個 397 位的合數 $n$ ，但是所有小於 307 的基數 $a$ 都是 $n$ 是質數的「強偽證」。不出所料，這個大合數通過了 Maple 程序的 `isprime()` 函數（被判定為質數）。這個函數通過檢測 $a=2,3,5,7,11$ 這幾種情況來進行質性檢定。

## 例子

假設我們需要檢定 $n=221$ 是否是一個質數。首先 $n-1=220=(2^{2})*55$ ，於是我們得到了 $s=2$ 和 $d=55$ . 我們隨機從 $[1,n-1]$ 中選取 $a$ ，假設 $a=174$ ，於是有：

因為我們得到了一個 -1，所以要麼 $n=221$ 確實是一個質數，要麼 $a=174$ 是一個「強偽證」。我們再選取 $a=137$ ，於是有：即 $a=137$ 是 $n=221$ 為合數的一個憑證，而 $a=174$ 是一個「強偽證」。選取特定的整數可以在一定範圍內確定（而非單純基於機率猜測）某個整數是質數還是合數。對於小於 $2^{32}$ 的情形，選取 $2,7,61$ 共三個憑據可以做到這一點；對於小於 $2^{64}$ 的情形，選取 $2,325,9375,28178,450775,9780504,1795265022$ 共七個憑據可以做到這一點。

## 算法複雜度

這一算法可以被表示成如下偽代碼：

```
**Input #1**: _n_ > 3, an odd integer to be tested for primality;
**Input #2**: _k_, a parameter that determines the accuracy of the test
**Output**: _composite_ if _n_ is composite, otherwise _probably prime_
```

```
write _n_ − 1 as 2r·_d_ with _d_ odd by factoring powers of 2 from _n_ − 1
WitnessLoop: **repeat** _k_ times:
   pick a random integer _a_ in the range [2, _n_ − 2]
   _x_ ← _a_d mod _n_
   **if** _x_ = 1 or _x_ = _n_ − 1 **then**
      **continue** WitnessLoop
   **repeat** _r_ − 1 times:
      _x_ ← _x_2 mod _n_
      **if** _x_ = _n_ − 1 **then**
         **continue** WitnessLoop
   **return** _composite_
**return** _probably prime_
```

使用模冪運算，這個算法的時間複雜度是 $O(k\log ^{3}n)$ ， $k$ 是我們測試的不同的 $a$ 的值。也就是說，這是一個高效的多項式時間算法。使用快速傅立葉轉換能夠將這個時間推進到 O (_k_ log2_n_ log log _n_ log log log _n_) = Õ(_k_ log2_n_).

如果我們加入最大公因數算法到上述算法中，我們有時候可以得到 $n$ 的因數，而不僅僅是證明 $n$ 是一個合數。例如，若 $n$ 是一個基於 $a$ 的可能的質數，但不是一個大機率質數，則 $\gcd((a^{d}{\bmod {n}})-1,n)$ 或 $\gcd((a^{2^{r}d}{\bmod {n}})-1,n)$ 將得到 $n$ 的因數。如果因式分解是必要的，這一 $GCDs$ 算法可以加入到上述的算法中，代價是增加了一些額外的運算時間。

例如，假設 $n=341$ ，則 $n-1=340=85*4$ . 於是 $2^{85}{\bmod {3}}41=32$ ，這也告訴我們 $n$ 不是一個大機率質數，即 $n$ 是一個合數。如果這個時候我們求最大公因數，我們可以得到一個 $n=341$ 的因數： $\gcd((2^{85}{\bmod {3}}41)-1,341)=31$ . 這時可行的，因為 $n=341$ 是一個基於 2 的偽質數，但不是一個「強偽質數」。

## 示例代碼

下面是根據以上定義而使用 Magma 語言編寫的「米勒 - 拉賓」檢定程序。

## 參見

* 盧卡斯 - 萊默檢定法
* 埃拉托斯特尼篩法
* 試除法
* 費馬質性檢定
* 產業等級質數

## 參考資料