斯特林数的四种求法

· · 个人记录

斯特林数是个比较复杂的东西，数学基础不吼的同学慎入。

斯特林数有四种O(nlogn)求法,你知道么?

我们知道 m^n=\sum\limits_{i=0}^m\begin{Bmatrix}n \\ i\end{Bmatrix}*i!*\dbinom{m}{i} 这个式子是能够二项式反演的。

而且众所周知，二项式反演是可以卷积加速的。

反演得到 \begin{Bmatrix}n \\ m\end{Bmatrix}m!=\sum\limits_{i=0}^m(-1)^{m-i}\dbinom{m}{i}i^n

\begin{Bmatrix}n \\ m\end{Bmatrix}m!=\sum\limits_{i=0}^m(-1)^{m-i}\dfrac{m!}{i!(m-i)!}i^n \begin{Bmatrix}n \\ m\end{Bmatrix}=\sum\limits_{i=0}^m\dfrac{(-1)^{m-i}}{(m-i)!}\dfrac{i^n}{i!}

设 g[x]=\begin{Bmatrix}n \\ x\end{Bmatrix};\ f[x]=\dfrac{x^n}{x!};\ h[x]=\dfrac{(-1)^x}{x!}

得到 g[m]=\sum\limits_{i=0}^mh[m-i]h[i] ，卷积即可。

Code:

#include<algorithm>
#include<cstdio>
#define mod 167772161
#define G 3
#define Maxn 200500
using namespace std;
int n,m,r[Maxn<<2];
long long f[Maxn<<2],g[Maxn<<2],invn,invG;
long long powM(long long a,long long t=mod-2)
{
  long long ans=1,buf=a;
  while(t){
    if(t&1)ans=(ans*buf)%mod;
    buf=(buf*buf)%mod;
    t>>=1;
  }return ans;
}
void NTT(long long *f,bool op)
{
  for (int i=0;i<n;i++)
    if (r[i]<i)swap(f[r[i]],f[i]);
  for (int len=1;len<n;len<<=1){
    int w=powM(op==1 ? G:invG,(mod-1)/len/2);
    for (int p=0;p<n;p+=len+len){
      long long buf=1;
      for (int i=p;i<p+len;i++){
        int sav=f[i+len]*buf%mod;
        f[i+len]=f[i]-sav;
        if (f[i+len]<0)f[i+len]+=mod;
        f[i]=f[i]+sav;
        if (f[i]>=mod)f[i]-=mod;
        buf=buf*w%mod;
        }//F(x)=FL(x^2)+x*FR(x^2)
         //F(W^k)=FL(w^k)+W^k*FR(w^k)
         //F(W^{k+n/2})=FL(w^k)-W^k*FR(w^k)
    }
  }
}
long long inv[Maxn];
int main()
{
  scanf("%d",&n);
  inv[0]=1;
  for (int i=1;i<=n;i++)
    inv[i]=inv[i-1]*powM(i)%mod;
  for (int i=0;i<=n;i++)
    f[i]=powM(i,n)*inv[i]%mod;
  for (int i=0;i<=n;i++)
    g[i]=powM(mod-1,i)*inv[i]%mod;
  n++;
  for(m=n+n,n=1;n<m;n<<=1);
  invn=powM(n);invG=powM(G);
  for (int i=0;i<n;i++)
    r[i]=(r[i>>1]>>1)|((i&1)?n>>1:0);
  NTT(f,1);
  NTT(g,1);
  for (int i=0;i<n;i++)f[i]=f[i]*g[i]%mod;
  NTT(f,0);
  for (int i=0;i<m/2;i++)
   printf("%lld ",f[i]*invn%mod);
  return 0;
}

前面提到了 x^{\overline n}=\sum\limits_{i=0}^n \begin{bmatrix}n\\i \end{bmatrix}x^i

将其视作关于 x 的多项式，显而易见求出 x^{\overline n} 的各项系数就求出了答案。

或者结合递推公式来理解。

每乘 $(x+(n-1))$ 一次就是批量转移一行，答案就是 $\prod\limits_{i=0}^{n-1}(x+i)=x^{\overline{n}}$ 。 分治 FFT 显然可以做到 $O(nlog^2n)$ ，但还不够优秀。 我们知道一个多项式 $F(x)$ 是可以快速求出 $F(x+c)$ 的，称之为多项式平移。 不了解这个技巧的可见前置芝士。 这样，利用 $x^{\overline{2n}}=x^{\overline n}*(x+n)^{\overline n}$ ，可以做到次数倍增。 要恰好倍增到某个次数，写法就和多项式次数倍增略有不同，要分奇偶讨论，具体见代码。 总复杂度 $T(n)=T(n/2)+O(n\log n)=O(n\log n)$ ,注意常数。 ```cpp #include<algorithm> #include<cstdio> #define mod 167772161 #define G 3 #define Maxn 270500 using namespace std; int r[Maxn<<2]; long long invn,invG; long long fac[Maxn],inv[Maxn]; long long powM(long long a,long long t=mod-2) { long long ans=1,buf=a; while(t){ if(t&1)ans=(ans*buf)%mod; buf=(buf*buf)%mod; t>>=1; }return ans; } void NTT(long long *f,bool op,int n) { for (int i=0;i<n;i++) if (r[i]<i)swap(f[r[i]],f[i]); for (int len=1;len<n;len<<=1){ int w=powM(op==1 ? G:invG,(mod-1)/len/2); for (int p=0;p<n;p+=len+len){ long long buf=1; for (int i=p;i<p+len;i++){ int sav=f[i+len]*buf%mod; f[i+len]=f[i]-sav; if (f[i+len]<0)f[i+len]+=mod; f[i]=f[i]+sav; if (f[i]>=mod)f[i]-=mod; buf=buf*w%mod; }//F(x)=FL(x^2)+x*FR(x^2) //F(W^k)=FL(w^k)+W^k*FR(w^k) //F(W^{k+n/2})=FL(w^k)-W^k*FR(w^k) } } } long long g[Maxn<<2]; //令f=f*g (mod x^lim) void times(long long *f,long long *gg,int len,int lim) { int m=len+len,n; for(int i=0;i<len;i++)g[i]=gg[i]; for(n=1;n<m;n<<=1);invn=powM(n); for(int i=len;i<n;i++)g[i]=0; for(int i=0;i<n;i++) r[i]=(r[i>>1]>>1)|((i&1)?n>>1:0); NTT(f,1,n);NTT(g,1,n); for(int i=0;i<n;++i)f[i]=(f[i]*g[i])%mod; NTT(f,0,n); for(int i=0;i<lim;++i)f[i]=(f[i]*invn)%mod; for(int i=lim;i<n;++i)f[i]=0; } long long p[Maxn]; //求出F(x+c) void fplus(long long *s,long long *f,int len,long long c) { for (int i=0;i<len;i++) p[len-i-1]=f[i]*fac[i]%mod; long long buf=1; for (int i=0;i<len;i++,buf=buf*c%mod) s[i]=buf*inv[i]%mod; times(p,s,len,len); for (int i=0;i<len;i++)s[len-i-1]=p[i]*inv[len-i-1]%mod; for (int i=len;i<len+len;i++)s[i]=0; } long long f[Maxn],s[Maxn]; void solve(long long *f,int n) { if (n==1){f[0]=0;f[1]=1;} else if (n&1){ solve(f,n-1);f[n]=0; //再乘上(x+n-1)就好了 for (int i=n;i>0;i--) f[i]=(f[i-1]+(n-1)*f[i])%mod; f[0]=f[0]*(n-1)%mod; }else { solve(f,n/2); //S(x)=F(x+n/2) fplus(s,f,n/2+1,n/2); times(f,s,n/2+1,n+1); } } int n,m; int main() { scanf("%d",&n);invG=powM(G); inv[1]=inv[0]=fac[0]=1; for (int i=1;i<=n;i++)fac[i]=fac[i-1]*i%mod; for (int i=2;i<=n;i++) inv[i]=inv[mod%i]*(mod-mod/i)%mod; for (int i=2;i<=n;i++)inv[i]=inv[i-1]*inv[i]%mod; solve(f,n); for (int i=0;i<=n;i++)printf("%lld ",f[i]); return 0; } ``` - [P5396 【模板】第二类斯特林数·列](https://www.luogu.org/problemnew/show/P5396) 没有什么优美的公式帮忙了……先来看看递推公式: $\begin{Bmatrix}n \\ m\end{Bmatrix}=\begin{Bmatrix}n-1 \\ m-1\end{Bmatrix}+m \begin{Bmatrix}n-1 \\ m\end{Bmatrix}

把答案设为多项式 H_k(x)=\sum\limits_{i=0}\begin{Bmatrix} i \\k \end{Bmatrix}x^i。

有 H_k(x)=\sum\limits_{i=0}\left(\begin{Bmatrix}i-1 \\ k-1\end{Bmatrix}+k \begin{Bmatrix}i-1 \\ k\end{Bmatrix}\right)x^i

=\sum\limits_{i=0}\begin{Bmatrix}i-1 \\ k-1\end{Bmatrix}x^i+k\sum\limits_{i=0}\begin{Bmatrix}i-1 \\ k\end{Bmatrix}x^i =x\sum\limits_{i=0}\begin{Bmatrix}i \\ k-1\end{Bmatrix}x^i+kx\sum\limits_{i=0}\begin{Bmatrix}i \\ k\end{Bmatrix}x^i=xH_{k-1}(x)+kxH_k(x)

即 H_k(x)=\dfrac{x}{1-kx}H_{k-1}(x)

边界是 H_0(x)=1。

则 H_k(x)=\prod\limits_{i=1}^k\dfrac{x}{1-ix}=x^k\left(\prod\limits_{i=1}^k(1-ix)\right)^{-1}

问题转化为计算 \prod\limits_{i=1}^k(1-ix)。

\prod\limits_{i=1}^k(1-ix)=\prod\limits_{i=1}^kx(\dfrac{1}{x}-i)=x^k\prod\limits_{i=1}^k(x^{-1}-i) \prod\limits_{i=1}^k(x-i)$ 翻转过来之后就是 $x^k\prod\limits_{i=1}^k(x^{-1}-i)

那么注意到\prod\limits_{i=1}^k(x-i)=\dfrac{x^{\underline{k+1}}}{x}我们只要求出x^{\underline{k+1}}就好。

我们还是考虑向上一题那样倍增 : x^{\underline{2n}}=x^{\underline n}*(x-n)^{\underline n}

复杂度还是 O(n\log n)。

其他细节都差不多，不过要加一个求逆，具体见代码。

#include<algorithm>
#include<cstdio>
#define mod 167772161
#define G 3
#define Maxn 270000
using namespace std;
int n,k,r[Maxn<<2];
long long invn,invG;
long long fac[Maxn],inv[Maxn];
long long powM(long long a,long long t=mod-2)
{
  long long ans=1,buf=a;
  while(t){
    if(t&1)ans=(ans*buf)%mod;
    buf=(buf*buf)%mod;
    t>>=1;
  }return ans;
}
void NTT(long long *f,bool op,int n)
{
  for (int i=0;i<n;i++)
    if (r[i]<i)swap(f[r[i]],f[i]);
  for (int len=1;len<n;len<<=1){
    int w=powM(op==1 ? G:invG,(mod-1)/len/2);
    for (int p=0;p<n;p+=len+len){
      long long buf=1;
      for (int i=p;i<p+len;i++){
        int sav=f[i+len]*buf%mod;
        f[i+len]=f[i]-sav;
        if (f[i+len]<0)f[i+len]+=mod;
        f[i]=f[i]+sav;
        if (f[i]>=mod)f[i]-=mod;
        buf=buf*w%mod;
        }//F(x)=FL(x^2)+x*FR(x^2)
         //F(W^k)=FL(w^k)+W^k*FR(w^k)
         //F(W^{k+n/2})=FL(w^k)-W^k*FR(w^k)
    }
  }
}
long long g[Maxn<<2];
void rev(long long *f,int len)
{
  for (int i=0;i<len;i++)g[i]=f[i];
  for (int i=0;i<len;i++)f[len-i-1]=g[i];
}
//令f=f*g (mod x^lim) 
void times(long long *f,long long *gg,int len,int lim)
{
  int m=len+len,n;
  for(int i=0;i<len;i++)g[i]=gg[i];
  for(n=1;n<m;n<<=1);invn=powM(n);
  for(int i=len;i<n;i++)g[i]=0;
  for(int i=0;i<n;i++)
    r[i]=(r[i>>1]>>1)|((i&1)?n>>1:0);
  NTT(f,1,n);NTT(g,1,n);
  for(int i=0;i<n;++i)f[i]=(f[i]*g[i])%mod;
  NTT(f,0,n);
  for(int i=0;i<lim;++i)f[i]=(f[i]*invn)%mod;
  for(int i=lim;i<n;++i)f[i]=0;
}
void Init(int lim)
{
  inv[1]=inv[0]=fac[0]=1;
  for (int i=1;i<=lim;i++)fac[i]=fac[i-1]*i%mod;
  for (int i=2;i<=lim;i++)
    inv[i]=inv[mod%i]*(mod-mod/i)%mod;
  for (int i=2;i<=lim;i++)inv[i]=inv[i-1]*inv[i]%mod;
}
long long p[Maxn<<2];
//求出F(x-c) 
void fminus(long long *s,long long *f,int len,int c)
{
  c=mod-c;
  for (int i=0;i<len;i++)
    p[len-i-1]=f[i]*fac[i]%mod;
  long long buf=1;
  for (int i=0;i<len;i++,buf=buf*c%mod)
    s[i]=buf*inv[i]%mod;
  times(p,s,len,len);
  for (int i=0;i<len;i++)s[len-i-1]=p[i]*inv[len-i-1]%mod;
  for (int i=len;i<len+len;i++)s[i]=0;
}
long long f[Maxn<<2],s[Maxn<<2];
void solve(long long *f,int n)
{
  if (n==1){f[0]=0;f[1]=1;}
  else if (n&1){
    solve(f,n-1);f[n]=0;
    //再乘上(x-n+1)就好了
    for (int i=n;i>0;i--)
      f[i]=(f[i-1]+(mod-n+1)*f[i])%mod;
    f[0]=f[0]*(mod-n+1)%mod;
  }else {
    solve(f,n/2);
    //S(x)=F(x+n/2)
    fminus(s,f,n/2+1,n/2);
    times(f,s,n/2+1,n+1);
  }
}
void invp(long long *f,int len)
{
  for (int i=0;i<k+1;i++)s[i]=p[i]=0;
  //注意清空 
  long long *r=s,*rr=p;
  int n=1;for(;n<len;n<<=1);
  rr[0]=powM(f[0]);
  for (int len=2;len<=n;len<<=1){
    for (int i=0;i<len;i++)
      r[i]=rr[i]*2%mod;
    times(rr,rr,len/2,len);
    times(rr,f,len,len);
    for (int i=0;i<len;i++)
      rr[i]=(r[i]-rr[i]+mod)%mod;
  }for (int i=0;i<len;i++)
    f[i]=rr[i];
}
int main()
{
  scanf("%d%d",&n,&k);
  if (k>n){
    for (int i=0;i<=n;i++)printf("0 ");
    return 0;
  }invG=powM(G);
  Init(k);solve(f,k+1);
  for (int i=0;i<k+1;i++)f[i]=f[i+1];
  rev(f,k+1);
  for (int i=n-k+1;i<k+1;i++)f[i]=0;
  for (int i=k+1;i<n-k+1;i++)f[i]=0;
  invp(f,n-k+1);
  for (int i=0;i<k;i++)printf("0 ");
  for (int i=0;i<n-k+1;i++)printf("%lld ",f[i]);
  return 0;
}

也可以使用一列第二类斯特林数的 EGF ：\frac{1}{k}(e^x-1)^k

需要多项式快速幂，常数较大，所以不推荐。

直接使用一列第一类斯特林数的 EGF : =\dfrac{\big(-\ln(1-x)\big)^k}{k!}。

多项式快速幂带走即可。

#include<algorithm>
#include<cstdio>
#define mod 167772161
#define G 3
#define Maxn 270500
using namespace std;
int n,k,r[Maxn<<2],invn,invG;
long long powM(long long a,long long t=mod-2)
{
  long long ans=1,buf=a;
  while(t){
    if(t&1)ans=(ans*buf)%mod;
    buf=(buf*buf)%mod;
    t>>=1;
  }return ans;
}
void NTT(long long *f,int n,short op)
{
  for (int i=0;i<n;i++)
   if (r[i]<i)swap(f[r[i]],f[i]);
  for (int p=2;p<=n;p<<=1){
    int len=p>>1,
        w=powM(op==1 ?G:invG,(mod-1)/p);
    for (int k=0;k<n;k+=p){
      long long buf=1;
      for (int i=k;i<k+len;i++){
        int sav=f[len+i]*buf%mod;
        f[len+i]=(f[i]-sav+mod)%mod;
        f[i]=(f[i]+sav)%mod;
        buf=buf*w%mod;
      }
    }
  }
}
long long _g[Maxn<<2];
void times(long long *f,long long *gg,int len1,int len2,int limit)
{
  int n=1;
  for(;n<len1+len2;n<<=1);
  long long *g=_g;
  for (int i=0;i<len2;i++)g[i]=gg[i];
  for (int i=len2;i<n;i++)g[i]=0;
  invn=powM(n);
  for(int i=0;i<n;i++)
   r[i]=(r[i>>1]>>1)|((i&1)?n>>1:0);
  NTT(f,n,1);NTT(g,n,1);
  for(int i=0;i<n;++i)f[i]=(f[i]*g[i])%mod;
  NTT(f,n,-1);
  for(int i=0;i<limit;++i)f[i]=f[i]*invn%mod;
  for(int i=limit;i<n;++i)f[i]=0;
}
long long _r[Maxn<<2],_rr[Maxn<<2];
void invp(long long *f,int m)
{
  long long *r=_r,*rr=_rr;
  int n=1;for(;n<m;n<<=1);
  rr[0]=powM(f[0]);
  for (int len=2;len<=n;len<<=1){
    for (int i=0;i<len;i++)
      r[i]=rr[i]*2%mod;
    times(rr,rr,len/2,len/2,len);
    times(rr,f,len,len,len);
    for (int i=0;i<len;i++)
      rr[i]=(r[i]-rr[i]+mod)%mod;
  }for (int i=0;i<m;i++)
    f[i]=rr[i];
  for (int i=0;i<n;i++)rr[i]=r[i]=0;
}
void dao(long long *f,int m)
{
  for (int i=1;i<m;i++)
   f[i-1]=f[i]*i%mod;
  f[m-1]=0;
}
void jifen(long long *f,int m)
{
  for (int i=m;i;i--)
   f[i]=f[i-1]*powM(i)%mod;
  f[0]=0;
}
long long _lns[Maxn<<2];
void lnp(long long *f,int m)
{
  long long *sav=_lns;
  for (int i=0;i<m;i++)sav[i]=f[i];
  invp(sav,m);dao(f,m);
  times(f,sav,m-1,m,m);
  jifen(f,m-1);
  for (int i=0;i<m;i++)sav[i]=0;
}
long long _xp[Maxn<<2],_xp2[Maxn<<2];
void exp(long long *f,int m)
{
  long long *s=_xp,*ss=_xp2;
  int n=1;for(;n<m;n<<=1);
  ss[0]=1;
  for (int len=2;len<=n;len<<=1){
    for (int i=0;i<len/2;i++)s[i]=ss[i];
    for (int i=len/2;i<len;i++)s[i]=0;
    lnp(s,len);
    for (int i=0;i<len;i++)
     s[i]=(f[i]-s[i]+mod)%mod;
    s[0]=(s[0]+1)%mod;
    times(ss,s,len/2,len,len);
  }for (int i=0;i<m;i++)f[i]=ss[i];
  for (int i=0;i<n;i++)s[i]=ss[i]=0;
}
long long fac[Maxn<<2],inv[Maxn<<2];
void Init(int lim)
{
  fac[0]=1;
  for (int i=1;i<=lim;i++)
    fac[i]=fac[i-1]*i%mod;
}
inline void print(long long *f,int len)
{
  for (int i=0;i<len;i++)
    printf("%lld ",f[i]);
  puts("");
}
long long a[Maxn<<2];
int main()
{
  invG=powM(G);
  scanf("%d%d",&n,&k);
  if (k>n){
    for(int i=0;i<=n;i++)printf("0 ");
    return 0;
  }n++;
  Init(n);
  for(int i=0;i<n-k;i++)a[i]=powM(i+1);
  lnp(a,n-k);
  for(int i=0;i<n-k;i++)a[i]=a[i]*k%mod;
  exp(a,n-k);
  for(int i=n-k-1;i>=0;i--)a[i+k]=a[i];
  for(int i=0;i<k;i++)a[i]=0;
  long long invk=1;
  for(int i=1;i<=k;i++)invk=invk*i%mod;
  invk=powM(invk);
  for(int i=0;i<n;i++)a[i]=a[i]*fac[i]%mod*invk%mod;
  print(a,n);
  return 0;
}

贝尓数 B[n]= 有 n 个不同元素的集合，将其分为任意个无序非空子集的方案数。

也就是 n 个球放进不限个数的若干个盒子里的方案数。

枚举盒子的个数就可以得到 {\rm Ans}=\sum\limits_{i=1}^n\begin{Bmatrix}n \\ i\end{Bmatrix} ,也就是一行第二类斯特林数的和。

这就得到了一个O(n\log n)求单个贝尔数的方法。

下面给出 O(n\log n) 求 B[1...n] 的做法。

单个集合的 EGF 为 \{0,1,1,1...\} 即 e^x-1。

使用 \exp 生成有标号集合，可得贝尓数的 EGF 为 \exp(e^x-1)。

也可以硬推一行第二类斯特林数的和。

B(x)=\sum\limits_{i=0}B[i]\dfrac{x^i}{i!}=\sum\limits_{i=0}\sum\limits_{j=0}\begin{Bmatrix}i \\ j\end{Bmatrix}\dfrac{x^i}{i!}=S_2(x,1)=\exp(e^x-1) 
#include<algorithm>
#include<cstdio>
#include<ctime>
#define ll long long
#define mod 998244353
#define G 3
#define Maxn 270000
using namespace std;
inline int read()
{
  register int X=0;
  register char ch=0;
  while(ch<48||ch>57)ch=getchar();
  while(ch>=48&&ch<=57)X=X*10+(ch^48),ch=getchar();
  return X;
}
ll powM(ll a,ll t=mod-2)
{
  ll ans=1,buf=a;
  while(t){
    if(t&1)ans=(ans*buf)%mod;
    buf=(buf*buf)%mod;
    t>>=1;
  }return ans;
}
int tr[Maxn<<1];
ll invG=powM(G);
void NTT(ll *f,short op,int n)
{
  for (int i=0;i<n;i++)
    if (i<tr[i])swap(f[i],f[tr[i]]);
  for(int p=2;p<=n;p<<=1){
    int len=p>>1,
        tG=powM(op==1 ? G:invG,(mod-1)/p);
    for (int k=0;k<n;k+=p){
      ll buf=1;
      for (int i=k;i<k+len;i++){
        int sav=f[len+i]*buf%mod;
        f[len+i]=(f[i]-sav+mod)%mod;
        f[i]=(f[i]+sav)%mod;
        buf=buf*tG%mod;
      }
    }
  }
}
ll _g1[Maxn<<1];
void times(ll *f,ll *g,int len1,int len2,int lim)
{
  int m=len1+len2-1,n;
  for(int i=0;i<len2;i++)_g1[i]=g[i];
  #define g _g1
  for(n=1;n<m;n<<=1);
  for(int i=len2;i<n;i++)g[i]=0;
  ll invn=powM(n);
  for(int i=0;i<n;i++)
    tr[i]=(tr[i>>1]>>1)|((i&1)?n>>1:0);
  NTT(f,1,n);NTT(g,1,n);
  for(int i=0;i<n;++i)f[i]=f[i]*g[i]%mod;
  NTT(f,-1,n);
  for(int i=0;i<lim;++i)
    f[i]=f[i]*invn%mod;
  for(int i=lim;i<n;++i)f[i]=0;
  #undef g
}
ll _w2[Maxn<<1],_r2[Maxn<<1];
void invp(ll *f,int m)
{
  int n;for (n=1;n<m;n<<=1);
  #define w _w2
  #define r _r2
  w[0]=powM(f[0]);
  for (int len=2;len<=n;len<<=1){
    for (int i=0;i<(len>>1);i++)
      r[i]=(w[i]<<1)%mod;
    times(w,w,len>>1,len>>1,len);
    times(w,f,len,len,len);
    for (int i=0;i<len;i++)
      w[i]=(r[i]-w[i]+mod)%mod;
  }for (int i=0;i<m;i++)f[i]=w[i];
  for (int i=0;i<n;i++)w[i]=r[i]=0;
  #undef w
  #undef r
}
void dao(ll *f,int m)
{
  for (int i=1;i<m;i++)
   f[i-1]=f[i]*i%mod;
  f[m-1]=0;
}
void jifen(ll *f,int m)
{
  for (int i=m;i;i--)
    f[i]=f[i-1]*powM(i)%mod;
  f[0]=0;
}
ll _s3[Maxn<<2];
void lnp(ll *f,int m)
{
  ll *sav=_s3;
  for (int i=0;i<m;i++)sav[i]=f[i];
  invp(sav,m);dao(f,m);
  times(f,sav,m-1,m,m);
  jifen(f,m-1);
  for (int i=0;i<m;i++)sav[i]=0;
}
ll _xp[Maxn<<2],_xp2[Maxn<<2];
void exp(ll *f,int m)
{
  ll *s=_xp,*s2=_xp2;
  int n=1;for(;n<m;n<<=1);
  s2[0]=1;
  for (int len=2;len<=n;len<<=1){
    for (int i=0;i<(len>>1);i++)s[i]=s2[i];
    lnp(s,len);
    for (int i=0;i<len;i++)
      s[i]=(f[i]-s[i]+mod)%mod;
    s[0]=(s[0]+1)%mod;
    times(s2,s,len>>1,len,len);
  }for (int i=0;i<m;i++)f[i]=s2[i];
  for (int i=0;i<n;i++)s[i]=s2[i]=0;
}
ll fac[Maxn],inv[Maxn];
void Init(int lim)
{
  fac[0]=1;
  for (int i=1;i<=lim;i++)
    fac[i]=fac[i-1]*i%mod;
  inv[lim]=powM(fac[lim]);
  for (int i=lim;i;i--)
    inv[i-1]=inv[i]*i%mod;
}
int T,n,m[1050];
ll F[Maxn<<1];
int main()
{
  T=read();
  for (int i=1;i<=T;i++)
    n=max(n,m[i]=read());
  Init(n);n++;
  for (int i=1;i<n;i++)F[i]=inv[i];
  exp(F,n);
  for (int i=1;i<=T;i++)
    printf("%lld\n",fac[m[i]]*F[m[i]]%mod);
  return 0;
}