Skip to content

Commit b79d09a

Browse files
adamant-pwnadamant-pwn
committed
Factor out primality + factorize
1 parent 080445b commit b79d09a

File tree

5 files changed

+83
-73
lines changed

5 files changed

+83
-73
lines changed

‎cp-algo/math/factorize.hpp

Lines changed: 47 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,47 @@
1+
#ifndef CP_ALGO_MATH_FACTORIZE_HPP
2+
#define CP_ALGO_MATH_FACTORIZE_HPP
3+
#include "primality.hpp"
4+
#include "../random/rng.hpp"
5+
namespace cp_algo::math {
6+
// https://en.wikipedia.org/wiki/Pollard%27s_rho_algorithm
7+
void factorize(uint64_t m, std::vector<int64_t> &res) {
8+
if(m % 2 == 0) {
9+
factorize(m / 2, res);
10+
res.push_back(2);
11+
} else if(is_prime(m)) {
12+
res.push_back(m);
13+
} else if(m > 1) {
14+
using base = dynamic_modint;
15+
base::with_mod(m, [&]() {
16+
base t = random::rng();
17+
auto f = [&](auto x) {
18+
return x * x + t;
19+
};
20+
base x, y;
21+
base g = 1;
22+
while(g == 1) {
23+
for(int i = 0; i < 64; i++) {
24+
x = f(x);
25+
y = f(f(y));
26+
if(x == y) [[unlikely]] {
27+
t = random::rng();
28+
x = y = 0;
29+
} else {
30+
base t = g * (x - y);
31+
g = t == 0 ? g : t;
32+
}
33+
}
34+
g = std::gcd(g.getr(), m);
35+
}
36+
factorize(g.getr(), res);
37+
factorize(m / g.getr(), res);
38+
});
39+
}
40+
}
41+
std::vector<int64_t> factorize(int64_t m) {
42+
std::vector<int64_t> res;
43+
factorize(m, res);
44+
return res;
45+
}
46+
}
47+
#endif // CP_ALGO_MATH_FACTORIZE_HPP

‎cp-algo/math/number_theory.hpp

Lines changed: 1 addition & 69 deletions
Original file line numberDiff line numberDiff line change
@@ -2,14 +2,12 @@
22
#define CP_ALGO_MATH_NUMBER_THEORY_HPP
33
#include "../random/rng.hpp"
44
#include "affine.hpp"
5-
#include "modint.hpp"
5+
#include "factorize.hpp"
66
#include <algorithm>
77
#include <optional>
88
#include <vector>
99
#include <bit>
1010
namespace cp_algo::math {
11-
std::vector<int64_t> factorize(int64_t m);
12-
1311
int64_t euler_phi(int64_t m) {
1412
auto primes = factorize(m);
1513
std::ranges::sort(primes);
@@ -91,72 +89,6 @@ namespace cp_algo::math {
9189
}
9290
}
9391
}
94-
// https://en.wikipedia.org/wiki/Miller–Rabin_primality_test
95-
bool is_prime(uint64_t m) {
96-
if(m == 1 || m % 2 == 0) {
97-
return m == 2;
98-
}
99-
// m - 1 = 2^s * d
100-
int s = std::countr_zero(m - 1);
101-
auto d = (m - 1) >> s;
102-
using base = dynamic_modint;
103-
auto test = [&](base x) {
104-
x = bpow(x, d);
105-
if(std::abs(x.rem()) <= 1) {
106-
return true;
107-
}
108-
for(int i = 1; i < s && x != -1; i++) {
109-
x *= x;
110-
}
111-
return x == -1;
112-
};
113-
return base::with_mod(m, [&](){
114-
// Works for all m < 2^64: https://miller-rabin.appspot.com
115-
return std::ranges::all_of(std::array{
116-
2, 325, 9375, 28178, 450775, 9780504, 1795265022
117-
}, test);
118-
});
119-
}
120-
// https://en.wikipedia.org/wiki/Pollard%27s_rho_algorithm
121-
void factorize(uint64_t m, std::vector<int64_t> &res) {
122-
if(m % 2 == 0) {
123-
factorize(m / 2, res);
124-
res.push_back(2);
125-
} else if(is_prime(m)) {
126-
res.push_back(m);
127-
} else if(m > 1) {
128-
using base = dynamic_modint;
129-
base::with_mod(m, [&]() {
130-
base t = random::rng();
131-
auto f = [&](auto x) {
132-
return x * x + t;
133-
};
134-
base x, y;
135-
base g = 1;
136-
while(g == 1) {
137-
for(int i = 0; i < 64; i++) {
138-
x = f(x);
139-
y = f(f(y));
140-
if(x == y) [[unlikely]] {
141-
t = random::rng();
142-
x = y = 0;
143-
} else {
144-
base t = g * (x - y);
145-
g = t == 0 ? g : t;
146-
}
147-
}
148-
g = std::gcd(g.getr(), m);
149-
}
150-
factorize(g.getr(), res);
151-
factorize(m / g.getr(), res);
152-
});
153-
}
154-
}
155-
std::vector<int64_t> factorize(int64_t m) {
156-
std::vector<int64_t> res;
157-
factorize(m, res);
158-
return res;
159-
}
16092
int64_t primitive_root(int64_t p) {
16193
using base = dynamic_modint;
16294
return base::with_mod(p, [p](){

‎cp-algo/math/primality.hpp

Lines changed: 33 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,33 @@
1+
#ifndef CP_ALGO_MATH_PRIMALITY_HPP
2+
#define CP_ALGO_MATH_PRIMALITY_HPP
3+
#include "modint.hpp"
4+
#include <algorithm>
5+
namespace cp_algo::math {
6+
// https://en.wikipedia.org/wiki/Miller–Rabin_primality_test
7+
bool is_prime(uint64_t m) {
8+
if(m == 1 || m % 2 == 0) {
9+
return m == 2;
10+
}
11+
// m - 1 = 2^s * d
12+
int s = std::countr_zero(m - 1);
13+
auto d = (m - 1) >> s;
14+
using base = dynamic_modint;
15+
auto test = [&](base x) {
16+
x = bpow(x, d);
17+
if(std::abs(x.rem()) <= 1) {
18+
return true;
19+
}
20+
for(int i = 1; i < s && x != -1; i++) {
21+
x *= x;
22+
}
23+
return x == -1;
24+
};
25+
return base::with_mod(m, [&](){
26+
// Works for all m < 2^64: https://miller-rabin.appspot.com
27+
return std::ranges::all_of(std::array{
28+
2, 325, 9375, 28178, 450775, 9780504, 1795265022
29+
}, test);
30+
});
31+
}
32+
}
33+
#endif // CP_ALGO_MATH_PRIMALITY_HPP

‎verify/number_theory/factorize.test.cpp

Lines changed: 1 addition & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -1,8 +1,7 @@
11
// @brief Factorize
22
#define PROBLEM "https://judge.yosupo.jp/problem/factorize"
33
#pragma GCC optimize("Ofast,unroll-loops")
4-
#pragma GCC target("tune=native")
5-
#include "cp-algo/math/number_theory.hpp"
4+
#include "cp-algo/math/factorize.hpp"
65
#include <bits/stdc++.h>
76

87
using namespace std;

‎verify/number_theory/primality.test.cpp

Lines changed: 1 addition & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -1,8 +1,7 @@
11
// @brief Primality Test
22
#define PROBLEM "https://judge.yosupo.jp/problem/primality_test"
33
#pragma GCC optimize("Ofast,unroll-loops")
4-
#pragma GCC target("tune=native")
5-
#include "cp-algo/math/number_theory.hpp"
4+
#include "cp-algo/math/primality.hpp"
65
#include <bits/stdc++.h>
76

87
using namespace std;

0 commit comments

Comments
 (0)