Miller-Rabin 与 Pollard-Rho 学习笔记
跟着大佬脚步学习:D
Miller-Rabin
题目
给你一个数字 \((z\le10^{18})\)，判断它是不是质数。
公式
费马小定理：（ \(p\) 为质数）
\[a^{p-1}\equiv1\pmod p
\]
通俗地讲是 \(a^{p-1}\bmod p=1\)，但是逆定理不一定成立（大概率成立），即若 \(a^{p-1}\equiv1\pmod p\) 成立， \(p\) 也不一定是质数，不过如果不成立， \(p\) 一定不是质数。
二次探测定理：若 \(p\) 为质数（ \(2\) 除外）且 \(a^2\equiv1\pmod p\)，则 \(a\equiv 1\) 或 \((p-1)\pmod p\) 。
算法
特判 \(z\) 是不是偶数，直接再见。
现在 \(z\) 是奇数了， \((z-1)\) 为偶数，所以令 \(z-1=2^q+p\)，然后多次进行下面算法流程：
- 随机选择一个 \(a\in[1,z-1]\)，若 \(a^{z-1}\not\equiv1\pmod p\) ，则它是合数，再见。（费马小定理）
- 否则枚举所有的 \(b\in [0,q-1]\)，如果 \(\Large a^{2^b}\bmod z\not=1\) 但是 \(\Large a^{2^{b-1}\times p}\bmod z\not=1\) 或 \(\Large(z-1)\) 那么它是合数，再见。 （二次探测定理）
重复若干遍后如果还是没有再见，那么它就是素数了。
时间复杂度是 \(O(k\log^2n)\) 的， \(k\) 为重复次数。但是如果数据小，效率会低于 \(O(\sqrt n)\) 的暴力，所以酌情选择。
代码
注意要用龟速乘，怕溢出。（用了 LH 大奆的 \(O(1)\) '鬼'速乘）
Drand() 是小常数生成随机数（其实 rand() 更快，但是一般我们都是用 rand()*rand()+rand() 之类的以提高范围和随机性，这就会变很慢）。
#include<bits/stdc++.h>
#define rep(i,x,y) for(int i=x,WALL=y;i<=WALL;++i)
#define lon long long
using namespace std;
lon rad=1;
int rd(){
int shu=0;char ch=getchar();
while( !isdigit(ch) )ch=getchar();
while( isdigit(ch) )shu=(shu<<1)+(shu<<3)+(ch^48),ch=getchar();
return shu;
}
lon Drand(){
unsigned lon x=rad;
x^=x<<13,x^=x>>17,x^=x<<5;
rad=x/4;return rad;
}
lon Dmul(lon p,lon q,lon mo) {
lon shul=p*(q>>32ll)%mo*(1ll<<25)%mo;
lon shur=p*(q & ( (1ll<<32)-1 ) )%mo;
return (shul+shur)%mo;
}
lon Dpow(lon p,lon q,lon mo){
lon tot=1;
while(q){
if(q&1)tot=Dmul(tot,p,mo)%mo;
p=Dmul(p,p,mo)%mo,q=q>>1;
}
return tot;
}
bool rabin(lon z){
if(z<=1)return 0;
if(z==2)return 1;
if(!(z&1))return 0;
lon q=0,p=z-1;
while(p&1)p=p/2,q++;
int chktime=10;
while(chktime--){
lon x=Drand()%(z-1)+1;
if(Dpow(x,z-1,z)^1)return 0;
lon las=1;
rep(i,0,q-1){
lon tmp=Dpow(x,(1ll<<i)*p,z);
if(tmp==1&&las^1&&las^(z-1))return 0;
las=tmp;
}
if(las^1&&las^(z-1))return 0;
}
return 1;
}
int main(){
int T=rd();
while(T--){
lon n=rd();
puts(rabin(n)?"Prime":"notPrime");
}
return 0;
}
Pollard-Rho
贴个代码，待填
Luogu模板题，求出一个数的最大质因数
#include<bits/stdc++.h>
#define rep(i,x,y) for(int i=x,WALL=y;i<=WALL;++i)
#define lon long long
#define abs(x) ((x)<0?-(x):(x))
using namespace std;
lon rad=1,ans;
lon rd(){
lon shu=0;char ch=getchar();
while( !isdigit(ch) )ch=getchar();
while( isdigit(ch) )shu=(shu<<1)+(shu<<3)+(ch^48),ch=getchar();
return shu;
}
lon Drand(){
unsigned lon x=rad;
x^=x<<13,x^=x>>17,x^=x<<5;
rad=x/4;return rad;
}
lon Dmul(lon p,lon q,lon mo){
p=p%mo,q=q%mo;
return (p*q-(lon)((long double)p/mo*q)*mo+mo)%mo;
}
lon Dpow(lon p,lon q,lon mo){
lon tot=1;
while(q){
if(q&1)tot=Dmul(tot,p,mo)%mo;
p=Dmul(p,p,mo)%mo,q=q>>1;
}
return tot;
}
bool mrabin(lon z){
if(z<=1)return 0;
if(z==2)return 1;
if(!(z&1))return 0;
lon q=0,p=z-1;//z=2^q+p
while(!(p&1))p=p/2,q++;
int chktime=10;
while(chktime--){
lon x=Drand()%(z-1)+1;
if(Dpow(x,z-1,z)^1)return 0;
lon las=1;
rep(i,0,q-1){
lon tmp=Dpow(x,(1ll<<i)*p,z);
if(tmp==1&&las^1&&las^(z-1))return 0;
las=tmp;
}
if(las^1&&las^(z-1))return 0;
}
return 1;
}
lon Dgcd(lon p,lon q){
lon tmp;
while(q)tmp=p,p=q,q=tmp%q;
return p;
}
lon prho(lon z){
int len=1;
lon now=0,tot=1,lead,c=Drand()%(z-1)+1;
while(1){
lead=now;
rep(i,1,len){
now=(Dmul(now,now,z)+c)%z;
tot=Dmul(tot,abs(lead-now),z);
if( !(i%127) ){
lon gcd=Dgcd(tot,z);
if(gcd>1)return gcd;
}
}
lon gcd=Dgcd(tot,z);
if(gcd>1)return gcd;
len=len<<1;
}
}
void work(lon z){
if(z<=ans||z<2)return;
if( mrabin(z) ){ans=max(ans,z);return;}
lon x=z;
while(x==z)x=prho(z);
while( !(z%x) )z=z/x;
work(x),work(z);
}
int main(){
int T=rd();
while(T--){
lon n=rd();
if( mrabin(n) )puts("Prime");
else ans=0,work(n),printf("%lld\n",ans);
}
return 0;
}
