因为他我学了龟速乘
Millar-robin 米勒罗宾
这个小东西是用来素数判定的,且听我细细道来。
前置知识
-
肥妈小定理
\[x^{p-1} \equiv 1 \pmod{p} \]又名费马小定理:
当一个数 \(x\) 不是一个质数 \(p\) 的倍数时有:证明:
对于一个序列
\[b = \left \{1,2,3....p-1\right \} \]令
\[a = \left \{ x,2x,3x...(p-1)x \right \} \]那么 \(a\) 序列对 \(p\) 取模后恰好为 \(p-1\) 的一个排列。
证明:
假设存在 \(s,t\) 使得 \(sx \equiv tx \pmod p\)
那么 \((s-t)x \equiv 0 \pmod p\)
又有 \(x \ne 0 \pmod p \wedge s-t < p\)
所以 \(s-t=0\)
所以不存在两个位置模数相同,又有互质,所以模数恰好是一个 \(p-1\) 的排列。(分割线)
所以有
\[(p-1)! \equiv (p-1)! \times x^{p-1} \pmod p \]那么根据
\[\gcd((p-1)!,p) = 1 \]所以有
\[x^{p-1} \equiv 1 \pmod p \] -
二次探测定理
\[x^2 \equiv 1 \pmod p \]
对于质数 \(p\) 和一个 \(x \in \left [1,p-1 \right ]\)则有
\[x_{1} = 1,x_{2} = {p-1} \]证明
有
\[x^2 - 1 \equiv 0 \pmod p \]\[(x+1)(x-1) \equiv 0 \pmod p \]故当
\[x+1 \equiv 0 \pmod p \]有
\[x = p - 1 \]当
\[x - 1 \equiv 0 \pmod p \]有
\[x = 1 \]
进入正题
如何检查 \(n\) 是否是质数?
首先忽略 \(2\) 的情况,那么 \(n\) 必定为奇,那么 \(n-1\) 必定为偶。
令 \(n-1 = 2^k \times m\)
米勒罗宾是通过费马小定理来测试 \(n\) 是否为素数的。
由费马小定理有这样的式子,若 $n $ 为质数那么 \(x^{n-1} \equiv 1 \pmod n\),但是这只是 \(n\) 是素数的必要条件,所以这是一个 \(\mathbf{测试}\) 算法。
考虑我们怎么能得到这样的式子,接下来默认 \(y \in [1,n-1]\) ,因为用到了二次探测。
若 \(y^{n-1} \equiv y^{2^km} \equiv1 \pmod n\)
那么 \(y^{2^{k-1}m} = 1\) 或 \(n-1\)
若 \(y^{2^{k-1}m} = 1\) 则 \(y^{2^{k-2}m} = 1\) 或 \(n-1\)
如此反复直到最后 \(y^m = 1\) 或 \(y^m = n - 1\)。
若 \(n\) 是质数,那么 \(y^m = 1\) 或者 存在 \(g\in[0,k-1]\) 使得 \(y^{2^gm} = n - 1\) 符合该条件则通过测试
这个测试,合数不一定通不过,素数肯定通过,所以通过的不一定是素数,通不过的一定是合数 我在说什么。
我们选用素数来测试 \(n\),算法竞赛中选 \(30\) 内的就足以让 unsigned long long 范围内的数合法。选的数越多,这个算法越准。至于更深的探究,我也不会了。
例题 sp288。
code :
#include<bits/stdc++.h>
#define rep(i,a,b) for(int i=(a);i<=(b);++i)
#define Rep(i,a,b) for(int i=(a);i<(b);++i)
#define rrep(i,a,b) for(int i=(a);i>=(b);--i)
using namespace std;
template <typename T>
inline void read(T &x){
x=0;char ch=getchar();bool f=0;
while(ch<'0'||ch>'9'){if(ch=='-')f=1;ch=getchar();}
while(ch>='0'&&ch<='9')x=(x<<1)+(x<<3)+(ch^48),ch=getchar();
if(f)x=-x;
}
template <typename T,typename ...Args>
inline void read(T &tmp,Args &...tmps){read(tmp);read(tmps...);}
typedef long long ll;
ll pr[10] = {2,3,5,7,11,13,17,19,23,29};
inline void add(ll &x,ll y,ll z){
x += y;
if(x >= z)x -= z;
}
ll gui_mul(ll x,ll y,ll z){
ll ans = 0;
while(y){
if(y & 1)ans = (ans + x) % z;
x = (x + x) % z;
y >>= 1;
}
return ans;
}
ll ksm(ll x,ll y,ll z){
ll res = 1;
while(y){
if(y & 1)res = gui_mul(res,x,z);
x = gui_mul(x,x,z);
y >>= 1;
}
return res;
}
bool millar_rabin(ll p){
if(p < 2)return 0;
if(p != 2 && p % 2 == 0)return 0;
ll s = p - 1;
while(!((s) & 1))s >>= 1;
rep(i,0,9){
if(p == pr[i])return 1;
ll t = s,m = ksm(pr[i],s,p);
while(t != (p - 1) && m != 1 && m != p - 1){
m = gui_mul(m,m,p);
t <<= 1;
}
if(m != p - 1 && !(t & 1))return 0;
}
return 1;
}
signed main(){
int t;
read(t);
while(t--){
ll x;
read(x);
puts(millar_rabin(x) ? "YES" : "NO");
}
}
标签:ch,Millar,pmod,定理,while,罗宾,ll,equiv
From: https://www.cnblogs.com/wsxxs/p/16646162.html