0%

板子 Pollard-Rho & Miller-Rabin

对于大数字的玄学复杂度质因数分解方法 Pollard-Rho,以及快速的素数判定方法 Millar-Rabin。

P4718 【模板】Pollard-Rho算法

前置知识

快速乘

1
2
3
inline LL mul(LL a, LL b, LL p) {
return (a * b - (LL)((long double)a * b / p) * p + p) % p;
}

快速幂

GCD

费马小定理

Miller-Rabin

很好的讲解及代码:OI Wiki

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
bool Miller_Rabin(LL n) {
if (n < 3) return n == 2;
LL a = n - 1, b = 0;
while (!(a & 1)) a >>= 1, b++;
for (int i = 1, j; i <= 8; ++i) {
LL u = rand() % (n - 2) + 2, v = power(u, a, n);
if (v == 1 || v == n - 1) continue;
for (j = 1; j <= b; ++j) {
v = mul(v, v, n);
if (v == n - 1) break;
}
if (j == b + 1) return false;
}
return true;
}

Pollard-Rho

很好的讲解及代码:OI WikiLinearODE's blog

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
inline LL Pollard_Rho(LL n) {
LL s = 0, t = 0, val = 1, c = rand() % (n - 1) + 1;
for (int goal = 1; ; goal <<= 1) {
s = t, val = 1;
for (int step = 1; step <= goal; ++step) {
t = next_rand(t, c, n);
val = mul(val, t - s, n);
if (step % 127 == 0) {
LL d = gcd(val, n);
if (d > 1) return d;
}
}
LL d = gcd(val, n);
if (d > 1) return d;
}
}

模板题代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
#include <ctime>
#include <cstdio>
#include <cctype>
#include <cstdlib>
#include <algorithm>

using namespace std;

namespace RenaMoe {

template <typename TT> inline void read(TT &x) {
x = 0; bool f = false; char ch = getchar();
while (!isdigit(ch)) f |= ch == '-', ch = getchar();
while (isdigit(ch)) x = x * 10 + ch - '0', ch = getchar();
x = f ? -x : x;
}

typedef long long LL;
typedef long double LD;

LL T, n, maxans;

inline LL mul(LL a, LL b, LL p) {
return (a * b - (LL)((LD)a * b / p) * p + p) % p;
}

inline LL power(LL a, LL b, LL p) {
LL res = 1;
while (b) {
if (b & 1) res = mul(res, a, p);
a = mul(a, a, p), b >>= 1;
}
return res;
}

inline LL next_rand(LL x, LL c, LL p) {
return (mul(x, x, p) + c) % p;
}

LL gcd(LL a, LL b) {
return b ? gcd(b, a % b) : a;
}

bool Miller_Rabin(LL n) {
if (n < 3) return n == 2;
LL a = n - 1, b = 0;
while (!(a & 1)) a >>= 1, b++;
for (int i = 1, j; i <= 8; ++i) {
LL u = rand() % (n - 2) + 2, v = power(u, a, n);
if (v == 1 || v == n - 1) continue;
for (j = 1; j <= b; ++j) {
v = mul(v, v, n);
if (v == n - 1) break;
}
if (j == b + 1) return false;
}
return true;
}

inline LL Pollard_Rho(LL n) {
LL s = 0, t = 0, val = 1, c = rand() % (n - 1) + 1;
for (int goal = 1; ; goal <<= 1) {
s = t, val = 1;
for (int step = 1; step <= goal; ++step) {
t = next_rand(t, c, n);
val = mul(val, t - s, n);
if (step % 127 == 0) {
LL d = gcd(val, n);
if (d > 1) return d;
}
}
LL d = gcd(val, n);
if (d > 1) return d;
}
}

inline void fac(LL n) {
if (n <= maxans || n < 2) return;
if (Miller_Rabin(n)) {
maxans = max(maxans, n);
return;
}
LL p = n;
while (p >= n) p = Pollard_Rho(n);
while (!(n % p)) n /= p;
fac(n), fac(p);
}

inline void main() {
srand(time(NULL));
read(T);
while (T--) {
read(n);
maxans = 0;
fac(n);
if (maxans == n) puts("Prime");
else printf("%lld\n", maxans);
}
}
}

int main() {
RenaMoe::main();
return 0;
}