分解大质数模板(复杂度小于sqrt(n))

//POJ 1811
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <vector>
#include <time.h>

using namespace std;

typedef __int64 lld;

lld ran() {
    return rand() << 16 | rand();
}

lld gcd(lld a, lld b) {
    return !b ? a : gcd(b, a % b);
}
inline void add(lld &x, lld ad, lld mod) {
    x += ad;
    if (x >= mod) x -= mod;
}
lld mul_mod(lld a, lld b, lld mod) {
    lld ret = 0;
    while (b) {
        if (b & 1) {
            add(ret, a, mod);
        }
        b >>= 1; add(a, a, mod);
    }
    return ret;
}
lld pow_mod(lld x, lld n, lld mod) {
    lld ret = 1 % mod;
    while (n) {
        if (n & 1) {
            ret = mul_mod(ret, x, mod);
        }
        n >>= 1; x = mul_mod(x, x, mod);
    }
    return ret;
}
bool test(lld n, lld b) {
    lld m = n - 1;
    int counter = 0;
    while (~m & 1) {
        m >>= 1;
        counter ++;
    }
    lld ret = pow_mod(b, m, n);
    if (ret == 1 || ret == n - 1) {
        return true;
    }
    counter --;
    while (counter >= 0) {
        ret = mul_mod(ret, ret, n);
        if (ret == n - 1) {
            return true;
        }
        counter --;
    }
    return false;
}
const int BASE[12] = {2,3,5,7,11,13,17,19,23,29,31,37};
bool is_prime(lld n) {
    if (n < 2) {
        return false;
    }
    if (n < 4) {
        return true;
    }
    if (n == 3215031751LL) {
        return false;
    }
    for (int i = 0; i < 12 && BASE[i] < n; i++) {
        if (!test(n, BASE[i])) {
            return false;
        }
    }
    return true;
}
lld pollard_rho(lld n, lld seed) {
    lld x, y, head = 1, tail = 2;
    x = y = ran() % (n - 1) + 1;
    while (true) {
        x = mul_mod(x, x, n);
        add(x, seed, n);
        if (x == y) {
            return n;
        }
        lld d = gcd(x > y ? x - y : y - x, n);
        if (1 < d && d < n) {
            return d;
        }
        head ++;
        if (head == tail) {
            y = x;
            tail <<= 1;
        }
    }
}
vector <lld> divisors;
void factorize(lld n) {
    if (n > 1) {
        if (is_prime(n)) {
            divisors.push_back(n);
        }else {
            lld d = n;
            while (d >= n) {
                d = pollard_rho(n, ran() % (n - 1) + 1);
            }
            factorize(n / d);
            factorize(d);
        }
    }
}

int main() {
    //srand(time(NULL));
    int T;
    scanf("%d", &T);
    for (int cas = 1; cas <= T; cas++) {
        lld x;
        scanf("%I64d", &x);
        if (is_prime(x)) {
            printf("Prime
");
        }else {
            divisors.clear();
            factorize(x);
            sort(divisors.begin(), divisors.end());
            printf("%I64d
", divisors[0]);
        }
    }
    return 0;
}

  

原文地址:https://www.cnblogs.com/get-an-AC-everyday/p/4764141.html