## Miller_rabin 素數測試 一種用來判斷素數的演算法。 ### 前置芝士 #### 威爾遜定理 若 $p$ 為素數，$(p-1)! \equiv -1 (\mod p)$。 證明： 充分性證明： 如果 $p$ 不是素數，那麼他的因數必定存在於$ 1,2,3,\dots,p−1$ 之中，所 ...
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;
}
