Miller_rabin 素数测试 学习笔记
Miller_rabin 素数测试
一种用来判断素数的算法。
前置芝士
威尔逊定理
若 \(p\) 为素数，\((p-1)! \equiv -1 (\mod p)\)。
证明：
充分性证明：
如果 \(p\) 不是素数，那么他的因数必定存在于$ 1,2,3,\dots,p−1$ 之中，所以 \(\gcd((p-1)!,p)\)，那么 \((p-1)! \not\equiv -1\)。
必要性证明：
首先，我们知道 $$p-1\equiv -1 (\mod p)$$
那么我们只需要证明 \((p-2)! \equiv 1 (\mod p)\) 就可以了。
设 \(A=\{2,3\dots,p-2\}\)
对于 \(x\in A\)，肯定存在一个 \(a \in A\)，使得 \(ij\equiv 1(\mod p)\)
当 \(x=a\) 时，考虑二次剩余 \(x^2 \equiv 1(\mod p)\)。
移项就可以得到 \((x+1)(x-1) \equiv 0 (\mod p)\)，
那么只有两个解，\(x \equiv 1(\mod p)\),\(x \equiv p-1(\mod p)\)，不成立。
所以其他情况都可以找到对应的 \(a\)。
所以 \((p-2)!\equiv 1(\mod p)\)，\((p-1)!\equiv p-1(\mod p)\)。
费马小定理
若 \(p\in\mathbb{P},\gcd(a,p)=1\)，则 \(a^{p-1}\equiv 1(\mod p)\)，
证明：
因为 \(1,2,3,\dots ,p-1\) 是 \(p\) 的完全剩余系，\(a,p\) 互质，
所以 \(a,2* a,3* a,4* a, \dots (p-1)* a\) 也是 \(\mod p\) 的完全剩余系。
所以 \(1 * 2 * 3 * \dots * (p-1) \equiv a * 2a * 3a * \dots * (p-1)a (\mod p)\)。
就是 \((p-1)! \equiv (p-1)!a^{p-1} (\mod p)\)。
由威尔逊定理得出，\((p-1)! \equiv -1(\mod p)\)，
两边同时约去 \((p-1)!\)，
所以可以得出 \(a^{p-1} \equiv 1(\mod p)\)。
二次探测定理
若 \(p\) 为素数，\(x^2 \equiv 1(\mod p)\)，那么\(x \equiv \pm 1(\mod p)\)。
证明：
原式移项就可以得到 \((x+1)(x-1) \equiv 0 (\mod p)\)，
那么只有两个解，\(x \equiv 1(\mod p)\),\(x \equiv p-1(\mod p)\)。
但是这些又和 Miller_rabin 有什么关系呢？
我们把 \(p-1\) 分解为 \(2^k* t\)，当 \(p\) 是素数时，则有 \(a^{2^k* t} \equiv 1(\mod p)\)。
然后随机出一个数 \(a\)，可以算出 \(a^t\)，然后再自乘，就可以得到 \(a^{2^k* t}\) ，也就是 \(a^{p-1}\)。
但如果我们在自乘之后 \(\equiv 1(\mod p)\)，而之前 \(\not\equiv 1(\mod p)\)，那么就违背了二次探测定理，则 \(p\) 不是素数。
else
\(Zwad\) 告诉我，若 \(p\) 通过一次测试，则p不是素数的概率仅为25%。
那么用 \(2,325,9375,28178,450775,9780504,1795265022\) 这几个数进行判断,在 $long long $ 范围内能保证正确。
Code
#include<bits/stdc++.h>
#define int __int128
#define N_4 10004
#define N_5 100005
#define N_6 1000006
#define Mod 1000000007
#define FOR(i,j,k) for(long long i=j;i<=k;++i)
#define ROF(i,j,k) for(long long i=j;i>=k;--i)
using namespace std;
inline int read(){
int x=0,f=1;
char ch;
ch=getchar();
while(ch<'0'||ch>'9'){
if(ch=='-') f=-f;
ch=getchar();
}
while(ch>='0'&&ch<='9'){
x=x*10+(ch-'0');
ch=getchar();
}
return x*f;
}
int T,n,tot;
int gg[6]={0,2,7,61,63,97};
bool isprime[500005];
int prime[500005],an[500005];
inline void Euler(int n){
isprime[1]=true;isprime[0]=true;
for(register int i=2;i<=n;i++){
if(!isprime[i])
prime[++tot]=i;
an[i]=tot;
for(register int j=1;j<=tot&&prime[j]*i<=n;j++){
isprime[i*prime[j]]=true;
if(i%prime[j]==0)
break;
}
}
return;
}
inline int ksm(int a,int b,int mod){
int ans=1,p=a,g=b;
while(g){
if(g&1) ans=(ans*p)%mod;
p=(p*p)%mod;
g>>=1;
}
return ans;
}
int mul(int a,int b,int p){return (a*b-(int)((__float128)a/p*b)*p+p)%p;}
inline bool miller_rabin(int P){
if(P==1) return 0;
int t=P-1,k=0;
while(!(t&1)) ++k,t>>=1;
for(register int i=1;i<=5;++i){
if(P==gg[i]) return true;
int a=ksm(gg[i],t,P),nxt=a;
for(register int j=1;j<=k;++j){
nxt=mul(nxt,nxt,P);
if(nxt==1&&a!=1&&a!=P-1) return false;
a=nxt;
}
if(a!=1) return false;
}
return true;
}
signed main(){
T=read();
Euler(500000);
while(T--){
n=read();
if(n<500000){
if(isprime[n]) cout<<"NO"<<endl;
else cout<<"YES"<<endl;
}
else{
if(!miller_rabin(n)) cout<<"NO"<<endl;
else cout<<"YES"<<endl;
}
}
return 0;
}