CSP-S 难度的数学方法讲解：分解质因数的运用、欧拉函数、逆元、组合数

一、分解质因数的运用

我们知道，根据唯一分解定理，对于一个数 $a\in \mathbf{Z}$ ，满足存在唯一的 $$a=\prod _{i=1}^{N_a}{p}_i^{{\alpha }_i}$$ 其中 $N_a$ 为 $a$ 的质因数个数， $\forall 1\le i\le N_a,p_i$ 都是质数。这个性质不仅可以用来分解质因数，也可以用于解决其他所有与质因数有关的问题。由质因数的性质，我们可以推导出几个重要推论。

1. 一个数的因数个数

结论：假设一个数的因数个数为 $\text{fac}(x)$ ，则 $$\text{fac}(x)=\prod _{i=1}^{N_x}({\alpha }_i+1)$$ 。

证明：考虑将这个数的每个因数分别分解质因数。可以知道，这个数的质因数集合一定是 $x$ 的质因数集合的子集。用 $x$ 的质因数分解表示 $x$ 的所有因数为：$1,p_1,\cdots ,p_1^{{\alpha }_1},p_2,\cdots ,p_1{p}_2,p_1^2{p}_2,\cdots$ 。因此，根据乘法原理，对于 $x$ 的每个质因子 $p_i$，可以取 $0\sim {\alpha }_i$ 个与所取其他质因子相乘组成 $x$ 的一个因数，即 $x$ 的质因子个数 $\text{fac}(x)=\prod _{i=1}^{N_x}({\alpha }_i+1)$ 。

求一个数的因数个数，代码如下。

template<typename _Tp>
inline long long countFactors(const _Tp& x, const int mod = 998644353)
{
	vector<int> a(sqrt(x) + 1);
	for (int i = 2; i <= sqrt(x); ++i)
		if (x % i == 0)
			while (x % i == 0)
			{
				a[i]++;
				x /= i;
			}
	long long ans = 1;
	for (int i = 2; i <= sqrt(x); ++i)
		ans = (ans * a[i]) % mod;
	return ans;
}

时间复杂度： $O(\sqrt{x}\cdot \log x)$
推论： $$\text{fac}(n!)=\sum _{T=1}^{N_{n!}}\sum _{i=1}^{{\alpha }_{p_T}}\lfloor \frac{n!}{p_T^i}\rfloor$$ 。

2. 一个数的因数之和

经过上面的证明，我们应该不难得出这个结论：设一个数的因数之和为 $\text{sumfac}(x)$ ，则 $$\text{sumfac}(x)=\prod _{i=1}^{N_x}\sum _{j=0}^{{\alpha }_i}{p}_i^j$$

证明：因为每个因数可以取 $0\sim {\alpha }_i$ 个，所以分别选取 $0\sim {\alpha }_i$ 个相乘即可。该公式还可以根据等比数列的求和公式进一步优化为 $\text{sumfac}(x)=\prod _{i=1}^{N_x}\frac{p_i^{{\alpha }_i+1}-p_i}{p_i-1}$ 。

求一个数的因数之和，代码如下。

template<typename _Tp>
inline long long countFactors(const _Tp& x, const int mod = 998644353)
{
	vector<int> a(sqrt(x) + 1);
	for (int i = 2; i <= sqrt(x); ++i)
		if (x % i == 0)
			while (x % i == 0)
			{
				a[i]++;
				x /= i;
			}
	long long ans = 1;
	for (int i = 2; i <= sqrt(x); ++i)
	{
		long long fac = 0;
		for (int j = 0; j < a[i]; ++j)
			fac += pow(p[i], j);
		ans = ans * fac % mod;
	}
		
	return ans;
}

时间复杂度： $O(x)$ ，并可以根据等比数列求和公式和逆元继续优化至 $O(\sqrt{x}\cdot \log x)$ 。

3. 其他运用

例题： 给定 $n\in \mathbf{Z}$ ，求 $(x,y)$ 的组数。其中 ， $(x,y)$ 满足 $$\frac{1}{x}+\frac{1}{y}=\frac{1}{n!}$$ 。
方程变形为： $\frac{x+y}{xy}=\frac{1}{n!}$ ，分离变量得 $y=\frac{xn!}{x-n!}=n!+\frac{n{!}^2}{x-n!}$ ，因此我们只需要求出 $n{!}^2$ 的因数数量即可。

#include <bits/stdc++.h>

using namespace std;

const int MOD = 1e9 + 7;

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

    vector<int> primes;
    vector<int> vis(n + 1);
    for (int i = 2; i <= n; ++i)
    {
        if (!vis[i]) primes.push_back(i);
        for (int j = 0; primes[j] * i <= n; ++j)
        {
            vis[primes[j] * i] = 1;
            if (i % primes[j] == 0) break;
        }
    }
    long long ans = 1;
    for (int i = 0; i < primes.size(); ++i)
    {
        int x = n;
        long long res = 0;
        while (x)
        {
            res += x / primes[i];
            x /= primes[i];
        }
        ans = (ans * (2 * res + 1)) % MOD;
    }

    cout << ans << endl;

    return 0;
}

时间复杂度： $O(\phi (n))$ 。

二、欧拉函数

欧拉函数是一个积性函数，定义如下：
$$\phi (x)=\sum _{i=1}^x[\gcd (i,x)=1]$$
也就是说，欧拉函数等于小于一个数并且与这个数互质的数的个数。
使用数学方法计算 $\phi$ 的公式：
$$\phi (x)=x\cdot \prod _{i=1}^{N_x}\frac{p_i-1}{p_i}$$

1. 使用欧拉筛求 $\phi$

因此， $\phi$ 可以使用欧拉筛筛出：

inline void getPhi(const int& n)
{
	vector<int> primes(n + 1);
	vector<int> vis(n + 1);
	vector<double> phi(n + 1);
	phi[1] = 1;
	for (int i = 2; i <= n; ++i)
	{
		phi[i] *= i;
		if (!vis[i]) 
		{
			primes.push_back(i);
		}
		for (int j = 0; primes[i] * j <= n; ++j)
		{
			phi[primes[j] * i] *= (i - 1) / i;
			vis[primes[j] * i] = 1;
			if (i % primes[j] == 0) break;
		}
	}
}

时间复杂度： $O(n)$ 。

2. 求一个数的 $\phi$

若只需要求出一个数的 $\phi$ ，也可以使用以下方法：

inline int getPhi(const long long& n)
{
	long long ans = n;
	for (int i = 2; i * i <= n; ++i)
	{
		if (n % i == 0)
		{
			while (n % i == 0) n /= i;
			ans -= ans / i;
		}
	}
	return n > 1 ? ans - ans / n : ans;
}

时间复杂度： $O(\sqrt{n}\cdot \log n)$ 。

三、逆元

在使用模运算的时候，由于模运算不存在分配率，所以我们不能使用除法去进行取模，即
$$\frac{a}{b}\text{ mod }p\ne \frac{a\text{ mod }p}{b\text{ mod }p}$$ 。
而逆元的定义是对于 $a,b,p\in \mathbf{Z}$ ，存在一个整数 $x$ ，使得 $\frac{a}{x}\equiv bx(\mathrm{mod}p)$ 。具体可以使用 Fermat 小定理求出。
Fermat 小定理：对于整数 $n$ 和质数 $p$ ， $p|n^p-n$ 恒成立。
变形，得 $n^p\equiv n(\mathrm{mod}p)$ ， $n^{p-2}\equiv n^{-1}(\mathrm{mod}p)$ 。因此， $n^{p-2}$ 即为 $n$ 的逆元。

inline int getInv(long long x, long long p)
{
	const auto quickPow = [&](long long a, long long b, long long p) 
	{
		long long res = 1;
		for (; b; b >>= 1)
		{
			if (b & 1) res = (res * a) % p;
			a = (a * a) % p;
		}
		return res;
	};
	return quickPow(x, p - 2, p);
}

时间复杂度： $O(\log p)$ 。
条件： $p$ 必须是质数。

四、组合数

组合数的定义： ${\mathrm{C}}_m^n$ 表示从 $m$ 个物品里面拿出 $n$ 个物品进行排列的总方案数。

1. 杨辉三角求组合数

杨辉三角示例如下：
$$\begin{array}{ccccccccc}& & & & 1& & & & \\ & & & 1& & 1& & & \\ & & 1& & 2& & 1& & \\ & 1& & 3& & 3& & 1& \\ 1& & 4& & 6& & 4& & 1\\ & & & & \dots & & & & \end{array}$$
其中，第 $m$ 行第 $n$ 列的数 $(n\le m+1)$ 即为 ${\mathrm{C}}_m^n$ 。
状态转移方程： ${\mathrm{C}}_m^n={\mathrm{C}}_{m-1}^n+{\mathrm{C}}_{m-1}^{n-1}$ 。

template<typename _Tp>
inline int getC(_Tp& C)
{
	C[0][0] = 1;
	for (int i = 1; i <= 1000; ++i)
	{
		for (int j = 0; j <= i; ++j)
		{
			C[i][j] = C[i - 1][j] + C[i - 1][j - 1];
		}
	}
}

时间复杂度： $O(n^2)$
适用范围：$n\le 10^4$

2. 组合数公式求组合数

组合数公式如下：
$${\mathrm{C}}_m^n=\frac{n!}{m!(n-m)!}$$
又因为 $n!,m!,(n-m)!$ 可以通过分解质因数和快速幂计算，所以比杨辉三角更加高效。

#include <bits/stdc++.h>

using namespace std;

const int MOD = 1e9 + 7;

int main()
{
    ios::sync_with_stdio(0);
    int N;
    cin >> N;

    const auto quickPow = [&](long long a, long long b, long long p) -> long long
    {
        long long res = 1;
        for (; b; b >>= 1)
        {
            if (b & 1) res = (res * a) % p;
            a = (a * a) % p;
        }
        return res;
    };
    vector<long long> fact(1e5 + 1);
    vector<long long> inv(1e5 + 1);
    fact[0] = inv[0] = 1;
    for (int i = 1; i <= 1e5; ++i)
    {
        fact[i] = fact[i - 1] * i % MOD;
        inv[i] = inv[i - 1] * quickPow(i, MOD - 2, MOD) % MOD;
    }

    while (N--)
    {
        long long a, b;
        cin >> a >> b;
        cout << ((fact[a] * inv[b]) % MOD * inv[a - b]) % MOD << endl;
    }

    return 0;
}

时间复杂度： $O(n\log n)$ 。适用于 $n\le 10^5$ 的题目。

3. 使用 $\text{Lucas}$ 定理求组合数

$\text{Lucas}$ 定理描述如下：
$${\mathrm{C}}_m^n\equiv {\mathrm{C}}_{n/p}^{m/p}\cdot {\mathrm{C}}_{m\text{ mod }p}^{n\text{ mod }p}(\mathrm{mod}p)$$
使用逆元和快速幂可以高效求出 ${\text{C}}_m^n$ 。

#include <bits/stdc++.h>

using namespace std;

int main()
{
    int N;
    cin >> N;

    const auto quickPow = [&](long long a, long long b, long long p) -> long long
    {
        long long res = 1;
        for (; b; b >>= 1)
        {
            if (b & 1) res = (res * a) % p;
            a = (a * a) % p;
        }
        return res;
    };
    const auto C = [&](long long n, long long m, long long p) -> long long
    {
        if (n < m) return 0;
        if (m == 0) return 1;
        m = (m > n - m ? n - m : m);
        long long up = 1, down = 1;
        for (int i = 1; i <= m; ++i)
            down = (down * i) % p;
        for (int i = n - m + 1; i <= n; ++i)
            up = (up * i) % p;
        return up * quickPow(down, p - 2, p) % p;
    };
    function<long long(long long, long long, long long)> Lucas = [&](long long n, long long m, long long p) -> long long
    {
        if (m == 0) return 1;
        return Lucas(n / p, m / p, p) * C(n % p, m % p, p) % p;
    };

    while (N--)
    {
        long long n, m, p;
        cin >> n >> m >> p;
        cout << Lucas(n, m, p) << endl;
    }
    return 0;
}

注：本文中所有代码未经调试，如有错误请在评论区中指出，谢谢。