一道好题,第一次用生成函数做题。感谢赛珂狼教我这个做法。
首先我们显然可以把题目中的限制转化成一个二分图的模型:左边有$n$个点,右边有$m$个点,如果在棋盘$(i,j)$这个点上放了炮,那么我们把左边第$i$个点向右边第$j$个点连边,那么最终这个图上左边每个点的度数都是$2$,右边每个点的度数都小于等于$2$。求合法图的个数。
可以发现,这个二分图是由一些环和一些链组成的,每一条链都是属于右边的点数比属于左边的点数多一(我们把属于右边的单独的一个点也算作一条链)。并且这样链的条数恰好为$m-n$。
对于每条链或者每一个环,一旦分配给它的点固定后,其内部的方案数就和外部无关了,只和该链(环)的大小有关。
我们先计算有关于链的部分。
我们定义函数$g_n$表示左边有$n-1$个点,右边有$n$个点的一条链内部排列的方案数。那就有:
$$g = <0, 1, 1, 6, dots, frac{n!(n-1)!}{2}>$$
简述上述式子的由来,左边$n-1$个点有排列$(n-1)!$,右边$n$个点有排列$n!$,乘起来之后左右对称要除以$2$。我们构造一个生成函数把$g$挂上去,表示有一条链的情况:
$$G(x) = sum_{n} g_n frac{x^n}{n!(n-1)!}$$
如果有两条链我们就把它自己卷一下,并且用$g^{[2]}$来表示新得到的数列:
$$G^2(x) = sum_{n} left( sum_{i = 0}^{n} inom{n}{i} inom{n - 2}{i - 1} g_i g_{n - i}
ight) frac{x^n}{n!(n-2)!}$$
$$G^2(x) =sum_{n} g^{[2]}_n frac{x^n}{n!(n-2)!}$$
注意到这个生成函数的形式和$G(x)$相比有些许变动,即在$x^n$下的系数有改动,这是因为在两条链的时候右边的点数会比左边的点数多二,而不是一。于是更一般的,如果有$k$条链,就可以写成以下形式:
$$G^k(x) = sum_{n} g^{[k]}_n frac{x^n}{n!(n-k)!}$$
其中数列$g^{[m-n]}$就是我们想要的有恰好$m-n$条链的时候的方案数。但是细心的读者可以发现这里有一个小问题,就是我们在做卷积的过程中,我们每卷一次就会增加一条链,这相当于枚举了新加入的链的位置,于是在每一个合法的情况中,这$m-n$条链的加入顺序不同会被算成不同的方案,于是最终我们需要的应该是数列$frac{g^{[m-n]}}{(m-n)!}$。
理解这一点很重要,我们在计算环的时候同样需要用到这一点,请大家自行揣摩。
我们知道了$G(x)$,于是我们只要求出$G^{m-n}(x)$就行了,直接快速幂可以$O(mlog^2m)$,如果比较厉害写$O(mlogm)$的多项式快速幂也是可以的。
接下来我们算环的部分。
我们定义函数$f_n$表示左右边都有$n$个点的环内部排列的方案数,那就有:
$$f = <0, 0, 1, 6, dots ,frac{(n!)^2}{2n}>$$
简单说明上述式子的由来,左右边分别有圆排列$(n-1)!$,合并起来有$n$中方案,同时镜像对称要除以$2$。同时,这里有特殊情况$f_1=0$,因为你不能有一个两个点的环。这里有读者会发现,我们定义了$f_0 = 0$,这是因为我们接下来做卷积时必须保证每次选择环的大小时不能用一个空环(也就是会导致实际环的个数会少一)。
同样我们构造了一个生成函数,表示一个环的情况:
$$F(x) = sum_{n} f_n frac{x^n}{(n!)^2}$$
与链同理,如果我们用了$k$个环,就有:
$$F^k(x) = sum_{n} f^{[k]}_n frac{x^n}{(n!)^2}$$
与链同理,每一种用了$k$个环的合法情况会被重复计算到$k!$次,我们要除掉它,同时我们需要枚举我们用的环的个数$k$,于是我们需要的是:
$$H(x) = sum_{k} frac{F^k(x)}{k!} = e^{F(x)}$$
虽然我不会多项式$exp$,但是这里的形式比较简单我们可以手推。
我们需要先求出$F(x)$的一个比较简单的形式。根据$f$我们可以写出:
$$F(x) = -frac{x}{2} + frac{1}{2} sum_n frac{x^n}{n} = -frac{x}{2} - frac{1}{2} ln(1 - x)$$
于是有:
$$H(x) = sum_n h_n frac{x^n}{(n!)^2} = e^{F(x)} = e^{-frac{1}{2} x} e^{-frac{1}{2} ln(1-x)} = e^{-frac{1}{2}x} (1-x)^{-frac{1}{2}}$$
前半部分用$e$直接展开:
$$e^{-frac{1}{2}x} = sum_n (-frac{1}{2})^n frac{x^n}{n!}$$
后半部分用广义二项式展开:
$$(1-x)^{-frac{1}{2}} = sum_n frac{(-frac{1}{2})^{underline{n}}}{n!} (-1)^n x^n$$
把这两个卷起来就得到了$H(x)$。
最后我们需要把链和环拼起来,这一步做卷积或手动展开都可以了。
总复杂度是$O(mlogmlog(m-n))$,或者$O(mlogm)$。于是就和$m-n$无关了。
#include <cstdio> #include <vector> #include <algorithm> using namespace std; typedef vector<int> Poly; const int N = (int)1e5 + 5; const int MOD = 998244353; const int _IV2 = (MOD - 1) / 2; int n, m; namespace { int fac[N], ifac[N]; int Add(int a, int b) { return (a += b) >= MOD? a - MOD : a; } int Sub(int a, int b) { return (a -= b) < 0? a + MOD : a; } int Mul(int a, int b) { return (long long)a * b % MOD; } int Pow(int a, int b) { int r = 1; for (; b; b >>= 1, a = Mul(a, a)) if (b & 1) r = Mul(r, a); return r; } void Math_Init() { fac[0] = ifac[0] = 1; for (int i = 1; i < N; ++i) fac[i] = Mul(fac[i - 1], i); ifac[N - 1] = Pow(fac[N - 1], MOD - 2); for (int i = N - 2; ~i; --i) ifac[i] = Mul(ifac[i + 1], i + 1); } } namespace PO { int rev[(int)1e6 + 5]; void Dft(Poly &a, int le) { for (int i = 0; i < le; ++i) if (i < rev[i]) swap(a[i], a[rev[i]]); for (int i = 1; i < le; i <<= 1) { int wn = Pow(3, (MOD - 1) / (i << 1)); for (int j = 0; j < le; j += i << 1) { int w = 1, x, y; for (int k = 0; k < i; ++k, w = Mul(w, wn)) { x = a[j + k]; y = Mul(a[j + k + i], w); a[j + k] = Add(x, y); a[j + k + i] = Sub(x, y); } } } } void Idft(Poly &a, int le) { reverse(a.begin() + 1, a.end()); Dft(a, le); int iv = Pow(le, MOD - 2); for (int i = 0; i < le; ++i) a[i] = ::Mul(a[i], iv); } Poly Mul(Poly a, Poly b, int lim) { int n = a.size(), m = b.size(), l; for (l = 1; l < n + m - 1; l <<= 1); for (int i = 0; i < l; ++i) rev[i] = (rev[i >> 1] >> 1) | (i & 1? l >> 1 : 0); a.resize(l), Dft(a, l); b.resize(l), Dft(b, l); for (int i = 0; i < l; ++i) a[i] = ::Mul(a[i], b[i]); Idft(a, l), a.resize(lim); return a; } Poly Pow(Poly a, int b, int lim) { Poly r(lim); for (r[0] = 1; b; b >>= 1) { if (b & 1) r = Mul(r, a, lim); a = Mul(a, a, lim); } return r; } } int main() { Math_Init(); scanf("%d%d", &n, &m); Poly F0(m + 1), F1(m + 1), fm(m + 1); for (int i = 0, dw = 1; i <= m; ++i) { F0[i] = Mul(dw, ifac[i]); if (i & 1) F0[i] = Sub(0, F0[i]); dw = Mul(dw, Sub(_IV2, i)); } for (int i = 0, pw = 1; i <= m; ++i) { F1[i] = Mul(pw, ifac[i]); pw = Mul(pw, _IV2); } Poly expF = PO::Mul(F0, F1, m + 1); for (int i = 0; i <= m; ++i) { fm[i] = Mul(expF[i], Mul(fac[i], fac[i])); } Poly G(m + 1), gk(m + 1); G[1] = 1; for (int i = 2; i <= m; ++i) { G[i] = (MOD + 1) / 2; } G = PO::Pow(G, m - n, m + 1); for (int i = m - n; i <= m; ++i) { gk[i] = Mul(Mul(G[i], Mul(fac[i], fac[i - m + n])), ifac[m - n]); } int ans = 0; for (int i = 0; i <= n; ++i) { int ch1 = Mul(fac[m], Mul(ifac[i], ifac[m - i])); int ch2 = Mul(fac[n], Mul(ifac[i], ifac[n - i])); ans = Add(ans, Mul(Mul(fm[i], gk[m - i]), Mul(ch1, ch2))); } printf("%d ", ans); return 0; }
这里简单讲一下为什么有$sum_{n=1} frac{x^n}{n} = -ln(1-x)$。
根据定义有:$ln'(x) = frac{1}{x}$
根据导数的复合:$(ln(1-x))' = ln'(1-x)(1-x)' = frac{-1}{1-x} = -sum_n x^n$
对两边积分:$ln(1-x) = -sum_{n=1} frac{x^n}{n}$