组合数学
排列组合
排列数：\(A_n^m = n \times (n-1) \times \cdots \times (n-m+1) = \dfrac{n!}{(n-m)!}\)
组合数：\(C_n^m = \dfrac{A_n^m}{m!} = \dfrac{n!}{m!(n-m)!}\)
组合数与杨辉三角：\(C_n^m = C_{n-1}^{m-1} + C_{n-1,m}\)，特殊地，\(C_0^0 = 1, C_i^0 = i\)
组合数常用等式：\(C_n^m = C_n^{n-m}, \ \sum \limits_{i=0}^n C_n^i = 2^n\)
当 \(n,m\) 较小（允许 \(O(nm)\) 做法），模数不一定是质数时，可以使杨辉三角递推计算组合数。
当 \(n,m\) 较大但比模数小（只允许 \(O(n)\) 或 \(O(n \log n)\) 做法）且模数是质数时，可以使用杨辉三角递推计算组合数。
当 \(m\) 很小时，可以用 \(\dfrac{n \times (n-1) \times \cdots \times (n-m+1)}{m!}\) 暴力乘除（如果需要取模且模数是质数，除法用乘以逆元替代）。
习题：P2822 [NOIP 2016 提高组] 组合数问题
解题思路
\(t\) 组数据，\(k\) 都是相同的，所以可以预处理出所有 \(n,m \le 2000\) 组合数 \(\bmod k\) 的结果，本题的 \(k\) 不一定是质数，所以用杨辉三角递推预处理。
为了快速求出一组数据的答案，可以令 \(a_{i,j}\) 表示 \(C_i^j \bmod k\) 是否为 \(0\)，特殊地，如果 \(i \lt j\)，则 \(a_{i,j}=0\)。
对 \(a\) 数组做二维前缀和，单组数据求解的时间复杂度就是 \(O(1)\) 了。
总的时间复杂度为 \(O(nm+t)\)。
参考代码
#include <cstdio>
const int N = 2005;
int c[N][N], a[N][N];
int main()
{
int t, k; scanf("%d%d", &t, &k);
for (int i = 1; i <= 2000; i++) {
c[i][0] = 1; a[i][0] = a[i - 1][0];
for (int j = 1; j < i; j++) {
c[i][j] = (c[i - 1][j] + c[i - 1][j - 1]) % k;
a[i][j] = a[i - 1][j] + a[i][j - 1] - a[i - 1][j - 1] + (c[i][j] == 0);
}
c[i][i] = 1;
for (int j = i; j <= 2000; j++) a[i][j] = a[i][j - 1];
}
for (int i = 1; i <= t; i++) {
int n, m; scanf("%d%d", &n, &m);
printf("%d\n", a[n][m]);
}
return 0;
}
习题：B3717 组合数问题
解题思路
因为模数 \(p\) 是质数，且 \(n,m\) 均比模数小，可以用 \(C_n^m = \dfrac{n!}{m!(n-m)!}\) 求解。
预处理阶乘数组，\(f_i\) 表示 \(i! \bmod p\)，阶乘逆元数组，\(g_i\) 表示 \(i!\) 的乘法逆元。
预处理阶乘逆元时，可以先求 \(f_n\) 的逆元，再倒着求每个阶乘的逆元，\(g_i = g_{i+1} \times (i+1) \bmod p\)，这样求所有逆元的时间复杂度是线性的。
每次询问的答案就是 \(f_n \times g_m \times g_{n-m}\)，注意取模。
总时间复杂度为 \(O(N+T)\)。
参考代码
#include <cstdio>
const int N = 5000005;
const int MOD = 998244353;
int f[N], g[N];
int qpow(int x, int y) {
int res = 1;
while (y > 0) {
if (y & 1) res = 1ll * res * x % MOD;
x = 1ll * x * x % MOD;
y >>= 1;
}
return res;
}
void init(int n) {
f[0] = g[0] = 1;
for (int i = 1; i <= n; i++) f[i] = 1ll * f[i - 1] * i % MOD;
g[n] = qpow(f[n], MOD - 2);
for (int i = n - 1; i > 0; i--) g[i] = 1ll * g[i + 1] * (i + 1) % MOD;
}
int main()
{
int t, maxn;
scanf("%d%d", &t, &maxn);
init(maxn);
int ans = 0;
for (int i = 1; i <= t; i++) {
int n, m; scanf("%d%d", &n, &m);
int res = 1ll * f[n] * g[m] % MOD * g[n - m] % MOD;
ans ^= res;
}
printf("%d\n", ans);
return 0;
}
二项式定理
\((x+y)^n = \overbrace{(x+y) \times (x+y) \times \cdots \times (x+y)}^\text{共 n 项} = C_n^0 x^n y^0 + C_n^1 x^{n-1}y^1 + \cdots + C_n^r x^{n-r}y^r + \cdots + C_n^n x^0 y^n\)，对于 \(C_n^r x^{n-r} y^r\) 这一项来说，相当于有 \(r\) 个 \((x+y)\) 提供了 \(x\)，有 \(n-r\) 个 \((x+y)\) 提供了 \(y\)。
\((1+(-1)^n) = C_n^0 (-1)^0 + C_n^1 (-1)^1 + \cdots = -C_n^0 + C_n^1 - C_n^2 + C_n^3 - \cdots = 0\)，所以 \(\sum C_n^{奇数} = \sum C_n^{偶数} = 2^{n-1}\)。
习题：P1313 [NOIP 2011 提高组] 计算系数
解题思路
把括号展开再合并同类项，\(x^ny^m\) 这一项相当于有 \(n\) 个括号选了 \(ax\) 这一项，其它括号选了 \(by\) 这一项，因此其系数为 \(C_k^n \times a^n \times b^m\)。
本题这个组合数，可以用杨辉三角递推，也可以选择预处理阶乘和阶乘逆元来做。
参考代码
#include <cstdio>
#include <vector>
using std::vector;
const int MOD = 10007;
int main()
{
int a, b, k, n, m;
scanf("%d%d%d%d%d", &a, &b, &k, &n, &m);
int ans = 1;
for (int i = 1; i <= n; i++) ans = 1ll * ans * a % MOD;
for (int i = 1; i <= m; i++) ans = 1ll * ans * b % MOD;
vector<vector<int>> c(k + 1, vector<int>(n + 1));
for (int i = 1; i <= k; i++) {
c[i][0] = 1;
if (i <= n) c[i][i] = 1;
for (int j = 1; j < i && j <= n; j++)
c[i][j] = (c[i - 1][j] + c[i - 1][j - 1]) % MOD;
}
ans = 1ll * ans * c[k][n] % MOD;
printf("%d\n", ans);
return 0;
}
习题：P5390 [Cnoi2019] 数学作业
解题思路
对二进制下每一位分别讨论。
枚举当前位，假设有 \(k\) 个数的二进制在当前位上取 \(1\)。
显然当前位为 \(0\) 的数取不取无所谓，当且仅当选取的子集中有奇数个 \(1\) 的数才对答案有贡献。
\((C_k^1 + C_k^3 + \cdots) \times 2^{n-k} = 2^{k-1+n-k} = 2^{n-1}\)，也就是说只要这些数在第 \(i\) 位上有 \(1\)，就会产生 \(2^{n-i} \times 2^i\) 的贡献。
参考代码
#include <cstdio>
#include <vector>
using std::vector;
const int MOD = 998244353;
void solve() {
int n; scanf("%d", &n);
vector<int> v(n + 1);
int num = 0, p = 1;
for (int i = 1; i <= n; i++) {
scanf("%d", &v[i]); num |= v[i];
if (i > 1) p = 1ll * p * 2 % MOD;
}
int i = 0, ans = 0;
while (num > 0) {
if (num & 1) ans += (1 << i);
num >>= 1; i++;
}
ans = 1ll * ans * p % MOD;
printf("%d\n", ans);
}
int main()
{
int t; scanf("%d", &t);
for (int i = 1; i <= t; i++) solve();
return 0;
}
习题：CF895C Square Subsets
解题思路
一个正整数是完全平方数，当且仅当它的所有质因子的幂次都是偶数。
当选取一个子集并将其元素相乘时，乘积中某个质因子 \(p\) 的幂次，是子集中每个含有的 \(p\) 的幂次之和。
为了让最终乘积是完全平方数，需要让每一个质因子的总幂次为偶数。这等价于，对于每一个质因子 \(p\)，其幂次的奇偶性之和必须是偶数。
这个问题可以转化为一个关于模 2 的代数问题，由于数组元素最大只到 70，只需要考虑 70 以内的质数。这些质数共有 19 个，分别为 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67。
可以为每个质数对应一个“质因子奇偶表示” \(\text{mask}\)，可以看成一个长度为 19 的二进制数。如果某个质数在其质因数分解中幂次为奇数，则对应的那个二进制位为 1，否则为 0。
例如，\(12 = 2^2 \times 3^1\)，70 以内的质数是 \(2, 3, \dots\)。12 的质因子中，2 的幂次是偶数，3 的幂次是奇数。所以 \(\text{mask}(12) = 2 = (10)_2\)。
选择一个子集 \(\{ a_{i_1}, a_{i_2}, \dots \}\)，其乘积是完全平方数，等价于它们对应的 \(\text{mask}\) 异或和为 0，即 \(\text{mask}(a_{i_1}) \oplus \text{mask}(a_{i_2}) \oplus \cdots = 0\)。
问题就变成了从给定的一堆 \(\text{mask}\) 中，选出若干个，使得它们的异或和为 0，求有多少种选法。
由于 \(a_i\) 的值域很小，可以先统计每个数字（1 到 70）在输入数组中出现的次数，记为 \(c_x\)，然后通过动态规划来解决。
定义 \(f_{i,j}\) 表示只考虑数字 \(1\) 到 \(i\) 时，选出一个子集，使得它们的 \(\text{mask}\) 异或和为 \(j\) 的方案数。
当考虑数字 \(i\) 时（假设它在原数组中出现了 \(c_i\) 次），对于 \(f_{i-1}\) 的每个状态 \(f_{i-1,k}\)，可以决定如何从 \(c_i\) 个 \(i\) 中进行选择：
- 如果 \(c_i = 0\)，无法选择数字 \(i\)，所以 \(f_{i,j}=f_{i-1,j}\)。
- 如果 \(c_i \gt 0\)，可以选择奇数个 \(i\) 或者偶数个 \(i\)。
- 选择偶数个 \(i\)：从 \(c_i\) 个 \(i\) 中选出偶数个（包括 0 个）的方案数是 \(C_{c_i}^0 + C_{c_i}^2 + \cdots = 2^{c_i-1}\)，选择偶数个 \(i\)，它们的 \(\text{mask}\) 异或起来为 0，不会改变之前状态的异或和。
- 选择奇数个 \(i\)：从 \(c_i\) 个 \(i\) 中选出奇数个的方案数是 \(C_{c_i,1} + C_{c_i,3} + \cdots = 2^{c_i-1}\)，选择奇数个 \(i\)，它们的 \(\text{mask}\) 异或起来结果为 \(\text{mask}(i)\)，会使之前的异或和 \(k\) 变为 \(k \oplus \text{mask}(i)\)。
因此，\(f_{i,j}\) 的方案数由两部分构成：
- 从 \(f_{i-1,j}\) 转移而来：需要从 \(c_i\) 个 \(i\) 中选择偶数个，方案数为 \(f_{i-1,j} \times 2^{c_i-1}\)。
- 从 \(f_{i-1,j \oplus \text{mask}(i)}\) 转移而来：需要从 \(c_i\) 个 \(i\) 中选择奇数个，方案数为 \(f_{i-1, j \oplus \text{mask}_i} \times 2^{c_i-1}\)。
初始化 \(f_{0,0}=1\)，表示不从任何数字中选择（即空集），异或和为 0，方案数为 1。
最终要求的是 \(f_{70,0}\)，这表示从所有数字中选择，使得 \(\text{mask}\) 异或和为 0 的方案数。这个结果包含了空集（因为 \(f_{0,0}\) 初始化为 1），而题目要求非空子集，所以最终答案还要减 1。
时间复杂度为 \(O(\max \{a\} \cdot 2^L)\)，其中 \(\max \{ a \}\) 是数组中的最大值（\(\le 70\)），\(L\) 是数据范围内的质数个数（\(\le 19\)）。
参考代码
#include <cstdio>
#include <algorithm>
using namespace std;
const int N = 100005;
const int A = 71; // 数组元素的最大值+1
const int MOD = 1000000007;
int a[N]; // 输入数组
int p[A], len; // p: 存储小于70的质数, len: 质数个数
int cnt[A]; // cnt[i]: 数字i在输入中出现的次数
int mask[N]; // mask[i]: 数字i的质因子幂次奇偶性
// dp[cur/pre][j]: DP数组, j为质因子奇偶性(mask)
// dp[i][j]表示考虑完前i种数，mask异或和为j的方案数
int dp[2][1 << 19];
int p2[N]; // p2[i] = 2^i % MOD
bool isprime[A];
int main()
{
// 1. 预处理：筛出70以内的所有质数
for (int i = 2; i < A; i++) isprime[i] = true;
for (int i = 2; i < A; i++) {
if (isprime[i]) {
p[len++] = i; // len = 19
for (int j = i * i; j < A; j += i) isprime[j] = false;
}
}
int n; scanf("%d", &n);
// 预处理2^i
p2[0] = 1; for (int i = 1; i <= n; i++) p2[i] = p2[i - 1] * 2 % MOD;
// 统计每个数字出现的次数
int maxa = 0;
for (int i = 1; i <= n; i++) {
scanf("%d", &a[i]); cnt[a[i]]++; maxa = max(maxa, a[i]);
}
// 2. 预处理：为1到maxa的每个数计算其质因子奇偶性mask
for (int i = 1; i <= maxa; i++) {
int x = i;
for (int j = 0; j < len; j++) {
while (x % p[j] == 0) {
mask[i] ^= (1 << j); // 对应质因子幂次+1，奇偶性翻转
x /= p[j];
}
}
}
// 3. 动态规划
dp[0][0] = 1; // base case: 不选任何数，异或和为0，方案数为1 (空集)
int cur = 0;
// 遍历1到maxa的所有数
for (int i = 1; i <= maxa; i++) {
cur = i % 2; int pre = 1 - cur; // 滚动数组
// 如果数字i没有出现过，直接继承上一轮的结果
if (cnt[i] == 0) {
for (int j = 0; j < (1 << len); j++) dp[cur][j] = dp[pre][j];
continue;
}
// 如果数字i出现过cnt[i]次
for (int j = 0; j < (1 << len); j++) dp[cur][j] = 0; // 清空当前dp数组
// 遍历上一轮的所有mask状态
for (int j = 0; j < (1 << len); j++) {
if (dp[pre][j] == 0) continue; // 优化：如果上一轮该状态方案数为0，则跳过
// 状态转移方程: dp[i][j] = (dp[i-1][j] * ways_even) + (dp[i-1][j^mask(i)] * ways_odd)
// ways_even = ways_odd = 2^(cnt[i]-1)
// a. 从cnt[i]个i中选偶数个，mask不变，贡献给dp[cur][j]
dp[cur][j] = (dp[cur][j] + 1ll * p2[cnt[i] - 1] * dp[pre][j]) % MOD;
// b. 从cnt[i]个i中选奇数个，mask变为j^mask[i]，贡献给dp[cur][j^mask[i]]
dp[cur][j ^ mask[i]] = (dp[cur][j ^ mask[i]] + 1ll * p2[cnt[i] - 1] * dp[pre][j]) % MOD;
}
}
// 4. 输出结果
// dp[cur][0]是包括空集在内的方案数，题目要求非空子集，所以减1
printf("%d\n", (dp[cur][0] + MOD - 1) % MOD);
return 0;
}
组合数常见问题
插板法
\(n\) 个相同的球放进 \(m\) 个不同的盒子（盒子允许空）：\(C_{n+m-1}^n\)，将 \(n\) 个球与 \(m-1\) 个分隔板看成总共 \(n+m-1\) 个物品，其中 \(n\) 个是球。
\(n\) 个相同的球放进 \(m\) 个不同的盒子（盒子不能空）：\(C_{n-1}^{m-1}\)，在 \(n-1\) 个缝隙中选择 \(m-1\) 个放分隔板。
从 \(n\) 个不同元素中选 \(m\) 个互不相邻的元素：\(C_{n-m+1}^m\)，在 \(n-m\) 个元素形成的 \(n-m+1\) 个缝隙（包括头尾）中插进去 \(m\) 个元素。
格线游走
从 \((0,0)\) 出发，每步往右或往上走一个单位，走到 \((n,m)\) 的路径数：\(C_{n+m}^n\)，一共 \(n+m\) 步，其中 \(n\) 步往下。
习题：P14254 分割（divide）
解题思路
核心在于理解并简化等式 \(S_1 = \bigcap \limits_{i=2}^{k+1} S_i\)。
设 \(C = \{ b_1, \dots, b_k \}\) 为选中的节点集合。对于 \(i \in [1,k]\)，子树 \(i\) 的深度集合为 \(S_i\) 是一个连续的整数区间，即 \(S_i = [d_{b_i}, h_{b_i}]\)，其中 \(d_{b_i}\) 是节点 \(b_i\) 的深度，而 \(h_{b_i}\) 是 \(b_i\) 子树中节点的最大深度。题目要求 \(d_{b_1} \le d_{b_2} \le \cdots \le d_{b_k}\)；\(S_1\) 是后面若干集合的交集意味着 \(S_1 \subseteq S_2\)，说明 \([d_{b_1}, h_{b_1}] \le [d_{b_2}, h_{b_2}]\)，这要求 \(d_{b_2} \le d_{b_1}\)。结合两者，可得 \(d_{b_1} = d_{b_2}\)，依此类推，所有被选节点的深度必须相同，设这个共同深度为 \(D\)。
条件 \(S_1 \subseteq S_i\) 对于所有 \(i \in [2,k]\) 都成立，这意味着 \([d_{b_1}, h_{b_1}] \subseteq [d_{b_i}, h_{b_i}]\)，因为深度都为 \(D\)，这简化为 \(h_{b_1} \le h_{b_i}\)。这说明，\(b_1\) 必须是所选节点集合 \(C\) 中，\(h\) 值（子树最大深度）最小的节点之一。
\(S_{k+1}\) 是包含根节点 1 的剩余部分的深度集合。设 \(U_D\) 是深度为 \(D\) 的所有节点的集合，一个深度 \(d \gt D\) 要想出现在 \(S_{k+1}\) 中，当且仅当在 \(U_D\) 中存在一个未被选择的节点 \(u \notin C\)，且 \(u\) 的子树中包含深度为 \(d\) 的节点（即 \(h_u \ge D\)）。定义一个深度 \(d\) 被集合 \(C\) 覆盖，当且仅当所有未被选择的节点 \(u \in U_D \setminus C\) 的子树最大深度都小于 \(d\)（即 \(\forall u \in U_D \setminus C, h_u \lt d\)），反之，则称 \(d\) 未被覆盖。于是，条件 \(d \in S_{k+1}\) 等价于 \(d\) 未被覆盖。
令 \(h_{\min} = \min \limits_{c \in C} h_c\)，而 \(C_{\min} = \{ c \in C \mid h_c = h_{\min} \}\)。
- 情况 A：\(|C_{\min}| \gt 1\)（被选节点中，拥有最小 \(h\) 值的节点不唯一），此时 \(b_1\) 可以是 \(C_{\min}\) 中任意一个，\(S_1 = [D, h_{\min}]\)，要使等式成立，要求对于所有 \(d \in [D, h_{\min}]\)，都有 \(d \in S_{k+1}\)，即所有这些深度都未被覆盖。
- 情况 B：\(|C_{\min}| = 1\)（拥有最小 \(h\) 值的节点唯一），此时 \(b_1\) 必须是这个唯一的节点，设 \(h_2 = \min \limits_{c \in C \setminus C_{\min}} h_c\)（集合中第二小的 \(h\) 值）。等式要求成立：对所有 \(d \in [D, h_{\min}]\)，深度 \(d\) 未被覆盖；且对所有 \(d \in (h_{\min}, h_2]\)，深度 \(d\) 都被覆盖。
基于以上分析，可以设计一个组合计数的算法。
- 预处理：遍历整棵树，计算出每个节点的深度和子树最大深度，同时预处理计算排列组合需要的阶乘及其逆元。
- 按深度枚举：遍历所有可能的深度（从 2 到树的最大深度），对于每个深度 \(D\)，将只在这一层的节点中选取 \(k\) 个。
- 分组与计数
- 对于固定的深度 \(D\)，获取该层所有节点 \(U_D\)，将这些节点按它们的 \(h\) 值进行分组。
- 为了方便计算，可以将这些组按 \(h\) 值排序，设共有 \(p\) 个组 \(G_1, \dots, G_p\)，对应的 \(h\) 值为 \(H_1 \lt H_2 \cdots \lt H_p\)，各组大小为 \(m_1, \dots, m_p\)。
- 遍历每个组 \(G_j\)，假设选出的集合 \(C\) 的 \(h_{\min}\) 值就等于 \(H_j\)。
- 计算情况 A 的贡献
- 条件：\(h_{\min} = H_j\) 且 \(|C_{\min}| \ge 2\)。
- “未被覆盖”条件要求不能选走所有 \(h \gt H_j\) 的节点，设 \(h \gt H_j\) 的节点总数为 \(N_j\)，则要求 \(N_j \gt k\)。
- 贡献的序列总数为 \(m_j \times (A_{N_j-1}^{k-1} - A_{M_{\gt j}}^{k-1})\)，其中 \(M_{\gt j}\) 是 \(h \gt H_j\) 的节点总数，\(A\) 是排列数。
- 计算情况 B 的贡献
- 条件：\(h_{\min} = H_j\) 且 \(|C_{\min}| = 1\)。
- “被覆盖”条件要求必须选走所有 \(h \gt H_j\) 的节点，即 \(k-1 = M_{\gt j}\)。
- “未被覆盖”条件要求 \(m_j \gt 1\)。
- 若满足 \(k - 1 = M_{\gt j}\) 且 \(m_j \gt 1\)，则贡献的序列总数为 \(m_j \times (k-1)!\)。
- 将所有深度、所有分组计算出的贡献累加，即为最终答案。
参考代码
#include <cstdio>
#include <algorithm>
#include <vector>
using namespace std;
// 常量定义
const int N = 1e6 + 5;
const int MOD = 998244353;
// 全局变量
int d[N], h[N], maxd, f[N], g[N]; // d:深度, h:子树最大深度, maxd:全局最大深度, f:阶乘, g:阶乘逆元
vector<int> tr[N], vts[N]; // tr:树的邻接表, vts:按深度存储节点
// DFS遍历，计算每个节点的深度d和子树最大深度h
void dfs(int u, int depth) {
d[u] = depth;
maxd = max(maxd, depth);
vts[depth].push_back(u);
h[u] = depth; // 初始化h值为自身深度
for (int v : tr[u]) {
dfs(v, depth + 1);
h[u] = max(h[u], h[v]); // 用子节点的h值更新父节点的h值
}
}
// 快速幂，用于计算模逆元
int qpow(int x, int y) {
int res = 1;
while (y > 0) {
if (y & 1) res = 1ll * res * x % MOD;
x = 1ll * x * x % MOD;
y >>= 1;
}
return res;
}
// 预处理阶乘和阶乘逆元
void init(int n) {
f[0] = g[0] = 1;
for (int i = 1; i <= n; i++) f[i] = 1ll * f[i - 1] * i % MOD;
g[n] = qpow(f[n], MOD - 2);
for (int i = n - 1; i >= 1; i--) g[i] = 1ll * g[i + 1] * (i + 1) % MOD;
}
// 计算排列数 A(n, m) = n! / (n-m)!
int A(int n, int m) {
if (n < 0 || m < 0 || n - m < 0) return 0;
return 1ll * f[n] * g[n - m] % MOD;
}
int main()
{
int n, k;
scanf("%d%d", &n, &k);
init(n);
d[1] = 1;
for (int i = 2; i <= n; i++) {
int p; scanf("%d", &p);
tr[p].push_back(i);
}
dfs(1, 1);
int ans = 0;
// 主循环：遍历所有可能的深度D
for (int i = 2; i <= maxd; i++) {
if (vts[i].size() < k) continue; // 如果该层节点数不足k，则跳过
// 收集当前深度所有节点的h值并排序，用于分组
vector<int> hi;
for (int vtx : vts[i]) {
hi.push_back(h[vtx]);
}
sort(hi.begin(), hi.end());
int s = hi.size();
int r = s - 1, l = s - 1;
// 从大到小遍历所有h值分组
while (r >= 0) {
while (l >= 0 && hi[l] == hi[r]) l--; // 找到当前h值分组的左边界
// cnt: 当前h值分组的节点数 (m_j)
int cnt = r - l;
// suf: h值大于等于当前h值的节点总数 (N_j)
int suf = s - 1 - l;
// s-1-r: h值严格大于当前h值的节点总数 (M_{>j})
int suf_greater = s - 1 - r;
// 情况A: C_min大小大于1。要求h_min=H_j, |C_min|>=2
// 条件: suf > k (确保未覆盖), cnt >= 2 (确保能选出至少2个)
if (cnt >= 2 && suf > k) {
// 计算 A(N_j-1, k-1) - A(M_{>j}, k-1)
int tmp = (A(suf - 1, k - 1) + MOD - A(suf_greater, k - 1)) % MOD;
// 乘以 m_j, 得到总贡献
ans = (ans + 1ll * cnt * tmp % MOD) % MOD;
}
// 情况B: C_min大小等于1。要求h_min=H_j, |C_min|=1
// 条件: M_{>j} == k-1 (确保覆盖), m_j > 1 (确保未覆盖)
if (suf_greater == k - 1 && cnt > 1) {
// 贡献为 m_j * (k-1)!
ans = (ans + 1ll * cnt * f[k - 1] % MOD) % MOD;
}
r = l; // 移动到下一个分组
}
}
printf("%d\n", ans);
return 0;
}
错排
求有多少 \(1 \sim n\) 的排列 \(p\) 满足对于任意 $i = 1, \dots, n $，有 \(p_i \ne i\)。
方法 1（递推）：\(D_i = (D_{i-1} + D_{i-2}) \times (i-1)\)，其中 \(D_0 = 1, \ D_1 = 0\)
第一种情况，前 \(i-1\) 位形成错排，第 \(i\) 位跟前 \(i-1\) 位里的任意一位交换。
第二种情况，有 \(i-2\) 位形成错排，第 \(i\) 位跟前 \(i-1\) 位里唯一 \(p_i = i\) 的交换。
方法 2（容斥）：\(D_n = \sum \limits_{k=0}^n (-1)^k \times C_n^k \times (n-k)!\)
可以继续简化为 \(D_n = \sum \limits_{k=0}^n (-1)^k \times \dfrac{n!}{k!(n-k)!} \times (n-k)! = n! \times \sum \limits_{k=0}^n \dfrac{(-1)^k}{k!}\)。
其中 \(n!\) 可以从 \((n-1)!\) 递推，后边的求和项也是会比算 \(D_{n-1}\) 时多一项，所以这个 \(D_n\) 的表达式也可以递推。
习题：P4071 [SDOI2016] 排列计数
解题思路
从 \(n\) 个位置中选出 \(m\) 个位置固定 \(a_i = i\)，其它位置进行错排。
因此答案为 \(C_n^m \times D_{n-m}\)（\(D\) 表示错排数）。
组合数使用预处理阶乘和阶乘逆元的方法求。
参考代码
#include <cstdio>
const int N = 1000005;
const int MOD = 1000000007;
int d[N], f[N], g[N];
int qpow(int x, int y) {
int res = 1;
while (y > 0) {
if (y & 1) res = 1ll * res * x % MOD;
x = 1ll * x * x % MOD;
y >>= 1;
}
return res;
}
void init(int n) {
f[0] = 1; g[0] = 1;
for (int i = 1; i <= n; i++) f[i] = 1ll * f[i - 1] * i % MOD;
g[n] = qpow(f[n], MOD - 2);
for (int i = n - 1; i > 0; i--) g[i] = 1ll * g[i + 1] * (i + 1) % MOD;
d[0] = 1; d[1] = 0;
for (int i = 2; i <= n; i++) d[i] = 1ll * (i - 1) * (d[i - 1] + d[i - 2]) % MOD;
}
int C(int n, int m) {
return 1ll * f[n] * g[m] % MOD * g[n - m] % MOD;
}
int main()
{
int t; scanf("%d", &t);
init(1000000);
for (int i = 1; i <= t; i++) {
int n, m; scanf("%d%d", &n, &m);
int ans = m > n ? 0 : 1ll * C(n, m) * d[n - m] % MOD;
printf("%d\n", ans);
}
return 0;
}
卢卡斯定理
引理 1：\(C_{p}^x \equiv 0 \pmod{p}, \ 0 \lt x \lt p\)。因为 \(C_p^x = \dfrac{p!}{x!(p-x)!} = \dfrac{p(p-1)!}{x(x-1)!(p-x)!} = \dfrac{p}{x} C_{p-1}^{x-1}\)，所以 \(C_p^x \equiv p \cdot inv_{x} \cdot C_{p-1}^{x-1} \equiv 0 \pmod{p}\)。
引理 2：\((1+x)^p \equiv 1 + x^p \pmod{p}\)。二项式定理 \((1+x)^p = \sum \limits_{i=0}^p C_p^i x^i\)，由引理 1 得只剩 \(i=0,p\) 两项。
大组合数取模问题：给定整数 \(n,m,p\) 的值，求出 $C_n^m \pmod{p} $ 的值。其中 \(1 \le m \le n \le 10^{18}, \ 1 \le p \le 10^5\)，保证 \(p\) 为质数。
卢卡斯（Lucas）定理：\(C_n^m \equiv C_{n/p}^{m/p} C_{n \bmod p}^{m \bmod p} \pmod{p}\)，其中 \(p\) 为质数。\(n \bmod p\) 和 \(m \bmod p\) 一定是小于 \(p\) 的数，可以直接求解，\(C_{n/p}^{m/p}\) 可以继续用 Lucas 定理求解。
边界条件：当 \(m=0\) 时，结果为 \(1\)。
证明：令 \(n = ap + b, \ m = cp + d\)。
\((1+x)^n \equiv \sum \limits_{i=0}^n C_n^i x^i \pmod{p}\)，其中 \(x^m\) 的系数为 \(C_n^m\)。
\((1+x)^n \equiv (1+x)^{ap+b} \equiv (1+x^p)^a \cdot (1+x)^b \equiv \sum \limits_{i=0}^a C_a^i x^{ip} \cdot \sum \limits_{j=0}^{b} C_b^j x^b \pmod{p}\)，其中 \(x^m = x^{cp} \cdot x^d\) 的系数为 \(C_a^c C_b^d\)。
所以 \(C_n^m \equiv C_a^c C_b^d \pmod{p}\)，即 \(C_n^m \equiv C_{n/p}^{m/p} C_{n \bmod p}^{m \bmod p} \pmod{p}\)。
void init(int p) {
f[0] = 1; g[0] = 1;
for (int i = 1; i < p; i++)
f[i] = 1ll * f[i - 1] * i % p;
g[p - 1] = qpow(f[n - 1], p - 2, p);
for (int i = p - 2; i >= 1; i--)
g[i] = 1ll * g[i + 1] * (i + 1) % p;
}
int C(int n, int m, int p) {
if (n < m) return 0;
return 1ll * f[n] * g[m] % p * g[n - m] % p;
}
int lucas(ll n, ll m, int p) {
if (m == 0) return 1;
return 1ll * lucas(n / p, m / p, p) * C(n % p, m % p, p) % p;
}
时间复杂度：预处理 \(O(p)\)，计算一次组合数 \(O(\log_p n)\)。
习题：P3807 【模板】卢卡斯定理/Lucas 定理
解题思路
因为 \(n,m\) 有可能大于等于 \(p\)，且 \(p\) 为质数，所以使用 Lucas 定理。
需要注意每次 \(p\) 不同，所以需要重新推一下阶乘和阶乘逆元，并且是预处理到 \(p-1\)。
参考代码
#include <cstdio>
const int N = 100005;
int f[N], g[N];
int qpow(int x, int y, int p) {
int res = 1;
while (y > 0) {
if (y & 1) res = 1ll * res * x % p;
x = 1ll * x * x % p;
y >>= 1;
}
return res;
}
void init(int p) {
f[0] = 1; g[0] = 1;
for (int i = 1; i < p; i++)
f[i] = 1ll * f[i - 1] * i % p;
g[p - 1] = qpow(f[p - 1], p - 2, p);
for (int i = p - 2; i > 0; i--)
g[i] = 1ll * g[i + 1] * (i + 1) % p;
}
int C(int n, int m, int p) {
if (n < m) return 0;
return 1ll * f[n] * g[m] % p * g[n - m] % p;
}
int lucas(int n, int m, int p) {
if (m == 0) return 1;
return 1ll * lucas(n / p, m / p, p) * C(n % p, m % p, p) % p;
}
void solve() {
int n, m, p; scanf("%d%d%d", &n, &m, &p);
init(p);
printf("%d\n", lucas(n + m, n, p));
}
int main()
{
int t; scanf("%d", &t);
for (int i = 1; i <= t; i++) solve();
return 0;
}
习题：P2606 [ZJOI2010] 排列计数
解题思路
满足要求的排列相当于一棵二叉树，父亲节点的数小于儿子节点的数（小根堆）。
根肯定是最小数 \(1\)，然后把剩下 \(n-1\) 个数分成左右子树两部分，设左子树大小为 \(left\)，右子树大小为 \(right\)，则有 \(C(n-1,left)\) 种分法，问题就被拆成了大小为 \(left\) 和 \(right\) 的两个子问题。
因此存在一种 DP 方法，\(dp_i = dp_{left_i} \times dp_{i-1-left_i} \times C_{i-1}^{left_i}\)，其中 \(left_i\) 表示当大小为 \(i\) 时，其根的左子树的大小。
因为 \(n\) 有可能大于等于模数，且模数为质数，所以组合数使用 Lucas 定理求。
考虑 \(left\) 怎么求？一个新节点跟它的父节点同处于左子树或右子树，则有起始条件：\(2\) 在左子树，\(3\) 在右子树，后面的节点可以用父节点情况递推。
所以如果 \(i\) 在左子树，\(left_i = left_{i-1} + 1\)，否则 \(left_i = left_{i-1}\)。
参考代码
#include <cstdio>
#include <algorithm>
using std::min;
const int N = 1000005;
bool flag[N]; // true表示在左子树，false表示在右子树
int ans[N], sz_left[N], f[N], g[N];
int qpow(int x, int y, int mod) {
int res = 1;
while (y > 0) {
if (y & 1) res = 1ll * res * x % mod;
x = 1ll * x * x % mod;
y >>= 1;
}
return res;
}
void init(int n, int p) {
f[0] = g[0] = 1;
for (int i = 1; i < n; i++)
f[i] = 1ll * f[i - 1] * i % p;
g[n - 1] = qpow(f[n - 1], p - 2, p);
for (int i = n - 2; i > 0; i--)
g[i] = 1ll * g[i + 1] * (i + 1) % p;
}
int C(int n, int m, int p) {
if (n < m) return 0;
return 1ll * f[n] * g[m] % p * g[n - m] % p;
}
int lucas(int n, int m, int p) {
if (m == 0) return 1;
return 1ll * lucas(n / p, m / p, p) * C(n % p, m % p, p) % p;
}
int main()
{
int n, m; scanf("%d%d", &n, &m);
ans[0] = ans[1] = 1; flag[2] = 1; flag[3] = 0;
ans[2] = 1; ans[3] = 2;
sz_left[2] = 1; sz_left[3] = 1;
init(min(n, m), m);
for (int i = 4; i <= n; i++) {
flag[i] = flag[i / 2];
sz_left[i] = sz_left[i - 1] + flag[i];
ans[i] = 1ll * lucas(i - 1, sz_left[i], m) * ans[sz_left[i]] % m * ans[i - 1 - sz_left[i]] % m;
}
printf("%d\n", ans[n]);
return 0;
}
多重集排列
由 \(n_1\) 个 \(a_1\)，\(n_2\) 个 \(a_2\)，……，\(n_k\) 个 \(a_k\) 组成的多重集，全排列数量为：
- $ \dfrac{(n_1 + \cdots + n_k)!}{n_1! \times n_2! \times \cdots \times n_k!}$，任意排列，再把每种元素内部的顺序除掉
- \(C_{n_1+ \cdots + n_k}^{n_1} \times C_{n_2 + \cdots + n_k}^{n_2} \times \cdots \times C_{n_k}^{n_k}\)，从所有位置中先选出 \(n_1\) 个位置放 \(a_1\)，再在剩下的位置中选出 \(n_2\) 个位置放 \(a_2\)，……
直接计算时一般使用第一种，两种方法都可以与 DP 思路结合。
习题：P7961 [NOIP2021] 数列
解题思路
测试点 \(1 \sim 4\)
直接枚举每个 \(a_i\) 的取值，时间复杂度 \(O(m^n)\)。
参考代码
#include <cstdio>
const int MOD = 998244353;
const int S = 1000005;
const int M = 105;
int n, m, k, ans;
int cnt[S], v[M];
void dfs(int u, int s, int val) {
if (u == n + 1) {
if (cnt[s] <= k) ans = (ans + val) % MOD;
return;
}
for (int i = 0; i <= m; i++) {
// a[u]填i
dfs(u + 1, s + (1 << i), 1ll * val * v[i] % MOD);
}
}
int lowbit(int x) {
return x & -x;
}
int main()
{
for (int i = 1; i < S; i++) cnt[i] = cnt[i - lowbit(i)] + 1;
scanf("%d%d%d", &n, &m, &k);
for (int i = 0; i <= m; i++) scanf("%d", &v[i]);
dfs(1, 0, 1);
printf("%d\n", ans);
return 0;
}
测试点 \(5 \sim 10\)
设 \(dp_{i,j}\) 表示填前 \(i\) 个数，得到的 \(S = j\) 时，所有方案的权值和。假设 \(a_i\) 要填的是 \(k\)，则状态转移为 \(dp_{i,j} = \sum \limits_{k=0}^m dp_{i-1,j-2^k} \times v_k\)。时间复杂度为 \(O(n^2m2^m)\)。
参考代码
#include <cstdio>
const int MOD = 998244353;
const int S = 150005;
const int M = 105;
const int N = 35;
int cnt[S], v[M], dp[N][S];
int lowbit(int x) {
return x & -x;
}
int main()
{
for (int i = 1; i < S; i++) cnt[i] = cnt[i - lowbit(i)] + 1;
int n, m, k; scanf("%d%d%d", &n, &m, &k);
int maxs = (1 << m) * n;
for (int i = 0; i <= m; i++) scanf("%d", &v[i]);
dp[0][0] = 1;
for (int i = 1; i <= n; i++) {
for (int j = 0; j <= maxs; j++) {
for (int k = 0; k <= m; k++) {
if (j >= (1 << k)) {
dp[i][j] += 1ll * dp[i - 1][j - (1 << k)] * v[k] % MOD;
dp[i][j] %= MOD;
}
}
}
}
int ans = 0;
for (int j = n; j <= maxs; j++)
if (cnt[j] <= k) ans = (ans + dp[n][j]) % MOD;
printf("%d\n", ans);
return 0;
}
测试点 \(11 \sim 20\)
因为 \(S\) 的结果只和 \(a\) 中 \(0\) 到 \(m\) 的个数有关和具体顺序无关，而且如果按下标 \(1\) 到 \(n\) 的顺序来做的话需要记录 \(S\) 每一位有几个 \(1\)，状态数太多，所以转而按值域 \(0\) 到 \(m\) 考虑，每次考虑一种数值填在了 \(a\) 数组的几个位置上。
设 \(dp_{i,j,p,q}\) 表示填 \(i\) 这种数（依次枚举 \(0\) 到 \(m\)），在 \(a\) 中已经有 \(j\) 个位置填了数，\(S\) 的二进制低位中已经有了 \(p\) 个 \(1\)，\(S\) 的 \(2^i\) 这一位上已经从低位进位了 \(q\) 过来，这种状态下所有方案的权值之和。
如果考虑前边的状态转移过来（填表法）情况会有点多，从当前状态转出去（刷表法）会更好实现。
考虑填 \(x\) 个 \(i\)，则 \((i,j,p,q)\) 这个状态会转移到 \((i+1,j+x,p+(q+x) \bmod 2, \left\lfloor \dfrac{q+x}{2} \right\rfloor)\)，并且有 \(C_{n-j}^x\) 种填法，填这一次会让权值乘上 \(v_{i}^x\)。其中进位数 \(q\) 那一维不会超过 \(n/2\)（进位最多的情况是把 \(n\) 位都填成同一个数，往前进 \(n/2\) 位）。
最后看所有的 \(dp_{m+1,n,p,q}\)，注意此时这 \(q\) 个进位可能会让 \(S\) 的二进制中 \(1\) 的数量不止 \(p\)，\(S\) 的二进制中 \(1\) 的真正数量应该是 \(p+popcnt(q)\)，在这个结果 \(\le k\) 的情况下才能把 \(dp_{m+1,n,p,q}\) 算进答案。
状态树为 \(m \times n \times k \times n\)，转移的时间复杂度为 \(O(n)\)，总时间复杂度为 \(O(mn^3k)\)。
参考代码
#include <cstdio>
const int MOD = 998244353;
const int N = 35;
const int M = 105;
int dp[M][N][N][N], v[M], c[N][N], power[M][N];
int lowbit(int x) {
return x & -x;
}
int popcnt(int x) {
int res = 0;
while (x > 0) {
x -= lowbit(x);
res++;
}
return res;
}
int main()
{
int n, m, k; scanf("%d%d%d", &n, &m, &k);
c[0][0] = 1;
for (int i = 1; i <= n; i++) {
c[i][0] = c[i][i] = 1;
for (int j = 1; j < i; j++) c[i][j] = (c[i - 1][j] + c[i - 1][j - 1]) % MOD;
}
for (int i = 0; i <= m; i++) {
scanf("%d", &v[i]);
power[i][0] = 1;
for (int j = 1; j <= n; j++) power[i][j] = 1ll * power[i][j - 1] * v[i] % MOD;
}
dp[0][0][0][0] = 1;
for (int i = 0; i <= m; i++)
for (int j = 0; j <= n; j++)
for (int p = 0; p <= k; p++)
for (int q = 0; q <= n / 2; q++)
for (int x = 0; x <= n - j; x++) {
int add = 1ll * dp[i][j][p][q] * c[n - j][x] % MOD * power[i][x] % MOD;
(dp[i + 1][j + x][p + (q + x) % 2][(q + x) / 2] += add) %= MOD;
}
int ans = 0;
for (int p = 0; p <= k; p++)
for (int q = 0; q <= n / 2; q++) {
int cnt = p + popcnt(q);
if (cnt <= k) {
ans += dp[m + 1][n][p][q];
ans %= MOD;
}
}
printf("%d\n", ans);
return 0;
}
卡特兰数
\(H_n = C_{2n}^n - C_{2n}^{n-1} = \dfrac{C_{2n}^n}{n+1} = H_{n-1} \times \dfrac{4n-2}{n+1}\)。
合法括号序列匹配问题，总数减去非法序列数，找到非法括号匹配序列中第一个右括号比左括号多一个的地方，将这个前缀取反（左括号变右括号，右括号变左括号），就会变成一个有 \(n+1\) 个左括号，\(n-1\) 个右括号的括号序列。
并且每一个有 \(n+1\) 个左括号，\(n-1\) 个右括号的括号序列都可以通过找到第一个左括号比右括号多一个的地方，把这个前缀取反，对应回原来的一个非法括号匹配序列，说明非法括号匹配序列和有 \(n+1\) 个左括号，\(n-1\) 个右括号的括号序列一一对应。
如果用格线游走解释的话，非法路径可以通过找到第一个突破 \(y=x\) 这条线的地方，一定是到某个 \((x,x+1)\)，把前边这一段根据 \(y=x+1\) 这条线做对称变换，起点就会变成 \((x-(x+1),x+1-x)=(-1,1)\)，终点还是 \((n,n)\)，路径数就是 \(C_{2n}^{n-1}\)。
选择题：关于 Catalan 数 \(C_n = \dfrac{(2n)!}{(n+1)!n!}\)，下列说法中错误的是？
- A. \(C_n\) 表示有 n+1 个节点的不同形态的二叉树的个数
- B. \(C_n\) 表示含 n 对括号的合法括号序列的个数
- C. \(C_n\) 表示长度为 n 的入栈序列对应的合法出栈序列个数
- D. \(C_n\) 表示通过连接顶点而将 n+2 边的凸多边形分成三角形的方法个数
答案
A。
应为 \(n\) 个节点的不同形态的二叉树个数。
习题：P1044 [NOIP 2003 普及组] 栈
解题思路
合法出栈序列数问题等价于卡特兰数。
参考代码
#include <cstdio>
int main()
{
int n; scanf("%d", &n);
int ans = 1;
for (int i = 1; i <= n; i++) ans = 1ll * ans * (4 * i - 2) / (i + 1);
printf("%d\n", ans);
return 0;
}
习题：P1754 球迷购票问题
解题思路
本质上是在任何位置前 \(50\) 元的人数不能少于 \(100\) 元的人数，卡特兰数问题。
参考代码
#include <cstdio>
using ll = long long;
int main()
{
int n; scanf("%d", &n);
ll ans = 1;
for (int i = 1; i <= n; i++) ans = ans * (4 * i - 2) / (i + 1);
printf("%lld\n", ans);
return 0;
}
习题：P1375 小猫
解题思路
枚举第一个点连第几个点，把问题拆成两部分，就可以得到 \(H_n = \sum \limits_{i=1}^n H_{i-1} \times H_{n-i}\) 这个递推式。
因此本题是卡特兰数问题，可以用卡特兰数的各种公式解决。
参考代码
#include <cstdio>
const int MOD = 1000000007;
const int N = 200005;
int f[N], g[N];
int qpow(int x, int y) {
int res = 1;
while (y > 0) {
if (y & 1) res = 1ll * res * x % MOD;
x = 1ll * x * x % MOD;
y >>= 1;
}
return res;
}
void init(int n) {
f[0] = 1; g[0] = 1;
for (int i = 1; i <= n; i++) f[i] = 1ll * f[i - 1] * i % MOD;
g[n] = qpow(f[n], MOD - 2);
for (int i = n - 1; i > 0; i--) g[i] = 1ll * g[i + 1] * (i + 1) % MOD;
}
int C(int n, int m) {
return 1ll * f[n] * g[m] % MOD * g[n - m] % MOD;
}
int main()
{
int n; scanf("%d", &n);
init(n * 2);
int ans = C(2 * n, n);
ans = (ans + MOD - C(2 * n, n - 1)) % MOD;
printf("%d\n", ans);
return 0;
}
浙公网安备 33010602011771号