类似中国象棋，对于一个这样的二分图我们往往考察其每个连通块的构造。
分别用 \(x,y\) 计量左右部点的个数，两维均是 EGF。
一个连通块有三种情况：环、链、点，接下来我们分别来考虑它们。
对于环，它显然是偶环，那么左右部点的个数是一样的，并且正反各会被统计一次。而一边是排列，另一边是环排列，因此就是 \[ \frac 12\left(-\ln(1-xy)-xy\right) \]
对于链，右部点一定比左部点多一个，并且正反各会被统计一次。且两边都是排列，因此 \[ \frac 12 \frac{xy^2}{1-xy} \]
对于点，必然是一个右部点，因此 \[ y \]
现在将它们组合，那么答案就是 \[ \begin{aligned} n!m! [x^n y^m] F(x,y) &= n!m! [x^n y^m] \exp\left(\frac 12\left(-\ln(1-xy)-xy\right)+\frac 12 \frac{xy^2}{1-xy}+y\right) \\ &= n!m! [s^n t^{m-n}] \exp\left(\frac 12\left(-\ln(1-s)-s\right)+\frac 12 \frac{st}{1-s}+t\right) \\ &= n!m! [s^n t^{m-n}] \exp\left(\frac 12\left(-\ln(1-s)-s\right)\right) {\rm e}^t \sum\limits_{k\ge 0} \frac 1{k!} \left(\frac 12 \frac{st}{1-s}\right)^k \\ &= n!m^{\underline n} [s^n] \exp\left(\frac 12\left(-\ln(1-s)-s\right)\right) \sum\limits_{k\ge 0} \binom{m-n}k \left(\frac 12 \frac s{1-s}\right)^k \\ &= n!m^{\underline n} [s^n] \exp\left(\frac 12\left(-\ln(1-s)-s\right)\right) \left(1+\frac 12 \frac s{1-s}\right)^{m-n} \\ &= n!m^{\underline n} [s^n] \frac1{ \sqrt{1-s} } {\rm e}^{-s/2} \left(\frac{2-s}{2-2s}\right)^{m-n} \end{aligned} \]
而其三者均微分有限，故其乘积也微分有限，借助计算机便可得到递推式 \[ f_n = \frac1{2n}\left((m-2n-3) f_{n-1} + (3-n) f_{n-2} - \frac 12 f_{n-3}\right) \]
代码： 1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
using namespace std;
const int N = 2e6;
const int mod = 943718401;
inline int fpow(int a,int b)
{
int ret = 1;
for(;b;b >>= 1)
(b & 1) && (ret = (long long)ret * a % mod),a = (long long)a * a % mod;
return ret;
}
int n,k;
long long m;
int f[N + 5],ans;
int main()
{
scanf("%d%lld",&n,&m),k = (m - n) % mod;
f[0] = 1,
f[1] = 471859201LL * k % mod,
f[2] = (825753601LL * k % mod * k + 589824001LL * k + 707788801) % mod * 2 % mod,
f[3] = (924057601LL * k % mod * k % mod * k + 766771201LL * k % mod * k + 550502401LL * k + 786432001) % mod * 6 % mod;
for(register int i = 4;i <= n;++i)
f[i] = (k + 3LL * (i - 1)) * f[i - 1] % mod,
f[i] = (f[i] + (3LL - i + mod) * (i - 1) % mod * f[i - 2]) % mod,
f[i] = (f[i] + 471859200LL * (i - 1) % mod * (i - 2) % mod * f[i - 3]) % mod,
f[i] = 471859201LL * f[i] % mod;
ans = 1;
for(register int i = 1;i <= n;++i)
ans = (m - i + 1) % mod * ans % mod;
ans = (long long)ans * f[n] % mod;
printf("%d\n",ans);
}