米勒-拉賓質數判定法(英語:Miller–Rabin primality test)是一种質數判定法則,利用随机化算法判断一个数是合数还是可能是素数。1976年,卡内基梅隆大学的计算机系教授蓋瑞·米勒首先提出了基于广义黎曼猜想的确定性算法,由于广义黎曼猜想并没有被证明,於1980年,由以色列耶路撒冷希伯來大學的迈克尔·拉宾教授作出修改,提出了不依赖于该假设的随机化算法。
首先介绍一个相关的引理。我们发现
和
总是得到
,我们称
和
是
的“平凡平方根”,当
是素数且
时,不存在
的“非平凡平方根”。为了证明该引理,首先假设
是
的平方根之一,于是有
![{\displaystyle x^{2}\equiv 1{\pmod {p}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/f2404e3c31b0cfd08dda220c8baa8219babc5f97)
![{\displaystyle (x+1)(x-1)\equiv 0{\pmod {p}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/672a2fc37298d72e29c68b7013b486b6246c38db)
也就是说,素数
能够整除
或者
,根据欧几里得引理,
或者
能够被
整除,即
或
,
即
是
的平凡平方根。
现在假设
是一个素数,且
。于是
是一个偶数,可以被表示为
的形式,
和
都是正整数且
是奇数。对任意在
范围内的
,必须满足以下两种形式的一种:
![{\displaystyle a^{d}\equiv 1{\pmod {n}}{\text{ ① }}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/66ae04641adf19b3e11caaee8cdad13ef9165225)
![{\displaystyle a^{2^{r}d}\equiv -1{\pmod {n}}{\text{ ② }}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/a19331cde958ce0bea3027b18996863995c8b525)
其中
是某个满足
的整数。
因为由于 费马小定理 ,对于一个素数
,有
![{\displaystyle a^{n-1}\equiv 1{\pmod {n}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/3bfc8625369c7558deadd61823842db06983423d)
又由于上面的引理,如果我们不断对
取平方根后,我们总会得到
和
。如果我们得到了
,意味着②式成立,
是一个素数。如果我们从未得到
,那么通过这个过程我们已经取遍了所有
的幂次,即①式成立。
米勒-拉宾素性测试就是基于上述原理的逆否,也就是说,如果我们能找到这样一个
,使得对任意
以下两个式子均满足:
![{\displaystyle a^{d}\not \equiv 1{\pmod {n}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/0867bd3e86262241fae38a021ed361001e28f66a)
![{\displaystyle a^{2^{r}d}\not \equiv -1{\pmod {n}}}](https://wikimedia.org/api/rest_v1/media/math/render/svg/e20bbbf10f3ef34c4e8d4e398dd978250331af54)
那么
![{\displaystyle n}](https://wikimedia.org/api/rest_v1/media/math/render/svg/a601995d55609f2d9f5e233e36fbe9ea26011b3b)
就不是一个素数。这样的
![{\displaystyle a}](https://wikimedia.org/api/rest_v1/media/math/render/svg/ffd2487510aa438433a2579450ab2b3d557e5edc)
称为
![{\displaystyle n}](https://wikimedia.org/api/rest_v1/media/math/render/svg/a601995d55609f2d9f5e233e36fbe9ea26011b3b)
是合数的一个凭证(witness)。否则
![{\displaystyle a}](https://wikimedia.org/api/rest_v1/media/math/render/svg/ffd2487510aa438433a2579450ab2b3d557e5edc)
可能是一个证明
![{\displaystyle n}](https://wikimedia.org/api/rest_v1/media/math/render/svg/a601995d55609f2d9f5e233e36fbe9ea26011b3b)
是素数的“强伪证”(strong liar),即当
![{\displaystyle n}](https://wikimedia.org/api/rest_v1/media/math/render/svg/a601995d55609f2d9f5e233e36fbe9ea26011b3b)
确实是一个合数,但是对当前选取的
![{\displaystyle a}](https://wikimedia.org/api/rest_v1/media/math/render/svg/ffd2487510aa438433a2579450ab2b3d557e5edc)
来说上述两个式子均不满足,这时我们认为
![{\displaystyle n}](https://wikimedia.org/api/rest_v1/media/math/render/svg/a601995d55609f2d9f5e233e36fbe9ea26011b3b)
是基于
![{\displaystyle a}](https://wikimedia.org/api/rest_v1/media/math/render/svg/ffd2487510aa438433a2579450ab2b3d557e5edc)
的大概率素数。
每个奇合数
都有很多个对应的凭证
,但是要生成这些
并不容易。当前解决的办法是使用概率性的测试。我们从集合
中随机选择非零数
,然后检测当前的
是否是
为合数的一个凭证。如果
本身确实是一个合数,那么大部分被选择的
都应该是
的凭证,也即通过这个测试能检测出
是合数的可能性很大。然而,仍有极小概率的情况下我们找到的
是一个“强伪证”(反而表明了
可能是一个素数)。通过多次独立测试不同的
,我们能减少这种出错的概率。
对于测试一个大数是否是素数,常见的做法随机选取基数
,毕竟我们并不知道凭证和伪证在
这个区间如何分布。典型的例子是 Arnault 曾经给出了一个397位的合数
,但是所有小于307的基数
都是
是素数的“强伪证”。不出所料,这个大合数通过了 Maple 程序的isprime()
函数(被判定为素数)。这个函数通过检测
这几种情况来进行素性检验。
假设我们需要检验
是否是一个素数。首先
,于是我们得到了
和
.我们随机从
中选取
,假设
,于是有:
![{\displaystyle a^{2^{0}d}{\bmod {n}}=174^{55}{\bmod {2}}21=47\not =1,-1}](https://wikimedia.org/api/rest_v1/media/math/render/svg/19e9e9a882fac5f476786c36d8c3f4a6bcbaf360)
![{\displaystyle a^{2^{1}d}{\bmod {n}}=174^{110}{\bmod {2}}21=220=-1}](https://wikimedia.org/api/rest_v1/media/math/render/svg/8301ee7d2044aced56588bc94da1e173f8b65148)
因为我们得到了一个 -1,所以要么
![{\displaystyle n=221}](https://wikimedia.org/api/rest_v1/media/math/render/svg/8ddd2d6a8ddcb65d4097310c9aad8cbc5e7e609a)
确实是一个素数,要么
![{\displaystyle a=174}](https://wikimedia.org/api/rest_v1/media/math/render/svg/d2aa71c10280ec3e806c0d38327bafdf4f3cb2e4)
是一个“强伪证”。我们再选取
![{\displaystyle a=137}](https://wikimedia.org/api/rest_v1/media/math/render/svg/8f24189d651acb67bd6c8a374305b2a655432815)
,于是有:
![{\displaystyle a^{2^{0}d}{\bmod {n}}=137^{55}{\bmod {2}}21=188\not =1,-1}](https://wikimedia.org/api/rest_v1/media/math/render/svg/de7341b95c3e6c251dd3fecf82b518923981175f)
![{\displaystyle a^{2^{1}d}{\bmod {n}}=137^{110}{\bmod {2}}21=205\not =-1}](https://wikimedia.org/api/rest_v1/media/math/render/svg/afb910be8ecbc2d88d5972c6268a4ab1c5202727)
即
![{\displaystyle a=137}](https://wikimedia.org/api/rest_v1/media/math/render/svg/8f24189d651acb67bd6c8a374305b2a655432815)
是
![{\displaystyle n=221}](https://wikimedia.org/api/rest_v1/media/math/render/svg/8ddd2d6a8ddcb65d4097310c9aad8cbc5e7e609a)
为合数的一个凭证,而
![{\displaystyle a=174}](https://wikimedia.org/api/rest_v1/media/math/render/svg/d2aa71c10280ec3e806c0d38327bafdf4f3cb2e4)
是一个“强伪证”。
选取特定的整数可以在一定范围内确定(而非单纯基于概率猜测)某个整数是质数还是合数。对于小于
的情形,选取
共三个凭据可以做到这一点;对于小于
的情形,选取
共七个凭据可以做到这一点[1]。
算法复杂度[编辑]
这一算法可以被表示成如下伪代码:
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 ← ad mod n
if x = 1 or x = n − 1 then
continue WitnessLoop
repeat r − 1 times:
x ← x2 mod n
if x = n − 1 then
continue WitnessLoop
return composite
return probably prime
使用模幂运算,这个算法的时间复杂度是
,
是我们测试的不同的
的值。也就是说,这是一个高效的多项式时间算法。使用快速傅里叶变换能够将这个时间推进到 O(k log2n log log n log log log n) = Õ(k log2n).
如果我们加入最大公因数算法到上述算法中,我们有时候可以得到
的因数,而不仅仅是证明
是一个合数。例如,若
是一个基于
的可能的素数,但不是一个大概率素数,则
或
将得到
的因数。如果因式分解是必要的,这一
算法可以加入到上述的算法中,代价是增加了一些额外的运算时间。
例如,假设
,则
.于是
,这也告诉我们
不是一个大概率素数,即
是一个合数。如果这个时候我们求最大公因数,我们可以得到一个
的因数:
.这时可行的,因为
是一个基于2的伪素数,但不是一个“强伪素数”。
示例代码[编辑]
下面是根据以上定义而使用Magma语言编写的“米勒-拉宾”检验程序。
function ModPotenz(a,b,n)
erg:=1;
while (b ne 0) do
b, rest:=Quotrem(b,2);
if (rest ne 0) then erg:=((erg*a) mod n); end if;
a:= (a^2 mod n);
end while;
return erg;
end function;
;
function Miller_rabin(n)
if (n mod 2 ne 0) then
d:=n-1; s:=0;t:=50;
while (d mod 2) eq 0 do
d:=d div 2;
s:=s+1;
end while;
k:=0;
while(k lt t) do
a:=Random(1,n-1);
k:=k+1;
if GCD(a,n) ne 1 then
continue;
end if;
x:=ModPotenz(a,d,n);
if(x ne 1) then
for r in [0,s-1] do
x:=ModPotenz(a,2^r*d,n);
if(x eq (n-1)) then
return "probably prime";
end if;
end for;
elif (x eq 1) then
break;
end if;
end while;
end if;
return "composite";
end function;
參考資料[编辑]