大致题意: 计算(sum_{i=0}^nsum_{j=0}^iS(i,j)*2^j*(j!)),其中(S)为第二类斯特林数。
推式子
首先让我们来推一波式子:
因为当(i<j)时,(S(i,j)=0),所以,为了方便式子的化简,我们可以先将第二个(sum)的上限全部改成(n),即:
[sum_{i=0}^nsum_{j=0}^nS(i,j)*2^j*(j!)
]
这样一来,(sum_{j=0}^n2^j*(j!))这个式子就与(i)无关,则我们可以将其提前:
[sum_{j=0}^n2^j*(j!)sum_{i=0}^nS(i,j)
]
套用第二类斯特林数的通项公式(S(n,m)=sum_{i=0}^mfrac{(-1)^i*(m-i)^n}{i!(m-i)!})可以得到:
[sum_{j=0}^n2^j*(j!)sum_{i=0}^nsum_{k=0}^jfrac{(-1)^k*(j-k)^i}{k!(j-k)!}
]
接下来,我们可以考虑再将所有与(i)无关的项提前,得到下面这个式子:
[sum_{j=0}^n2^j*(j!)sum_{k=0}^jfrac{(-1)^k}{k!(j-k)!}*sum_{i=0}^n(j-k)^i
]
我们可以枚举(j),这样前面的(sum_{j=0}^n2^j*(j!))一项就可以轻松搞定。
那么接下来的问题就是如何对于给定的(j),求出下面这个式子的值:
[sum_{k=0}^jfrac{(-1)^k}{k!(j-k)!}*sum_{i=0}^n(j-k)^i
]
仔细观察,可以发现似乎这个式子前半部分的分母和后半部分都有((j-k))这个因式,则容易想到将它们移到一起,即:
[sum_{k=0}^jfrac{(-1)^k}{k!}*sum_{i=0}^nfrac{(j-k)^i}{(j-k)!}
]
这样一移的效果是显著的,因为此时这个式子的前半部分只与(k)有关,后半部分只与((j-k))有关。
则可以设(f)和(g)如下:
[f_x=frac{(-1)^x}{x!}
]
[g_x=sum_{i=0}^nfrac{x^i}{x!}=frac{sum_{i=0}^nx^i}{x!}=frac{frac{x^{n+1}-1}{x-1}}{x!}=frac{x^{n+1}-1}{(x-1)x!}
]
带入原式就可以得到:
[sum_{k=0}^jf_{k}*g_{j-k}
]
而这个式子实际上就相当于:
[(f*g)(j)
]
至此化简完毕,最终的式子就是:
[sum_{j=0}^n2^j*(j!)*(f*g)(j)
]
具体实现
我们可以先(O(n))((O(nlogn))?)预处理出(f)和(g)两个数组,然后(NTT)即可。
这应该还是比较简单的吧。
代码
#include<bits/stdc++.h>
#define Tp template<typename Ty>
#define Ts template<typename Ty,typename... Ar>
#define Reg register
#define RI Reg int
#define Con const
#define CI Con int&
#define I inline
#define W while
#define N 100000
#define X 998244353
#define Qinv(x) (Qpow(x,X-2))
#define swap(x,y) (x^=y^=x^=y)
#define Inc(x,y) ((x+=(y))>=X&&(x-=X))
#define XSum(x,y) ((x)+(y)>=X?(x)+(y)-X:(x)+(y))
#define XSub(x,y) ((x)-(y)<0?(x)-(y)+X:(x)-(y))
using namespace std;
int n,a[(N<<2)+5],b[(N<<2)+5],Fac[N+5],Inv[N+5];
I int Qpow(RI x,RI y) {RI res=1;W(y) y&1&&(res=1LL*res*x%X),x=1LL*x*x%X,y>>=1;return res;}
class NTT//NTT求模意义下的卷积
{
private:
static const int SZ=N,PR=3,IP=(X+1)/3;int P,L,R[(SZ<<2)+5];
I void Transform(int* s,CI op)
{
RI i,j,k,U,S,tx,ty;for(i=0;i^P;++i) i<R[i]&&swap(s[i],s[R[i]]);
for(i=1;i^P;i<<=1) for(U=Qpow(~op?PR:IP,(X-1)/(i<<1)),j=0;j^P;j+=(i<<1))
for(S=1,k=0;k^i;++k,S=1LL*S*U%X) tx=s[j+k],ty=1LL*S*s[i+j+k]%X,s[j+k]=XSum(tx,ty),s[i+j+k]=XSub(tx,ty);
}
public:
I void Solve(CI n,CI m,int* a,int* b)
{
RI i,t;P=1,L=0,memset(R,0,sizeof(R));
W(P<=n+m) P<<=1,++L;for(i=0;i^P;++i) R[i]=(R[i>>1]>>1)|((i&1)<<L-1);
for(Transform(a,1),Transform(b,1),i=0;i^P;++i) a[i]=1LL*a[i]*b[i]%X;
for(Transform(a,-1),t=Qinv(P),i=0;i<=n+m;++i) a[i]=1LL*a[i]*t%X;
}
}NTT;
int main()
{
RI i,t,ans=0;for(scanf("%d",&n),Fac[0]=1,i=1;i<=n;++i) Fac[i]=1LL*Fac[i-1]*i%X;//预处理阶乘
for(Inv[n]=Qpow(Fac[n],X-2),i=n-1;~i;--i) Inv[i]=1LL*Inv[i+1]*(i+1)%X;//预处理阶乘的逆元
for(t=1,i=0;i<=n;++i,t=1LL*t*(X-1)%X) a[i]=1LL*t*Inv[i]%X;//预处理第一个数组
for(b[0]=1,b[1]=n+1,i=2;i<=n;++i) b[i]=1LL*(Qpow(i,n+1)-1)*Qinv(i-1)%X*Inv[i]%X;//预处理第二个数组
for(NTT.Solve(n,n,a,b),t=1,i=0;i<=n;++i,(t<<=1)>=X&&(t-=X)) Inc(ans,1LL*a[i]*t%X*Fac[i]%X);//做一遍NTT,然后统计答案
return printf("%d",ans),0;//输出答案
}