Pollard-Rho
算法简介
Pollard-Rho 主要用于质因数分解，是一个随机化算法。
前置算法
伪随机函数
此算法用 \(f_x=f_{x-1}^2+c \bmod n\) 生成随机数用于寻找质因数。
（其中，\(f_1\) 和 \(c\) 可以随机生成， \(n\) 为正在分解的数）
但是，这个函数会产生环，所以叫伪随机函数。
判环
考虑用两个指针，一个每次走一步，另一个每次走两步。
如果两个指针的值在某一时刻相等，则找到了环。
因为两个指针位置的差每次增加1，所以不用担心找不到环。
算法主体
如果当前的数是质数，记录一下，直接返回。
判质数可用 Miller-Rabin 实现。
用上面的伪随机函数生成 \(f\) 数组。
根据生日悖论，随机两个数 \(x\) ， \(y\) ，使得 \(\gcd(|x-y|,n)>1\) 的概率比只随机一个数的概率要大。
所以，用判环使用的两个指针的值作为 \(x\) ， \(y\) ，
设 \(d=\gcd(|x-y|,n)\) 。
若 \(d>1\) 且 \(d<n\)，则递归处理 \(d\) 和 \(\frac{n}{d}\) 两个部分。
同时，如果找到了环，且还没有找到合法的因数，就重新随机 \(f_1\) 和 \(c\) 。
现在的代码
#include<iostream>
#include<cstdlib>
#include<cstdio>
#include<ctime>
#define int long long
#define LD long double
#define ull unsigned long long
using namespace std;
int T,n,ans;
int MUL(int a,int b,int p) //a*b%p
{
int x=(LD)a/p*b;
return ((ull)a*b-(ull)x*p+p)%p;
}
int POW(int a,int b,int p) //a^b%p
{
if(!b) return 1;
if(b==1) return a;
int sum=POW(a,b/2,p);
if(b%2) return MUL(MUL(sum,sum,p),a,p);
return MUL(sum,sum,p);
}
int MAX(int x,int y)
{
return x>y?x:y;
}
int ABS(int x)
{
return x>0?x:-x;
}
int gcd(int x,int y)
{
if(!y) return x;
return gcd(y,x%y);
}
int f(int x,int c,int p)
{
return (MUL(x,x,p)+c)%p;
}
bool MR(int x)
{
if(x==0||x==1) return false;
if(x==2) return true;
if(x%2==0) return false;
int p=x-1,q=0;
while(p%2==0) q++,p/=2;
for(int i=1;i<=10;i++)
{
int a=rand()%(x-1)+1;
if(POW(a,x-1,x)!=1) return false;
int lst=1;
for(int j=0;j<q;j++)
{
int t=POW(a,(1ll<<j)*p,x);
if(t==1&&lst!=1&&lst!=x-1) return false;
lst=t;
}
if(lst!=1&&lst!=x-1) return false;
}
return true;
}
int find(int x)
{
if(x%2==0) return 2;
if(MR(x)) return x;
int t=rand()%(x+1);
int a=t,b=t;
int c=rand()%(x+1);
while(1)
{
a=f(a,c,x),b=f(f(b,c,x),c,x);
int d=gcd(ABS(a-b),x);
if(d>1&&d<x) return d;
if(a==b) return find(x);
}
}
void PR(int x)
{
if(x<=1) return;
if(MR(x))
{
ans=MAX(ans,x);
return;
}
int y=find(x);
PR(y),PR(x/y);
}
signed main()
{
srand(time(0));
scanf("%lld",&T);
while(T--)
{
scanf("%lld",&n);
if(MR(n)) puts("Prime");
else
{
ans=0,PR(n);
printf("%lld\n",ans);
}
}
return 0;
}
黑科技 gcd
上面的代码因为效率过慢过不了洛谷上的 模板题。
观察到代码中会使用很多次 gcd 函数，所以考虑优化计算 gcd 的方法。
新算法
设现在要计算 \(a\) 和 \(b\) 两数的 gcd。
若 \(a\) 和 \(b\) 都是偶数，则递归计算 \(2 \times gcd(\frac{a}{2},\frac{b}{2})\) 。
若 \(a\) 是偶数，且 \(b\) 是奇数，则递归计算 \(gcd(\frac{a}{2},b)\) 。
若 \(a\) 是奇数，\(b\) 是偶数同理。
若 \(a\) 和 \(b\) 都是奇数，则递归计算 \(gcd(|a-b|,a)\) 。
证明：
设原来 \(a\) ， \(b\) 的 gcd 为 \(k\) 。
则 \(a=kA\) ， \(b=kB\) 。
\(|a-b|=k|A-B|\) 。
因为 \(A\) 与 \(B\) 必然互质，
所以，现在要证明的就是 \(A\) 与 \(|A-B|\) 互质。
为了方便，先假设 \(A>B\) ，去掉绝对值。相反的情况的证明也是类似的。
考虑反证法。
设 \(A\) 与 \(A-B\) 不互质，
则 \(A=qx\) ， \(A-B=qy\) 。
则 \(B=A-(A-B)=qx-qy=q(x-y)\) ，所以 \(A\) 与 \(B\) 不互质，不成立。
所以前面的递归是成立的。
把原算法中的 gcd 部分换成这种写法，并加上一些比较常见的卡常技巧即可轻松通过本题。
贴一个用迭代实现的板子。
//__builtin_ctzll(x) 函数用于求x在二进制中末尾0的个数
int gcd(int x,int y)
{
if(!x) return y;
if(!y) return x;
int t=__builtin_ctzll(x|y);
x>>=__builtin_ctzll(x);
do
{
y>>=__builtin_ctzll(y);
if(x>y) swap(x,y);
y-=x;
}while(y);
return x<<t;
}
