『学习笔记』Pollard-Rho 算法 题解
前言
说是学习笔记，其实窝并没有打算写太多（太麻烦了，而且我的理解还不是特别深，可能也写不清楚），所以打算大概写两句，然后贴个板子。
前置知识
- \(Miller-Rabin\) 素性测试。
- 倍增基础应用。
\(Miller-Rabin\) 素性测试
我们知道有费马小定理： \(a^{p - 1} \equiv 1\pmod p\)，\(p\) 是质数且 \(a\) 是小于 \(p\) 的正整数。
那么费马小定理的逆定理是否成立呢？
答案是否定的，\(2^{340} \equiv 1 \pmod{341}\)，但是 \(341 = 11 \times 13\)。
这样的数被称为伪素数。事实上，这样的数有很多。
有人会问：如果我选取小于 \(n\) 的所有质数当作底数都判断一遍，是否能让这个逆定理变得正确呢？
很遗憾，就算取再多的底数总是会有合数通过条件，例如 561，所以我们并不能用费马小定理的逆定理当作判断素数的方法。
这时 \(Miller\) 和 \(Rabbin\) 出现了，他们共同提出了二次探测定理。
这个定理是这样的：
假设有 \(x^2 \equiv 1 \pmod p\)，那么 \(x = 1\) 或 \(x = p - 1\)。
证明很简单，原式即为：
所以 \(x = 1\) 或 \(x = p - 1\)。
注意：\(x\) 等于上述两个数中的任意一个都表示它有可能是质数。
下面我们来演示一下上面的定理如何应用在费马素性测试上：
\(2^{340} = (2^{170})^2\)，把上面的定理套进去，并计算 \(2^{170}\ \ mod\ \ 341\) 的值，我们发现是 1。于是继续套用公式，计算 \(2^{85}\ \ mod\ \ 341\)，发现其等于 32。这时，341 的真正面目就暴露了出来，它是一个合数！
所以我们的计算过程就是不断地令指数除以 2，然后对模数取模进行判断。
当出现取模之后的结果不为 1 或 \(p - 1\) 时，证明它不是质数。如果指数已经是奇数了，那么它就是一个质数。
需要一提的是，单独使用 2 这个底数进行判断还是会有误差，例如 2047 就能通过以 2 为底数的 \(MR\) 素性测试。
上述就是 \(Miller-Rabin\) 素性测试的原理了。
\(Code:\)
namespace Miller_Rabin{
const int p[9] = {0, 2, 3, 5, 7, 11, 13, 17, 19}, totp = 8;
inline ll mul(ll a, ll b, ll mod) {return (a * b - (ll)((long double)a * b / mod) * mod + mod) % mod;}
inline ll qpow(ll a, ll b, ll mod){
ll res = 1; a %= mod;
while(b){
if(b & 1) res = mul(res, a, mod);
a = mul(a, a, mod), b >>= 1;
}
return res;
}
inline bool check(ll x, ll p){
if((x % p == 0) || qpow(p, x - 1, x) != 1) return 0;
ll k = x - 1;
while(!(k & 1)){
ll t = qpow(p, k >>= 1, x);
if(t != 1 && t != (x - 1)) return 0;
if(t == x - 1) return 1;
}
return 1;
}
inline bool MR(ll x){
if(x <= 1) return 0;
if(x <= 3) return 1;
if(!(x % 2) || !(x % 3)) return 0;
for(int i = 1; i <= totp; ++i){
if(x == p[i]) return 1;
if(!check(x, p[i])) return 0;
}
return 1;
}
}
using namespace Miller_Rabin;
倍增就不多说了。
\(Pollard-Rho\)
核心思想：每次猜一个数，判断是否是 \(n\) 的约数，并递归猜下去，如果是约数且是质数的话就更新答案。
考虑构造一个随机数列 \(x_i = x_{i - 1}^2 + c\)，通过列举不难发现，这个东西是一定会形成环的。
我们可以在上面倍增找点。当找到当前的环上的点，判断其是否是 \(n\) 的约数。
（不想写太多了，各位巨佬不喜勿喷）
代码
#include <bits/stdc++.h>
#define ll long long
using namespace std;
namespace IO{
inline ll read(){
ll x = 0;
char ch = getchar();
while(!isdigit(ch)) ch = getchar();
while(isdigit(ch)) x = (x << 3) + (x << 1) + ch - '0', ch = getchar();
return x;
}
template <typename T> inline void write(T x){
if(x > 9) write(x / 10);
putchar(x % 10 + '0');
}
}
using namespace IO;
namespace Miller_Rabin{
const int p[9] = {0, 2, 3, 5, 7, 11, 13, 17, 19}, totp = 8;
inline ll mul(ll a, ll b, ll mod) {return (a * b - (ll)((long double)a * b / mod) * mod + mod) % mod;}
inline ll qpow(ll a, ll b, ll mod){
ll res = 1; a %= mod;
while(b){
if(b & 1) res = mul(res, a, mod);
a = mul(a, a, mod), b >>= 1;
}
return res;
}
inline bool check(ll x, ll p){
if((x % p == 0) || qpow(p, x - 1, x) != 1) return 0;
ll k = x - 1;
while(!(k & 1)){
ll t = qpow(p, k >>= 1, x);
if(t != 1 && t != (x - 1)) return 0;
if(t == x - 1) return 1;
}
return 1;
}
inline bool MR(ll x){
if(x <= 1) return 0;
if(x <= 3) return 1;
if(!(x % 2) || !(x % 3)) return 0;
for(int i = 1; i <= totp; ++i){
if(x == p[i]) return 1;
if(!check(x, p[i])) return 0;
}
return 1;
}
}
using namespace Miller_Rabin;
namespace Pollard_Rho{
#define Rand(x) (1ll * rand() * rand() % (x) + 1)
inline ll add(ll x, ll mod) {return x >= mod ? x - mod : x;}
inline ll gcd(ll a, ll b) {return !b ? a : gcd(b, a % b);}
inline ll PR(ll n){
ll x0 = Rand(n - 1), x = x0, c = Rand(n), s = 1;
ll k = 1, step = 0, d;
while(true){
x = add(mul(x, x, n) + c, n);
s = mul(s, abs(x - x0), n);
if(!s || (x == x0)) return n;
if((++step) == k){
if((d = gcd(s, n)) != 1) return d;
x0 = x, k <<= 1;
}
}
}
inline void solve(ll x, ll &ans){
if((x == 1) || x < ans) return;
if(MR(x)) return ans = max(ans, x), void();
ll y = PR(x);
while(x == y) y = PR(x);
while(x % y == 0) x /= y;
solve(x, ans), solve(y, ans);
}
inline ll calc(ll n){
ll ans = 0;
solve(n, ans);
return ans;
}
}
using namespace Pollard_Rho;
ll T, n;
signed main(){
srand(20060707);
T = read();
while(T--){
ll n = read();
if(MR(n)) puts("Prime");
else write(calc(n)), puts("");
}
return 0;
}
