数论

约 713 字大约 2 分钟

2026-01-19

警告

本模板为了便于理解，包含了大量的详细注释以及规范的 std:: 命名空间前缀。在您将代码复制并提交至相关竞赛或评测平台前，请务必根据个人的编程习惯进行调整，重点建议删除冗余注释并精简代码风格。保留明显的模板特征极易被系统的查重算法或 AIGC 监测机制误标为“人工智能生成内容”，从而对您的评审进度或最终成绩产生不利影响。

i64 下的素数测试与因式分解（判定是否为素数）

Factorial.cc
using i64 = long long;
i64 mul(i64 a, i64 b, i64 m)
{
    return static_cast<__int128>(a) * b % m;
}
i64 power(i64 a, i64 b, i64 m)
{
    i64 res = 1 % m;
    for (; b; b >>= 1, a = mul(a, a, m))
        if (b & 1)
            res = mul(res, a, m);
    return res;
}
bool isprime(i64 n)
{
    if (n < 2)
        return false;
    static constexpr int A[] = {2, 3, 5, 7, 11, 13, 17, 19, 23};
    int s = __builtin_ctzll(n - 1);
    i64 d = (n - 1) >> s;
    for (auto a : A)
    {
        if (a == n)
            return true;
        i64 x = power(a, d, n);
        if (x == 1 || x == n - 1)
            continue;
        bool ok = false;
        for (int i = 0; i < s - 1; ++i)
        {
            x = mul(x, x, n);
            if (x == n - 1)
            {
                ok = true;
                break;
            }
        }
        if (!ok)
            return false;
    }
    return true;
}
std::vector<i64> factorize(i64 n)
{
    std::vector<i64> p;
    std::function<void(i64)> f = [&](i64 n)
    {
        if (n <= 10000)
        {
            for (int i = 2; i * i <= n; ++i)
                for (; n % i == 0; n /= i)
                    p.push_back(i);
            if (n > 1)
                p.push_back(n);
            return;
        }
        if (isprime(n))
        {
            p.push_back(n);
            return;
        }
        auto g = [&](i64 x)
        {
            return (mul(x, x, n) + 1) % n;
        };
        i64 x0 = 2;
        while (true)
        {
            i64 x = x0;
            i64 y = x0;
            i64 d = 1;
            i64 power = 1, lam = 0;
            i64 v = 1;
            while (d == 1)
            {
                y = g(y);
                ++lam;
                v = mul(v, std::abs(x - y), n);
                if (lam % 127 == 0)
                {
                    d = std::gcd(v, n);
                    v = 1;
                }
                if (power == lam)
                {
                    x = y;
                    power *= 2;
                    lam = 0;
                    d = std::gcd(v, n);
                    v = 1;
                }
            }
            if (d != n)
            {
                f(d);
                f(n / d);
                return;
            }
            ++x0;
        }
    };
    f(n);
    std::sort(p.begin(), p.end());
    return p;
}

快速幂

binpow.cc
using i64 = long long;
i64 binpow(i64 a, i64 b, i64 p)
{
    i64 res = 1;
    while (b > 0)
    {
        if (b & 1)
            res = res * a % p;
        a = a * a % p;
        b >>= 1;
    }
    return res;
}

欧拉筛素数

暴力枚举最多要枚举 $O(\sqrt{n})$ 次，欧拉筛可以根据倍数快速确定是不是素数，其本质是在 $O(n)$ 下生成 $1\sim n$ 的所有质数。

sieve.cc
std::vector<int> minp, primes;

void sieve(int n)
{
    minp.assign(n + 1, 0);
    primes.clear();

    for (int i = 2; i <= n; i++)
    {
        if (minp[i] == 0)
        {
            minp[i] = i;
            primes.push_back(i);
        }

        for (auto p : primes)
        {
            if (i * p > n)
            {
                break;
            }
            minp[i * p] = p;
            if (p == minp[i])
            {
                break;
            }
        }
    }
}

bool isprime(int n)
{
    return minp[n] == n;
}

模逆元

解决的问题是在模数下计算除法的时候将其转化为等价的乘法，即计算

\frac{a}{b}\equiv x ~\text{mod}~ p \Longrightarrow ab^{p-2}\equiv x ~\text{mod}~ p

inv.cc
using i64 = long long;
constexpr int p = 998244353;
i64 binpow(i64 a, i64 b, i64 p)
{
    i64 res = 1;
    while (b > 0)
    {
        if (b & 1)
            res = res * a % p;
        a = a * a % p;
        b >>= 1;
    }
    return res;
}
i64 inv(i64 a, i64 p)
{
    return binpow(a, p - 2, p);
}