前几天在BZ上的考试考到有关第二类斯特林数的东西
虽然说那道题目到最后并不需要用这个东西来化简把
不过抱着学习的态度还是学了学有关第二类斯特林数的东西
第二类斯特林数S(n,m)定义为把n个元素划分成m个无序集合的方案数
根据这个定义我们不难写出递推式
设状态S(i,j)，讨论第i个元素是否单独一个集合
若单独一个集合，则方案数等价于S(i-1,j-1)
若不是单独一个集合，则他可以在之前任意j个集合里，方案为S(i-1,j)*j
这样我们得到式子S(i,j)=S(i-1,j-1)+S(i-1,j)*j
递推的时间复杂度是O(n^2)的
根据定义我们也不难写出通项表达式
假设集合没有非空的限制，则答案显然是m^n
之后我们可以利用容斥原理，枚举至少有几个集合是空的，在套上容斥系数我们有(注意集合是无序的，所以最后要消序)
S(n,m)=1/m!*sigma((-1)^k*C(m,k)*(m-k)^n))
直接根据这个式子暴力求解某一行的斯特林数的复杂度是O(n^2)的
但显然这是一个卷积形式，利用FFT我们可以将时间复杂度优化成O(nlogn)
我们再来讨论一下贝尔数，贝尔数的定义是把n个元素划分成若干个无序集合的方案数
不难发现f(n)=sigma(S(n,m)),直接求解某一项贝尔数利用上面的通项公式我们可以做到O(nlogn)
但是我们可以做到更好，考虑第一个元素所在集合的大小，我们不难得到如下递推式
f(i)=sigma(C(i-1,k-1)*f(i-k))
我们对于这个式子利用CDQ+FFT即可在O(nlog^2n)的时间复杂度内求解前n项的贝尔数
但是我们还可以做到更好
我们考虑每个集合是由若干个元素组成的且集合非空，集合中元素无序
则可以得到集合的生成函数为g(x)=e^x-1
而贝尔数分解成若干个集合的划分方案,则我们有f=e^g
我们构造多项式g，之后多项式求exp即可
时间复杂度O(nlogn)
有意思的是我们可以考虑把n个元素划分成若干个有序集合的方案数，设为多项式f
考虑第一个集合的大小，我们有
f(i)=sigma(C(i,k)*f(i-k))
我们一样可以用CDQ+FFT来求解
但是我们利用之前在多项式求逆学习总结中提过的技巧进行化简，就可以化简成一个多项式求逆的式子
时间复杂度O(nlogn)
我们现在讨论一些跟斯特林数相关的化简技巧
其中最常用的莫过于i^k和S(k,j)之间的转化了
我们通过第二类斯特林数的通项公式不难发现S(k,j)可以转化成若干个带系数的i^k的和
而我们又存在公式i^k=sigma(S(k,j)*j!*C(i,j))
不难发现这里跟i有关的项其实只有C(i,j),我们就把幂通过第二类斯特林数转化成了组合数
更有趣的是，这个式子后面j!*C(i,j)实际上暴力计算的时间复杂度是O(j)的 QAQ
即计算i的下降阶乘幂
我们来看几道题目QAQ
BZ的题目由于某些奇怪的原因貌似不能放出来QAQ
首先是Crash的文明世界 cojs 1888
这是道很有趣的题目，首先我们如果对(i+1)^k暴力用二项式定理展开会得到一个时间复杂度为O(n*k^2)的DP方程
注意到这里状态是O(n*k)的，而转移是O(k)的
由于k很小，我们考虑运用第二类斯特林数将i^k转化为C(i,j)
由于C(i,j)=C(i-1,j)+C(i-1,j-1)，我们利用这个式子可以做到转移O(1)，这样时间复杂度就是O(n*k)的
具体做法是设f(i,j)表示第i个点的C(dis(i,o),j)的和
然后从底向上做一遍DP，然后从上向下在做一遍DP
求出第二类斯特林数还原答案即可
#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<cstdlib>
using namespace std;
typedef long long LL;
const int maxn=50010;
const int mod=10007;
int n,k,L,u,v;
int now,A,B,Q;
int h[maxn],cnt=0;
struct edge{
int to,next;
}G[maxn<<1];
int f[maxn][151];
int jc[151],inv[151];
int S[151][151];
void add(int x,int y){
++cnt;G[cnt].to=y;G[cnt].next=h[x];h[x]=cnt;
}
void Get_read(){
scanf("%d%d%d",&n,&k,&L);
scanf("%d%d%d%d",&now,&A,&B,&Q);
for(int i=1;i<n;++i){
now=(now*A+B)%Q;
int tmp=((i<L)?i:L);
u=i-now%tmp;v=i+1;
add(u,v);add(v,u);
}
}
int pow_mod(int v,int p){
int tmp=1;
while(p){
if(p&1)tmp=tmp*v%mod;
v=v*v%mod;p>>=1;
}return tmp;
}
void Get_S(){
S[0][0]=1;
for(int i=1;i<=k;++i){
S[i][0]=0;
for(int j=1;j<=i;++j)S[i][j]=(S[i-1][j-1]+S[i-1][j]*j)%mod;
}
jc[0]=1;
for(int i=1;i<=k;++i)jc[i]=jc[i-1]*i%mod;
inv[k]=pow_mod(jc[k],mod-2);
for(int i=k-1;i>=0;--i)inv[i]=inv[i+1]*(i+1)%mod;
}
void Get_DFS(int u,int fa){
f[u][0]=1;
for(int i=h[u];i;i=G[i].next){
int v=G[i].to;
if(v==fa)continue;
Get_DFS(v,u);
f[u][0]=(f[u][0]+f[v][0])%mod;
for(int j=1;j<=k;++j)f[u][j]=(f[u][j]+f[v][j]+f[v][j-1])%mod;
}return;
}
void Get_DP(int u,int fa){
for(int i=h[u];i;i=G[i].next){
int v=G[i].to;
if(v==fa)continue;
for(int j=k;j>=2;--j){
f[v][j]=(f[u][j]+f[u][j-1]-2*f[v][j-1]-f[v][j-2])%mod;
if(f[v][j]<0)f[v][j]+=mod;
}
f[v][1]=(f[u][0]+f[u][1]-2*f[v][0])%mod;
if(f[v][1]<0)f[v][1]+=mod;
f[v][0]=f[u][0];
Get_DP(v,u);
}return;
}
int main(){
freopen("civilization.in","r",stdin);
freopen("civilization.out","w",stdout);
Get_read();Get_S();
Get_DFS(1,-1);
Get_DP(1,-1);
for(int i=1;i<=n;++i){
int ans=0;
for(int j=0;j<=k;++j){
ans=ans+S[k][j]*jc[j]%mod*f[i][j]%mod;
ans%=mod;
}printf("%d
",ans);
}return 0;
}
cojs 2359 QAQ的图论题
QAQ 题解在blog里，写的非常详细 QAQ
具体做法是考虑每个点的贡献写出O(n)的计算式，之后把i^k用第二类斯特林数表示
化简之后计算的时间复杂度变成O(k)，然后利用FFT在O(klogk)的时间内求出我们所需要的斯特林数就可以了
#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<iostream>
#include<algorithm>
#define G 3
using namespace std;
typedef long long LL;
const int mod=998244353;
const int maxn=300010;
int n,k,N,len,ans,now;
int jc[maxn],inv[maxn],rev[maxn];
int A[maxn],B[maxn],w[maxn];
int pow_mod(int v,int p){
int tmp=1;
while(p){
if(p&1)tmp=1LL*tmp*v%mod;
v=1LL*v*v%mod;p>>=1;
}return tmp;
}
void FFT(int *A,int n,int flag){
for(int i=0;i<n;++i)if(i<rev[i])swap(A[i],A[rev[i]]);
for(int k=2;k<=n;k<<=1){
int mk=(k>>1);
int wn=pow_mod(G,flag==1?(mod-1)/k:(mod-1)-(mod-1)/k);
w[0]=1;
for(int i=1;i<mk;++i)w[i]=1LL*w[i-1]*wn%mod;
for(int i=0;i<n;i+=k){
for(int j=0;j<mk;++j){
int a=i+j,b=i+j+mk;
int x=A[a],y=1LL*A[b]*w[j]%mod;
A[a]=(x+y)%mod;
A[b]=x-y;if(A[b]<0)A[b]+=mod;
}
}
}
if(flag==-1){
int c=pow_mod(n,mod-2);
for(int i=0;i<n;++i)A[i]=1LL*A[i]*c%mod;
}return;
}
int main(){
freopen("QAQ_Graph.in","r",stdin);
freopen("QAQ_Graph.out","w",stdout);
scanf("%d%d",&n,&k);
jc[0]=1;
for(int i=1;i<=k;++i)jc[i]=1LL*jc[i-1]*i%mod;
inv[k]=pow_mod(jc[k],mod-2);
for(int i=k-1;i>=0;--i)inv[i]=1LL*inv[i+1]*(i+1)%mod;
for(int i=0;i<=k;++i){
if(i&1)A[i]=mod-inv[i];
else A[i]=inv[i];
}
for(int i=0;i<=k;++i)B[i]=1LL*pow_mod(i,k)*inv[i]%mod;
for(N=1,len=0;N<=k;N<<=1,len++);N<<=1,len++;
for(int i=0;i<N;++i)rev[i]=rev[i>>1]>>1|((i&1)<<(len-1));
FFT(A,N,1);FFT(B,N,1);
for(int i=0;i<N;++i)A[i]=1LL*A[i]*B[i]%mod;
FFT(A,N,-1);now=1;ans=0;
for(int i=0;i<=k;++i){
if(now==0)break;
ans=ans+1LL*A[i]*now%mod*pow_mod(2,n-i-1)%mod;
ans%=mod;now=1LL*now*(n-i-1)%mod;
}
ans=1LL*ans*n%mod;
ans=1LL*ans*pow_mod(2,(1LL*(n-1)*(n-2)/2)%(mod-1))%mod;
printf("%d
",ans);
return 0;
}
cojs 2272 HEOI2016 sum
QAQ 考试的时候没有做出来，所以要反反复复的多虐几遍 QAQ
第一发的做法是CDQ+FFT,第二发的做法是多项式求逆
而现在会了第二类斯特林数的一些性质，我们有了第三种做法
前两种做法都是设f(i)=sigma(S(i,j)*2^j*j!)
但是我们不妨转换思路,设f(i)=sigma(S(j,i))
这样我们求完之后对每一项乘以2^i*i!就可以了
考虑把S(j,i)用通项公式展开
对于S(n,m)我们有S(n,m)=sigma((-1)^k*C(m,k)*(m-k)^n)
不难发现跟n有关的项只有(m-k)^n，提出来求sigma可以用等比数列求和公式化简掉
之后带入原式我们发现这是一个裸的卷积形式，直接做FFT就可以了
#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<iostream>
#include<algorithm>
#define G 3
using namespace std;
typedef long long LL;
const int maxn=300010;
const int mod=998244353;
int n,N,len,ans,now;
int jc[maxn],inv[maxn];
int A[maxn],B[maxn],rev[maxn];
int w[maxn];
int pow_mod(int v,int p){
int tmp=1;
while(p){
if(p&1)tmp=1LL*tmp*v%mod;
v=1LL*v*v%mod;p>>=1;
}return tmp;
}
void FFT(int *A,int n,int flag){
for(int i=0;i<n;++i)if(i<rev[i])swap(A[i],A[rev[i]]);
for(int k=2;k<=n;k<<=1){
int mk=(k>>1);
int wn=pow_mod(G,flag==1?(mod-1)/k:(mod-1)-(mod-1)/k);
w[0]=1;
for(int i=1;i<mk;++i)w[i]=1LL*w[i-1]*wn%mod;
for(int i=0;i<n;i+=k){
for(int j=0;j<mk;++j){
int a=i+j,b=i+j+mk;
int x=A[a],y=1LL*A[b]*w[j]%mod;
A[a]=(x+y)%mod;
A[b]=x-y;if(A[b]<0)A[b]+=mod;
}
}
}
if(flag==-1){
int c=pow_mod(n,mod-2);
for(int i=0;i<n;++i)A[i]=1LL*A[i]*c%mod;
}return;
}
int main(){
freopen("heoi2016_sum.in","r",stdin);
freopen("heoi2016_sum.out","w",stdout);
scanf("%d",&n);
jc[0]=1;
for(int i=1;i<=n;++i)jc[i]=1LL*jc[i-1]*i%mod;
inv[n]=pow_mod(jc[n],mod-2);
for(int i=n-1;i>=0;--i)inv[i]=1LL*inv[i+1]*(i+1)%mod;
for(int i=0;i<=n;++i){
if(i&1)A[i]=mod-inv[i];
else A[i]=inv[i];
}B[1]=n;
for(int i=2;i<=n;++i){
B[i]=pow_mod(i,n+1)-i;
if(B[i]<0)B[i]+=mod;
B[i]=1LL*B[i]*pow_mod(i-1,mod-2)%mod;
B[i]=1LL*B[i]*inv[i]%mod;
}
for(N=1,len=0;N<=n;N<<=1,len++);N<<=1,len++;
for(int i=0;i<N;++i)rev[i]=rev[i>>1]>>1|((i&1)<<(len-1));
FFT(A,N,1);FFT(B,N,1);
for(int i=0;i<N;++i)A[i]=1LL*A[i]*B[i]%mod;
FFT(A,N,-1);ans=1;now=2;
for(int i=1;i<=n;++i){
ans=ans+1LL*A[i]*now%mod*jc[i]%mod;
ans%=mod;now=(now<<1)%mod;
}printf("%d
",ans);
return 0;
}
QAQ 第二类斯特林数还是很有趣的呢
准备在cojs补上一道求贝尔数的题目QAQ