Pollard's rho algorithm implementation

Question

I used std::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;
}

user555045 · Accepted Answer · 2024-08-18 14:28:06Z

I used std::multiset to 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 2 in your multiset.
a map is a less rare/unusual data structure than multiset.

// 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 use std::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 n
auto f = [n](const std::int64_t num) noexcept
{
    return (num * num + 1) % n;
};

This does not handle the case when num * num overflows, which it can easily do. You could fix this by restricting the input to fit in an int32_t (but still do this arithmetic in f with int64_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 get d == n when n is actually composite, but your code would decline to factor it.

"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? — ozan
– ozan, Commented Aug 18, 2024 at 15:11
I see, I also fixed adding factors of 2 you mentioned, looking for the other issues. Thanks. — ozan
– ozan, Commented Aug 18, 2024 at 16:06

Stack Exchange Network

Pollard's rho algorithm implementation

1 Answer 1

You must log in to answer this question.

Hot Network Questions

Pollard's rho algorithm implementation

1 Answer 1

You must log in to answer this question.

Related

Hot Network Questions