[模板] Miller_Rabin和Pollard_Rho
Miller_Rabin
快速($O(slogn)$，s为尝试次数)地测试一个数是否是质数
首先有费马小定理$a^{p-1}=1\ (mod\ p)$当p为质数时成立，所以可以随机选择a来以这个式子作为一定的判断依据，但并不是所有合数都不满足这个式子，甚至存在合数对所有的a都不满足这个式子
然后有二次探测定理：若p是质数，$a^2=1\ (mod\ p)$成立当且仅当$a=1\ (mod\ p)$或$a=p-1\ (mod\ p)$
证明：移项可得$(a-1)(a+1)=0\ (mod\ p)$，于是有$p|(a-1)(a+1)$，若p为质数，则只能有$p|(a-1)$或者$p(a+1)$
所以，如果有$a^2=1\ (mod\ p)$但a不等于1或者p-1的话，证明p一定是合数
然而当我们随便选一个a，很难保证它的平方就等于1...但我们还有费马小定理！
所以我们可以考虑把p-1表示成$s*2^t$的形式，然后从$a^s$一路平方上去，并进行判断
这个a可以考虑选择前10个质数，在1e18的范围内应该都是对的
Pollard_Rho
快速(期望$O(n^{1/4})$)地找到n的某个非平凡因子
首先，用Miller_Rabin来判断n是否是质数（这也是P_R和M_R唯一的关联）
然后我们考虑随机一个数x，然后计算$gcd(x,n)$，如果结果不是1或者n，就找到了n的一个非平凡因子
然而纯随机肯定是不行的..
考虑一种生成随机数的方式:$x_i=x_{i-1}^2+c\ (mod\ n)$，$c$和$x_0$为随机的常数，这样生成出的数一定会成一个"ρ"形，而且根据生日悖论(期望随机$\sqrt{n}$次就能得到在n的范围内相等的两个数)，环长和前面一段的长度都是期望$\sqrt{n}$的
更进一步的是，如果考虑a是n的一个最小质因数，那么$x_i$在模a意义下，会形成一个期望$\sqrt{a}$的ρ
这样的话，我随机在这个模n的ρ形上取两个数做gcd，期望只要$\sqrt{a}$也就是$n^{1/4}$次，就能找到一个n的因数
于是可以使用a走一步，b走两步的方法判圈，同时计算$gcd(a-b,n)$
但注意到一次gcd需要log的时间，所以我们可以采用每127（我也不知道这个玄学的数是哪来的）个计算一次gcd的方法来优化，就是说，每次连续算127个$a-b$并乘起来，然后算一下gcd，如果发现gcd大于1了，再重新算一下这127个看看是哪个大于1了
例题
luogu4718 (求n的最大质因子)
1 #include<bits/stdc++.h> 2 #include<tr1/unordered_map> 3 #define CLR(a,x) memset(a,x,sizeof(a)) 4 #define MP make_pair 5 #define fi first 6 #define se second 7 using namespace std; 8 typedef long long ll; 9 typedef unsigned long long ull; 10 typedef long double ld; 11 typedef pair<int,int> pa; 12 const int maxn=233; 13 14 inline ll rd(){ 15 ll x=0;char c=getchar();int neg=1; 16 while(c<'0'||c>'9'){if(c=='-') neg=-1;c=getchar();} 17 while(c>='0'&&c<='9') x=x*10+c-'0',c=getchar(); 18 return x*neg; 19 } 20 21 int pri[10]={2,3,5,7,11,13,17,19,23,29}; 22 23 inline ll mul(ll a,ll b,ll p){ 24 return ((a*b-(ll)((ld)a*b/p)*p)%p+p)%p; 25 } 26 27 inline ll fpow(ll a,ll b,ll p){ 28 ll r=1; 29 while(b){ 30 if(b&1) r=mul(r,a,p); 31 a=mul(a,a,p),b>>=1; 32 }return r; 33 } 34 35 inline bool judge(ll p){ 36 if(p==2) return 1; 37 else if(p==1||p%2==0) return 0; 38 ll s=0,t=p-1; 39 while(!(t&1)) s++,t>>=1; 40 for(int i=0;i<10&&pri[i]<p;i++){ 41 ll x=fpow(pri[i],t,p); 42 for(int j=1;j<=s;j++){ 43 ll r=mul(x,x,p); 44 if(r==1&&x!=1&&x!=p-1) return 0; 45 x=r; 46 } 47 if(x!=1) return 0; 48 }return 1; 49 } 50 51 inline ll gcd(ll a,ll b){return a?gcd(b%a,a):b;} 52 inline ll rand(ll r){ 53 return (((ll)rand()<<30)|rand())%(r-1)+1; 54 } 55 56 inline ll calc(ll x){ 57 if(judge(x)) return x; 58 while(1){ 59 ll a=rand(x-1),c=rand(x-1),b=(mul(a,a,x)+c)%x; 60 // printf("~%lld %lld %lld\n",a,b,c); 61 bool bl=0; 62 while(a!=b){ 63 ll s=1;ll aa=a,bb=b; 64 if(!bl){ 65 for(int i=1;i<=127&&aa!=bb;i++){ 66 s=mul(s,aa-bb,x); 67 aa=(mul(aa,aa,x)+c)%x; 68 bb=(mul(bb,bb,x)+c)%x,bb=(mul(bb,bb,x)+c)%x; 69 } 70 } 71 ll g=gcd(s+x,x); 72 if(g>1) bl=1; 73 else a=aa,b=bb; 74 if(bl){ 75 g=gcd(a-b+x,x); 76 if(g>1) return max(calc(g),calc(x/g)); 77 a=(mul(a,a,x)+c)%x; 78 b=(mul(b,b,x)+c)%x,b=(mul(b,b,x)+c)%x; 79 } 80 81 82 } 83 } 84 } 85 86 int main(){ 87 srand(1919810); 88 for(int T=rd();T;T--){ 89 ll x=rd(); 90 ll re=calc(x); 91 if(re==x) puts("Prime"); 92 else printf("%lld\n",re); 93 } 94 return 0; 95 }
