Miller-Rabin
事先声明,因为菜鸡Hastin知识水平有限就是菜,因此语言可能不是特别规范,仅供理解.
step 0
问一个数(p)是否为质数,(p<=10^{18}).
一个简单暴力的办法是(O( sqrt{n}))枚举约数.
然而显然会炸.
于是我们就有了Miller-Rabin.
讲了好多废话...
step 1
首先了解一下费马小定理:
若(p)为质数,则对于(ain[1,p-1])有(a^{p}equiv a(Mod) (p))
那么就有(a^{p-1} equiv 1( Mod) (p))
下面我们用数学归纳法来证明一下费马小定理:
显然(a=1)时结论成立,
若(a=n)时结论成立,
当(a=n+1)时,有
((n+1)^p)(=sum_{i=0}^{p}C_{p}^{i}n^{p-i})(二项式定理)
那么除了(n^{p})和(1)这两项外,
其它的都有一个系数(C_{p}^{i},iin[1,p-1]),所以都能被(p)整除.
而(n^pequiv n(Mod) (p)),
所以((n+1)^{p}equiv n+1(Mod) (p)),结论成立.
所以回到正题,如果对于一个数(p),存在(a in [1,p-1]),(a^{p-1} otequiv1(Mod) (p)),则(p)一定不是质数.
然而仍然有一些数能够逃掉这种检测,
于是就有了
step 2
二次探测!
对于一个质数(p),若(a in[1,p-1]),且(a^2equiv1(Mod) (p)),则(a=1)或(a=p-1).
证明:
若(a^2equiv1(Mod) (p)),则(a^2-1equiv 0(Mod) (p))
即(p|(a+1)(a-1)).
因为(p)为质数,且(ain [1,p-1]),所以只有当(a=1)或(a=p-1)时上式才成立.
所以反过来想当(a)等于其它数时,(p)就不是质数了.
step 3
下面再来讲一下具体的实现.
首先大的思路是费马小定理,在中间用二次探测判定.
对于要判定的数(p),
设(u=p-1),则根据费马小定理有(a^uequiv 1(Mod) (p)),(ain [1,p-1]).
然后把(u)写成(d*2^n)的形式,也就是((((d^2)^2)^2)^{2...})(n个2)
那么从(d^2)开始用二次探测判定,
最后再用费马小定理判定就行了.
而时间复杂度是(O(log_{ iny 2}) (p)).
并且一次的失误率是(frac{1}{4}),
那么(T)次测试的失误率就是(4^{-T}),时间复杂度为(O(Tlog_{ iny 2}p)),能够接受.
并且如果(n<2^{64}),只用选取(a=2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37)测试即可.
code:(似乎要开int128)
#include <iostream>
#include <cstdio>
#include <cstring>
#define int __int128
#define fre(x) freopen(x".in","r",stdin),freopen(x".out","w",stdout)
using namespace std;
inline int read(){
int sum=0,f=1;char ch=getchar();
while(ch>'9'||ch<'0'){if(ch=='-')f=-1;ch=getchar();}
while(ch>='0'&&ch<='9'){sum=sum*10+ch-'0';ch=getchar();}
return f*sum;
}
int p[101]={2,3,5,7,11,13,17,19,23,29,31,37};
int T,n,tot=12;
inline int fpow(int a,int b,int p){
int ret=1;
for(;b;a=a*a%p,b>>=1) if(b&1) ret=ret*a%p;
return ret;
}
inline bool Miller_Rabin(int n){
if(n<=2){
if(n==2) return 1;
return 0;
}
int u=n-1;
while(!(u%2)) u>>=1;
int d=u;
for(int i=0;i<tot&&p[i]<n;i++){
int x=p[i];u=d;x=fpow(x,u,n);
while(u<n){
int y=fpow(x,2,n);
if(y==1&&x!=1&&x!=n-1) return 0;
x=y;u<<=1;
}
if(x!=1) return 0;
}
return 1;
}
signed main(){
T=read();
while(T--){
n=read();
if(Miller_Rabin(n)) puts("Yes");
else puts("No");
}
return 0;
}