简述
枚举(2sim sqrt{n})中的整数判断(n)是否为质数的方法在(n)特别大的时候就失去了其价值,( exttt{Miller-Rabin})素数检测利用了一些质数的性质在一个错误概率非常小的前提下进行判断。
费马小定理 (( exttt{Fermat's little theorem}))
当(a)不被质数(p)整除时,(a^{p-1}equiv 1 (operatorname{mod} p)),其实质是欧拉定理的一个特殊形式。但其逆定理是不成立的,所以不能单纯地说满足这一性质的数即为质数,只能说是质数的概率更大。
如果我们在([2,n-1])的范围内多次随机得到的数都满足(a^{n-1}equiv 1 (operatorname{mod} n))时,(n)就有很大的概率是一个质数。
为了进一步减小出错的概率,还需要找到其它性质来帮忙。
卡迈克尔数 (( exttt{Carmichael number}))
我们将满足(forall aperp n,a^{n-1}equiv 1 (operatorname{mod} n))的合数(n)称为卡迈克尔数,也叫做费马伪素数。
一个整数(n)为卡迈克尔数当且仅当:
- (n)无平方因子
- 对于任一质因子(p)都有(p-1|n-1)
- 至少是三个质数的乘积
证明
[click]
用反证法,假设(n=p^{alpha}m(alpha ge 2,pperp m))且为卡迈克尔数,则根据中国剩余定理,(egin{cases}
xequiv 1+p (operatorname{mod} p^{alpha})\
xequiv 1 (operatorname{mod} m)
end{cases})存在一个解(a)。显然(aperp n)。
由于(n)为卡迈克尔数,则:
又由二项式定理:
则得到了(1equiv 1+p (operatorname{mod} p^{2})),矛盾,由此得出(n)无平方因子。
由于(n)无平方因子,则(pperp n/p)。质数的指数循环节易证为(p-1)。
如果(n)为两个质数的乘积,不妨设(n=pq,p<q)。
(0<p-1<q-1,p-1
otequiv 0 (operatorname{mod} q-1))与上面所证的(q-1|n-1)矛盾。
—————————————————————————————————————————————
显然其为奇数。
且若一个数(n)为卡迈克尔数,则(2^n-1)也为卡迈克尔数,因此卡迈克尔数是有无穷多个的。证明略去(不会证),详见 OEIS A006931 。
二次探测定理
若(p)为奇素数且(alpha le 1),则方程(x^2equiv 1 (operatorname{mod} p^{alpha}))在模意义下仅有两个解,即(x=pm 1)。若一个奇数存在非平凡的解,则说明其为合数。
证明
[click]
具体实现
由卡迈克尔数的性质可知,其必不为(p^{alpha})的形态。结合费马小定理和二次探测原理即可行之有效地排除伪质数进行质数检测。
其过程大致为:
设检测的数为(n),将(n-1)分解为(2^t u,2
otmid u)的形式。选取(s)作为检验的次数,在这(s)次中每次选择一个属于([2, n-1])的数(a),先算出(a^u)模(n)的值,再将它平方(t)次。若在平方的过程中发现了模意义下的1的非平凡平方根,即立刻判断为合数,或平方完得到(a^{n-1}
ot equiv 1 (operatorname{mod}))也说明这是合数;若一直平方到(a^{n-1})都没有出现上述情况,则此次判断中视为该数为质数。当然,通过一次判断并不代表这就是质数了,只能说是质数的可能性更大。
若(n)是一个奇合数,则在上述过程中能判断(n)为合数的数的个数至少为((n-1)/2);对于(forall n,s),( exttt{Miller-Rabin})检测出错的概率至多为(2^{-s})。证明略(太麻烦)。
所以只需要一个适合的判断次数以及一个靠谱的随机数生成方式,就可以将判断出错的概率控制在一个很小很小的范围内。经大量的研究表明,对几乎所有可以想象到的应用选取(s=50)应该都是足够的。
[click]
#include <cstdio>
#include <cctype>
#include <ctime>
#include <cstdlib>
typedef long long ll;
typedef long double ld;
ll random(ll n) {return (ll)rand()*rand()%n+1;}
ll mul(ll a,ll b,ll p) {return ((a*b-(ll)((ld)a/p*b)*p+p)%p+p)%p;}
//ll mul(ll a,ll n,ll p)
//{
// ll res=0;
// while(n)
// {
// if (n&1)
// res=(res+a)%p;
// a=(a+a)%p;
// n>>=1;
// }
// return res;
//}
ll power(ll a,ll n,ll p)
{
ll res=1;
while(n)
{
if (n&1)
res=mul(res, a, p);
a=mul(a, a, p);
n>>=1;
}
return res;
}
bool millerRabin(ll n)
{
if (n<3)
return n==2;
int t=0;
ll u=n-1;
while(!(u&1))
u>>=1,t++;
if (!t)
return 0;
for (int i=1;i<=50;i++)
{
ll a=power(random(n-1), u, n),pre=a;
if (a==1)
continue;
for (int j=1;j<=t;j++)
{
a=mul(a, a, n);
if (a==1)
{
if (pre!=1&&pre!=n-1)
return 0;
break;
}
pre=a;
}
if (a!=1)
return 0;
}
return 1;
}
int main()
{
srand((unsigned)time(0));
ll x;
while(scanf("%lld",&x)!=EOF)
printf("%c
",millerRabin(x)?'Y':'N');
return 0;
}
相关练习
luoguP3653 小清新数学题
参考
- CLRS solutions to introdution to algorithms-chapter_31
- OI WIKI
- 《算法导论(原书第3版)》 Thomas H.Cormen Charles E.Leiserson/Ronald L.Rivest Clifford Stein 著