「TJOI/HEOI2016」求和 - 指数型生成函数 + 多项式求逆 / 分治FFT | Bill Yang's Blog

「TJOI/HEOI2016」求和 - 指数型生成函数 + 多项式求逆 / 分治FFT

题目大意

    在2016年,佳媛姐姐刚刚学习了第二类斯特林数,非常开心。

    现在他想计算这样一个函数的值:

    $S(i,j)$表示第二类斯特林数,递推公式为:$S(i,j)=j\cdot S(i-1,j)+S(i-1,j-1),1\le j\le i-1$。

    边界条件为:$S(i, i) = 1(0 \leq i), \ S(i, 0) = 0(1 \leq i)$

    你能帮帮她吗?


题目分析

关于分治FFT直接看Candy?dalao的题解吧。
下面将多项式求逆的方法。

令$g(n)=\sum\limits_{i=1}^n {n \brace i}\cdot2^i\cdot i!$,那么答案就是其前缀和。
考虑$g(n)$的组合意义,代表将$n$个$1−n$的数分成若干集合,给每个集合黑白染色的方案数。

拆分组合数得到:

又出现了,指数型生成函数!
令:

那么有:

这里的$+1$是因为$G(x)$与$F(x)$卷积后第$0$次项(F从$1$开始卷积)没了,然而$G(x)$是有的,因此需要补$g(0)$回去,而根据组合意义(并非表达式)可以知道$g(0)=1$。

因此:

直接上多项式求逆即可。


代码

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
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
#include<bits/stdc++.h>
using namespace std;
typedef long long LL;
inline int Get_Int() {
int num=0,bj=1;
char x=getchar();
while(!isdigit(x)) {if(x=='-')bj=-1;x=getchar();}
while(isdigit(x)) {num=num*10+x-'0';x=getchar();}
return num*bj;
}
const int maxn=262144+5;
const LL mod=998244353,g=3;
LL Quick_Pow(LL a,LL b) {
LL sum=1;
for(; b; b>>=1,a=a*a%mod)if(b&1)sum=sum*a%mod;
return sum;
}
LL inv(LL x) {return Quick_Pow(x,mod-2);}
struct NumberTheoreticTransform {
int n,rev[maxn];
LL omega[maxn],iomega[maxn];
void init(int n) {
this->n=n;
int x=Quick_Pow(g,(mod-1)/n);
omega[0]=iomega[0]=1;
for(int i=1; i<n; i++)omega[i]=omega[i-1]*x%mod,iomega[i]=inv(omega[i]);
int k=log2(n);
for(int i=0; i<n; i++) {
int t=0;
for(int j=0; j<k; j++)if(i&(1<<j))t|=1<<(k-j-1);
rev[i]=t;
}
}
void transform(LL *a,LL *omega) {
for(int i=0; i<n; i++)if(i<rev[i])swap(a[i],a[rev[i]]);
for(int len=2; len<=n; len*=2) {
int mid=len>>1;
for(LL *p=a; p!=a+n; p+=len)
for(int i=0; i<mid; i++) {
LL t=omega[n/len*i]*p[mid+i]%mod;
p[mid+i]=(p[i]-t+mod)%mod;
p[i]=(p[i]+t)%mod;
}
}
}
void dft(LL *a) {transform(a,omega);}
void idft(LL *a) {
transform(a,iomega);
LL x=inv(n);
for(int i=0; i<n; i++)a[i]=a[i]*x%mod;
}
} ntt;
void polynomial_inverse(const LL* a,const int n,LL* b) {
if(n==1) {b[0]=inv(a[0]);return;}
polynomial_inverse(a,n>>1,b);
int p=n<<1;
static LL x[maxn];
copy(a,a+n,x),fill(x+n,x+p,0);
ntt.init(p),ntt.dft(x),ntt.dft(b);
for(int i=0; i<p; i++)b[i]=b[i]*((2-x[i]*b[i]%mod+mod)%mod)%mod;
ntt.idft(b),fill(b+n,b+p,0);
}
int n;
LL fac[maxn],invf[maxn],F[maxn],G[maxn];
int main() {
n=Get_Int();
int T=1;
while(T<=n)T<<=1;
fac[0]=1;
for(int i=1; i<=T; i++)fac[i]=fac[i-1]*i%mod;
invf[T]=inv(fac[T]);
for(int i=T; i>=1; i--)invf[i-1]=invf[i]*i%mod;
for(int i=1; i<=n; i++)F[i]=(mod-invf[i])*2%mod;
F[0]=1;
polynomial_inverse(F,T,G);
LL ans=0;
for(int i=n; i>=0; i--)ans=(ans+G[i]*fac[i]%mod)%mod;
printf("%lld\n",ans);
return 0;
}
本篇文章有用吗?有用还不快点击!
0%