POJ - 1811_Prime Test_miller-rabin模板&&polard_rho模板

题意

判断一个数是不是素数。如果不是，输出它的最小因子。

思路

miller_rabin 判素数
polard_rho 因数分解

链接

https://vjudge.net/contest/176263#problem/G

代码

#include<cstdio>
#include<iostream>
#include<cstring>
#include<cstdlib>

using namespace std;

typedef long long LL;

const int maxn = 1e4;
const int S = 20;

LL factor[maxn];
int tot;

LL muti_mod(LL a, LL b, LL c)
{
    a %= c;
    b %= c;
    LL res = 0;

    while(b)
    {
        if(b & 1)
        {
            res += a;
            if(res >= c) res -= c;
        }

        a <<= 1;
        if(a >= c) a -= c;
        b >>= 1;
    }

    return res;
}

LL pow_mod(LL x, LL n, LL mod)
{
    if(n == 1) return x % mod;

    int bit[90], k = 0;

    while(n)
    {
        bit[k ++] = n & 1;
        n >>= 1;
    }

    LL res = 1;
    for(k= k - 1; k>= 0; k --)
    {
        res = muti_mod(res, res, mod);
        if(bit[k] == 1) res = muti_mod(res, x, mod);
    }

    return res;
}

bool check(LL a, LL n, LL x, LL t)
{
    LL res = pow_mod(a, x, n), last = res;

    for(int i= 1; i<= t; i++)
    {
        res = muti_mod(res, res, n);
        if(res == 1 && last != 1 && last != n - 1) return true;
        last = res;
    }

    if(res != 1) return true;
    return false;
}

bool miller_rabin(LL n)
{
    LL x = n - 1, t = 0;
    while((x & 1) == 0) x >>= 1, t ++;

    bool flag = true;
    if(t >= 1 && (x & 1) == 1)
    {
        for(int k= 0; k< S; k++)
        {
            LL a = rand() % (n - 1) + 1;

            if(check(a, n, x, t))
            {
                flag = true;
                break;
            }

            flag = false;
        }
    }

    if(!flag || n == 2) return false;
    return true;
}

LL gcd(LL a, LL b)
{
    if(a == 0) return 1;
    if(a < 0) return gcd(-a, b);

    while(b)
    {
        LL t = a % b;
        a = b;
        b = t;
    }

    return a;
}

LL pollard_rho(LL x, LL c)
{
    LL i = 1, x0 = rand() % x, y = x0, k = 2;

    while(1)
    {
        i ++;
        x0 = (muti_mod(x0, x0, x) + c) % x;
        LL d = gcd(y - x0, x);

        if(d != 1 && d != x) return d;
        if(y == x0) return x;

        if(i == k) y = x0, k += k;
    }
}

void find_factor(LL n)
{
    if(!miller_rabin(n))
    {
        factor[tot ++] = n;
        return;
    }

    LL p = n;
    while(p >= n) p = pollard_rho(p, rand() % (n - 1) + 1);

    find_factor(p);
    find_factor(n / p);
}

int main()
{
    int T;
    scanf("%d", &T);

    while(T --)
    {
        LL n;
        cin >> n;

        if(!miller_rabin(n))
        {
            cout << "Prime\n";
            continue;
        }

        tot = 0;
        find_factor(n);

        LL res = factor[0];

        for(int i= 1; i< tot; i++)
            if(factor[i] < res) res = factor[i];
        cout << res << endl;
    }

    return 0;
}