LOJ #2058「TJOI / HEOI2016」求和

不错的推柿子题

LOJ #2058

题意：求$sumlimits_{i=0}^nsumlimits_{j=0}^nS(i,j)·2^j·j!$其中$ S(n,m)$是第二类斯特林数

---

$ Solution:$

首先考虑第二类斯特林数的意义：将$ n$个有标号元素放入$ m$个无标号集合（无空集）的方案数

我们枚举空集的数量容斥：$ S(n,m)=frac{1}{m!}sumlimits_{k=0}^m(-1)^kC_m^k(m-k)^n$

乘上$ frac{1}{m!}$是因为容斥的集合带标号而斯特林数本身不带标号

这样可以将原式展开得：

$ sumlimits_{i=0}^n sumlimits_{j=0}^n2^j sumlimits_{k=0}^j(-1)^kC_j^k(j-k)^i$     (消阶乘项)

把组合数展开得$ sumlimits_{i=0}^n sumlimits_{j=0}^n 2^j j! sumlimits_{k=0}^j frac{(-1)^k}{k!} frac{(j-k)^i}{(j-k)!}$

改变枚举顺序得$ sumlimits_{j=0}^n 2^j j! sumlimits_{k=0}^j frac{(-1)^k}{k!} frac{ sumlimits_{i=0}^n (j-k)^i}{(j-k)!}$

令$ A(x)= frac{(-1)^x}{x!}$，$ B(x)=frac{ sumlimits_{i=0}^n x^i}{x!}$

则原式为$ sumlimits_{j=0}^n 2^j j! sumlimits_{k=0}^jA(k)B(j-k)$

容易发现这是一个卷积形式，而函数$ A,B$均可以在$ O(n)$时间复杂度内完成

这样可以直接用$ NTT$优化，时间复杂度:$ O(n log   n)$

---

$ my code:$

#include<ctime>
#include<cmath>
#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<queue>
#include<vector>
#define p 998244353
#define rt register int
#define ll long long
using namespace std;
inline ll read(){
    ll x = 0; char zf = 1; char ch = getchar();
    while (ch != '-' && !isdigit(ch)) ch = getchar();
    if (ch == '-') zf = -1, ch = getchar();
    while (isdigit(ch)) x = x * 10 + ch - '0', ch = getchar(); return x * zf;
}
void write(ll y){if(y<0)putchar('-'),y=-y;if(y>9)write(y/10);putchar(y%10+48);}
void writeln(const ll y){write(y);putchar('
');}
int i,j,k,m,n,x,y,z,cnt;
int inv[100010],jc[100010],njc[100010];
int ksm(int x,int y){
    int ans=1;
    for(rt i=y;i;i>>=1,x=1ll*x*x%p)if(i&1)ans=1ll*ans*x%p;
    return ans;
}
vector<int>A,B,f,R;int lim=1;
int calc(int x,int L,int R){
    if(x==1)return R-L+1;
    return 1ll*(ksm(x,R+1)-ksm(x,L))*ksm(x-1,p-2)%p;
}
void init(int n){
    for(rt i=0;i<=1;i++)inv[i]=jc[i]=njc[i]=1;
    for(rt i=2;i<=n;i++){
        jc[i]=1ll*jc[i-1]*i%p;
        inv[i]=1ll*inv[p%i]*(p-p/i)%p;
        njc[i]=1ll*njc[i-1]*inv[i]%p;
    }    
    while(lim<=n+n)lim<<=1;
    A.resize(lim);B.resize(lim);f.resize(lim);
    A[0]=1;for(rt i=1,tag=-1;i<=n;i++,tag*=-1)A[i]=tag*njc[i];    
    B[0]=1;for(rt i=1;i<=n;i++)B[i]=1ll*njc[i]*calc(i,0,n)%p;
}
namespace poly{
    void getR(int n){
        R.resize(n);
        for(rt i=1;i<n;i++)R[i]=(R[i>>1]>>1)|(i&1)*(n>>1);
    }
    void NTT(int n,vector<int>&A,int fla){
        for(rt i=0;i<n;i++)if(i>R[i])swap(A[i],A[R[i]]);
        for(rt i=1;i<n;i<<=1){
            int w=ksm(3,(p-1)/2/i);
            for(rt j=0;j<n;j+=i<<1){
                int K=1;
                for(rt k=0;k<i;k++,K=1ll*K*w%p){
                    int x=A[j+k],y=1ll*K*A[i+j+k]%p;
                    A[j+k]=(x+y)%p,A[i+j+k]=(x-y)%p;
                }
            }
        }
        if(fla==-1){
            reverse(A.begin()+1,A.end());int invn=ksm(n,p-2);
            for(rt i=0;i<n;i++)A[i]=1ll*A[i]*invn%p;
        }
    }
}
using namespace poly;
int main(){
    n=read();init(n);
    int ans=0;getR(lim);
    NTT(lim,A,1);NTT(lim,B,1);
    for(rt i=0;i<lim;i++)f[i]=1ll*A[i]*B[i]%p;
    NTT(lim,f,-1);
    for(rt i=0;i<=n;i++)(ans+=1ll*ksm(2,i)*jc[i]%p*f[i]%p)%=p;
    cout<<(ans+p)%p;
    return 0;
}