组合数学

扩展卢卡斯定理

求 \(c_n^m \bmod p ^ k\)（\(p\) 为质数，\(p^k \leq 10^6\) ）。

---

\(c_n^m = \frac{n!}{m! \times (n - m)!}\)

因为 \(m!\) 和 \((n - m)!\) 不一定和 \(p ^ k\) 互质，所以不一定有逆元。

但如果可以把阶乘中所有质因数 \(p\) 单独提出来计算，剩下就可以用逆元计算。

---

\(g_n\) 表示 \(n!\) 含有几个质因数 \(p\)。

\(g_n = \lfloor n/p \rfloor + g_{\lfloor n / p \rfloor}\)，\(O(\log n)\) 复杂度计算。

---

因为组合数是整数，所以 \(g_n \geq g_m + g_{n - m}\)。

---

\(f_n\) 表示 \(n!\) 去掉所有质因数 \(p\) 后对 \(p^k\) 取模的结果。

\(f_n = f_{\lfloor n / p \rfloor} \times (\prod_{ i=0 }^{ p^k }i [i \not \equiv 0 \pmod p]) ^ {\lfloor n/p^k \rfloor} \times (\prod_{ i=0 }^{n \bmod p^k}i [i \not \equiv 0 \pmod p]) \bmod p^k\)。

第一项递归，后两项在 \(p^k\) 范围内，预处理即可。

---

\(C_n^m\)，不想写了。

点击查看代码

#include <cstdio>
#include <algorithm>
#include <vector>

long long qpow(long long a, long long n, long long mod = 1LL << 60)
{
	a %= mod; long long ans = 1, t = a;
	while (n > 0) ans = n & 1 ? ans * t % mod : ans, t = t * t % mod, n >>= 1;
	return ans;
}

const int maxn = 1e6 + 10;

long long f[maxn];

long long exgcd(long long a, long long b, long long &x, long long &y)
{
	if (b == 0) return x = 1, y = 0, a;
	long long d = exgcd(b, a % b, y, x);
	return y -= (a / b) * x, d;
}

long long calc(long long a, long long b, long long mod)
{
	long long x, y;
	long long t = exgcd(a, -mod, x, y);
	if (b % t != 0) return -1;
	long long d = std::abs(mod / t);
	return ((b / t * x) % d + d) % d;
}

long long excrt(int n, long long a[], long long b[])
{
	for (int i = 1; i < n; i++)
	{
		long long x = calc(a[i], b[i + 1] - b[i], a[i + 1]);
		a[i + 1] = a[i] * a[i + 1] / std::__gcd(a[i], a[i + 1]);
		b[i + 1] = (x * a[i] + b[i]) % a[i + 1];
	}
	return b[n];
}

std::pair<int, long long> exlucas(long long n, long long p, int k, long long pk)
{
	if (n < p) return {0, f[n]};
	auto t = exlucas(n / p, p, k, pk);
	return {t.first + n / p, t.second * qpow(f[pk - 1], n / pk, pk) % pk * f[n % pk] % pk};
}

int main()
{
	long long n, m, p;
	scanf("%lld %lld %lld", &n, &m, &p);
	
	std::vector<std::pair<long long, int>> fac;
	
	long long _p = p;
	for (long long i = 2; i * i <= p; i++)
	{
		if (p % i != 0) continue;
		int cnt = 0;
		while (p % i == 0) cnt++, p /= i;
		fac.push_back({i, cnt});
	}
	if (p > 1) fac.push_back({p, 1});
	p = _p;
	
	long long ans[70] = {0};
	long long mod[70] = {0};
	for (int i = 1; i <= (int)fac.size(); i++)
	{
		long long p = fac[i - 1].first;
		int k = fac[i - 1].second;
		mod[i] = qpow(p, k);
		f[0] = 1;
		for (int j = 1; j < mod[i]; j++) f[j] = j % p == 0 ? f[j - 1] : f[j - 1] * j % mod[i];
		auto x = exlucas(n, p, k, mod[i]), y = exlucas(m, p, k, mod[i]), z = exlucas(n - m, p, k, mod[i]);
		ans[i] = qpow(p, x.first - y.first - z.first, mod[i]) * x.second % mod[i] * calc(y.second, 1, mod[i]) % mod[i] * calc(z.second, 1, mod[i]) % mod[i];
	}
	
	printf("%lld", excrt(fac.size(), mod, ans));
	return 0;
}

计数数列

卡特兰数

小球与盒子

【日报】再谈小球与盒子 - 囧仙 - 洛谷专栏

P5824 十二重计数法（然而需要 FFT）

隔板法

第二类斯特林数与贝尔数

自然数分解

伯努利数