二项式反演
基本式子
但平时做题一般都是用的另一种形态:
还有一种:
其中第一种是恰好和至多间的转换,第二种是恰好和至少之间的转换
证明的话可以看 这篇博客,讲的非常详细
例题
[HDU1465]不容易系列之一
Description
求所有 \(a_i \not= i\) 的长度为 \(n\) 的排列个数。
Sol
设 \(f_i\) 为错排数恰好为 \(i\) 的排列个数,那么我们的答案就是 \(f_n\),但 \(f\) 并不好求,我们考虑用 \(g_i\) 表示错排数至多为 \(i\) 的方案数,那么我们有:
很显然,我们固定 \(i\) 个位置为错排,有 \(C_n^i\) 中选法,每种选法的方案数为 \(f_i\)。
那么容易得知 \(g_n=n!\),那么我们就可以设 \(g_i=i!\),那么二项式反演一下,即可得:
Code
#include<bits/stdc++.h>
#define int long long
using namespace std;
int Read() {
int x = 0, f = 1; char ch = getchar();
while(!isdigit(ch)) {if(ch == '-') f = -1; ch = getchar();}
while(isdigit(ch)) {x = (x << 3) + (x << 1) + ch - '0'; ch = getchar();}
return x * f;
}
int f[25], fac[25], C[25][25];
signed main() {
fac[0] = 1;
for(int i = 1; i <= 20; i++) fac[i] = fac[i - 1] * i;
C[1][1] = C[1][0] = 1;
for(int i = 2; i <= 20; i++) {
C[i][0] = 1;
for(int j = 1; j <= i; j++) C[i][j] = C[i - 1][j] + C[i - 1][j - 1];
}
for(int i = 1; i <= 20; i++) {
for(int j = 0; j <= i; j++)
if((i - j) % 2 == 0) f[i] += C[i][j] * fac[j];
else f[i] -= C[i][j] * fac[j];
}
int n;
while(cin >> n) cout << f[n] << endl;
return 0;
}
[bzoj3622]已经没有什么好害怕的了
Description
给定 \(2\) 个长为 \(n\) 的数组 \(a,b\),将 \(a,b\) 中数两两配对,求 \(a_i>b_i\) 的组数比 \(a_i <b_i\) 的组数恰好多 \(k\) 的方案数。
\(1\le n\le 2000\)
Sol
首先将恰好转化为至少,设 \(f_j\) 为 \(a_i >b_i\) 的组数恰好为 \(j\)的方案数,\(g_j\) 表示至少为 \(j\) 的方案数。
则
对于这个 \(C_i^k\) 的意义,我们可以理解为:我们选择某 \(i\) 个位置时,它可以由 \(C_i^k\) 个只选 \(k\) 个数的方案往后选而得出,故此处的“至少”是有重复的。
那么我们考虑求 \(g_i\),我们将 \(a,b\) 排序,设 \(g_{i,j}\) 表示前 \(i\) 个数恰好有 \(j\) 个 \(a_k>b_k\),记 \(cnt_i\) 为前 \(i\) 个数中最大的 \(a>b\) 的个数,那么我们有方程式:
注意此时的 \(g\) 的 \(n-i\) 个组合可以任意分配,所以 \(g_{n,i}\) 还需乘上 \((n-i)!\),得到二项式反演式子中的 \(g\),代表至少 \(i\) 个(可重)的方案数。
求出 \(g\) 后二项式反演即可。
Code
#include<bits/stdc++.h>
#define int long long
#define Mod 1000000009
using namespace std;
int Read() {
int x = 0, f = 1; char ch = getchar();
while(!isdigit(ch)) {if(ch == '-') f = -1; ch = getchar();}
while(isdigit(ch)) {x = (x << 3) + (x << 1) + ch - '0'; ch = getchar();}
return x * f;
}
int n, k, fac[2505], C[2505][2505], a[2505], b[2505], cnt[2505], f[2505][2505];
signed main() {
n = Read(), k = Read();
if((n + k) % 2 != 0) {
cout << 0;
return 0;
}
k = (n + k) / 2;
for(int i = 1; i <= n; i++) a[i] = Read();
for(int i = 1; i <= n; i++) b[i] = Read();
sort(a + 1, a + n + 1); sort(b + 1, b + n + 1);
int nw = 0;
for(int i = 1; i <= n; i++) {
while(nw + 1 <= n && a[i] > b[nw + 1]) ++nw;
cnt[i] = nw;
}
for(int i = 0; i <= n; i++) f[i][0] = 1;
for(int i = 1; i <= n; i++) {
for(int j = 1; j <= cnt[i]; j++) {
if(j) f[i][j] += f[i - 1][j - 1] * (cnt[i] - j + 1) % Mod;
f[i][j] += f[i - 1][j]; f[i][j] %= Mod;
}
}
fac[0] = 1;
for(int i = 1; i <= 2500; i++) fac[i] = fac[i - 1] * i % Mod;
C[0][0] = 1;
for(int i = 1; i <= 2500; i++) {
C[i][0] = 1;
for(int j = 1; j <= i; j++) C[i][j] = (C[i - 1][j - 1] + C[i - 1][j]) % Mod;
}
for(int i = 0; i <= n; i++) (f[n][i] *= fac[n - i]) %= Mod;
int ans = 0;
for(int i = k; i <= n; i++) {
if((i - k) & 1) ans -= C[i][k] * f[n][i] % Mod;
else ans += C[i][k] * f[n][i] % Mod;
ans = (ans % Mod + Mod) % Mod;
}
cout << ans << endl;
return 0;
}
[HAOI2018]染色
Description
\(n\le 10^7,m\le10^5,S\le150\)
Sol
转化题意,设 \(f_k\) 为至少出现了 \(k\) 种颜色的次数为 \(S\) 的方案数,然后再二项式反演回去就是答案了。
考虑求出 \(f_k\),首先我们钦定 \(k\) 中颜色只能选 \(S\) 种,剩下的颜色可以随便选,这样就符合这里 “至少” 的定义了
首先选 \(k\) 种颜色的方法为 \(C_m^k\),\(ks\) 种颜色放在 \(n\) 个空位里的方案数为 \(C_{n}^{kS}\)。
在我们选的 \(kS\) 个空位里,我们有 \(kS!\) 种排列方法,但是相同颜色交换不会对答案产生贡献,一种颜色有 \(S!\) 种情况重复,\(k\) 种颜色就是 \((S!)^k\) 种情况重复,所以排列的方法数为 \(\frac{(kS)!}{(S!)^k}\)。
最后还剩下 \(m-k\) 种颜色可以随便放到 \(n-kS\) 个空位内,方案数为 \((m-k)^{n-kS}\)
乘起来就是 \(f_k\)。
求出 \(f_k\) 后我们就可以进行二项式反演了。
拆开组合数, $$ans_i \times i! = \sum_{j=i}^n (-1)^{j-i} \frac{1}{(j-i)!} f_j\times j!$$
设 \(g_i=f_i\times i!, h_i=\frac{(-1)^i}{i!}\),那么翻转 \(g_i\) 做个 NTT 就可以得到 \(ans_i\times i!\),除以 \(i!\) 就是我们的答案了。
Code
#include<bits/stdc++.h>
#define int long long
using namespace std;
const int Mod = 1004535809;
inline int Read() {
int x = 0, f = 1; char ch = getchar();
while(!isdigit(ch)) {if(ch == '-') f = -1; ch = getchar();}
while(isdigit(ch)) {x = (x << 3) + (x << 1) + ch - '0'; ch = getchar();}
return x * f;
}
int fac[10000005], finv[10000005], w[500005], g[500005], h[500005];
int C(int n, int m) {return fac[n] * finv[m] % Mod * finv[n - m] % Mod;}
int n, m, s, num;
int qpow(int a, int b) {
int res = 1;
while(b) {
if(b & 1) res = res * a % Mod;
a = a * a % Mod;
b >>= 1;
}
return res;
}
int l, lim, rev[1000005];
void NTT(int *a, int inv) {
for(int i = 0; i < lim; i++)
if(i < rev[i]) swap(a[i], a[rev[i]]);
for(int i = 1; i < lim; i <<= 1) {
int x = qpow(3, (Mod - 1) / (i << 1));
for(int j = 0; j < lim; j += (i << 1)) {
int y = 1;
for(int k = 0; k < i; k++, y = y * x % Mod) {
int p = a[j + k], q = a[i + j + k] * y % Mod;
a[j + k] = (p + q) % Mod, a[i + j + k] = (p - q + Mod) % Mod;
}
}
}
if(inv == -1) {
int ny = qpow(lim, Mod - 2);
reverse(a + 1, a + lim);
for(int i = 0; i < lim; i++) a[i] = a[i] * ny % Mod;
}
}
signed main() {
n = Read(), m = Read(), s = Read();
num = min(m, n / s);
for(int i = 0; i <= m; i++) w[i] = Read();
fac[0] = 1;
for(int i = 1; i <= 10000000; i++) fac[i] = fac[i - 1] * i % Mod;
finv[10000000] = qpow(fac[10000000], Mod - 2);
for(int i = 9999999; i >= 0; i--) finv[i] = finv[i + 1] * (i + 1) % Mod;
for(int i = 0; i <= num; i++)
g[i] = C(m, i) * C(n, i * s) % Mod * fac[i * s] % Mod * qpow(finv[s], i) % Mod * qpow(m - i, n - i * s) % Mod * fac[i] % Mod;
for(l = 0, lim = 1; lim <= (num << 1); lim <<= 1, ++l);
for(int i = 0; i < lim; i++) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (l - 1));
for(int i = 0; i <= num; i++) h[i] = (i & 1) ? Mod - finv[i] : finv[i];
reverse(g, g + num + 1);
NTT(g, 1); NTT(h, 1);
for(int i = 0; i < lim; i++) g[i] = g[i] * h[i] % Mod;
NTT(g, -1);
reverse(g, g + num + 1);
int ans = 0;
for(int i = 0; i <= num; i++) (ans += w[i] * g[i] % Mod * finv[i] % Mod) %= Mod;
cout << ans << endl;
return 0;
}