I usedstd::multiset to determine how many times that factor appears.I'd like to know if it can be improved; I've focused on the algorithm rather than the code.
https://en.wikipedia.org/wiki/Pollard%27s_rho_algorithm
std::multiset<std::int64_t> pollards_rho(std::int64_t n){ std::multiset<std::int64_t> factors{}; if (n == 0 || n == 1) { return factors; } if (n < 0) { factors.insert(-1); n = -n; } // If n is even, add 2 and get rid of all even factors. if (n % 2 == 0) { factors.insert(2); while (n % 2 == 0) { n /= 2; } if (n == 1) { return factors; } } // Select an x_0 uniformly at random from [2, n - 1]. // x_0 is chosen as 2 for simplicity. // // Floyd's cycle-finding algorithm. // x => x_i // y => x_i+1 std::int64_t x{2}, y{x}, d{1}; // f(x_i+1) = x_i^2 + 1 mod n auto f = [n](const std::int64_t num) noexcept { return (num * num + 1) % n; }; // If we obtain a factor greater than 1, we are done. while (d == 1) { x = f(x); y = f(f(y)); d = std::gcd(std::abs(x - y), n); } /** * if (d == n) * { * return std::nullopt; * } * else * { * return d; * } */ // The algorithm normally returns either a non-trivial factor d or failure, // but we need to find all prime factors. if (d == n) { // If d and n are the same, we have a prime number. factors.insert(n); } else { // Look for other factors. factors.merge(pollards_rho(d)); factors.merge(pollards_rho(n / d)); } return factors;}
1 Answer1
I used
std::multisetto determine how many times that factor appears
It's an option. You could also use a plain old map to store a factor (key) and its exponent (value):
- In case we only care about factors and not their exponents, the exponents are easier to ignore that way.
- In case we do care about the exponents of factors, we just get them, instead of having to count how many times a factor appears in a multiset.
- "High" exponents don't cost as much duplication, but really an exponent can never be high, at worst you could have 63 copies of
2in your multiset. - a
mapis a less rare/unusual data structure thanmultiset.
// If n is even, add 2 and get rid of all even factors.if (n % 2 == 0){ factors.insert(2); while (n % 2 == 0) { n /= 2; } if (n == 1) { return factors; }}
This does not determine how many times the factor of 2 appears, it's only added once. By the way with C++20 you could usestd::countr_zero instead of a loop, that combines well with returning a map instead of a multiset.
// f(x_i+1) = x_i^2 + 1 mod nauto f = [n](const std::int64_t num) noexcept{ return (num * num + 1) % n;};
This does not handle the case whennum * num overflows, which it can easily do. You could fix this by restricting the input to fit in anint32_t (but still do this arithmetic inf withint64_t - the square of a 32-bit integer fits in 64 bits so it works), or implement a proper 64-bit modular multiplication, that's a bit tricky but possible.
// The algorithm normally returns either a non-trivial factor d or failure,// but we need to find all prime factors.
Careful: Pollard's Rho can also spuriously fail to find a factor when its input is not prime.
// If d and n are the same, we have a prime number.
So this is a bug. You can getd == n whenn is actually composite, but your code would decline to factor it.
- \$\begingroup\$"You can get d == n when n is actually composite, but your code would decline to factor it." Is this the limitation of the algorithm or did I implemented it badly?\$\endgroup\$ozan– ozan2024-08-18 15:11:28 +00:00CommentedAug 18, 2024 at 15:11
- 1\$\begingroup\$@ozan the algorithm inherently works that way\$\endgroup\$user555045– user5550452024-08-18 15:20:31 +00:00CommentedAug 18, 2024 at 15:20
- \$\begingroup\$I see, I also fixed adding factors of 2 you mentioned, looking for the other issues. Thanks.\$\endgroup\$ozan– ozan2024-08-18 16:06:42 +00:00CommentedAug 18, 2024 at 16:06
You mustlog in to answer this question.
Explore related questions
See similar questions with these tags.