第二类斯特林数.行
根据通项公式有:
拆开来就可以变成卷积的形式.
第一类斯特林数.行:
第一类斯特林数实际上是上升幂展开成普通幂的系数,那么只需要多项式平移就可以倍增求出
#include<bits/stdc++.h> #define vi vector<int> using namespace std; const int mod=167772161; int ksm(int x,int t){ int tot=1; while(t){ if(t&1) tot=1ll*tot*x%mod; x=1ll*x*x%mod; t/=2; } return tot; } const int N=1<<19;//1<<17 int where[N],limit; int poor[N*5],inv[N],n,m,fac[N],finv[N],h[N]; int*w[2][19],*now=poor; void ad(int&x,int y){x=(x+y>=mod)?(x+y-mod):(x+y);} void pre(){ inv[1]=1;for(int i=2;i<=n;i++) inv[i]=1ll*inv[mod%i]*(mod-mod/i)%mod; fac[0]=1;for(int i=1;i<=n;i++) fac[i]=1ll*fac[i-1]*i%mod; finv[n]=ksm(fac[n],mod-2);for(int i=n-1;i>=0;i--) finv[i]=1ll*finv[i+1]*(i+1)%mod; for(int u=0,l=2;u<19;u++,l<<=1){ w[1][u]=now;now+=l+1; w[0][u]=now;now+=l+1; w[1][u][0]=1;w[1][u][1]=ksm(3,(mod-1)/l); for(int i=2;i<=l;i++) w[1][u][i]=1ll*w[1][u][i-1]*w[1][u][1]%mod; for(int i=0;i<=l;i++) w[0][u][i]=w[1][u][l-i]; } } void prepare(int n){ int l=0;limit=1; while(limit<=n) limit<<=1,l++; for(int i=0;i<limit;i++) where[i]=(where[i>>1]>>1)|((i&1)<<(l-1)); } void DFT(vi&now,int op){ int**w1=w[op]; for(int i=0;i<limit;i++) if(i<where[i]) swap(now[i],now[where[i]]); for(int l=1,u=0,r=2;r<=limit;l=r,r<<=1,u++){ int*w2=w1[u]; for(int i=0;i<limit;i+=r){ for(int v=0;v<l;v++){ int&x=now[i+v],&y=now[i+l+v],b=1ll*y*w2[v]%mod; y=(x>=b)?(x-b):(x+mod-b); x=(x+b>=mod)?(x+b-mod):(x+b); } } } if(!op){ int tmp=ksm(limit,mod-2); for(int i=0;i<limit;i++) now[i]=1ll*now[i]*tmp%mod; } } vi operator*(vi f,vi g){ int tot=f.size()+g.size()-2; prepare(tot);f.resize(limit);g.resize(limit); DFT(f,1);DFT(g,1); for(int i=0;i<limit;i++) f[i]=1ll*f[i]*g[i]%mod; DFT(f,0);f.resize(tot+1); return f; } vi A,B; void gas(vi&f){ int n=f.size()-1,m=n/2; if(n==1){f[1]=1;return ;} vi f0;f0.resize(m+1);gas(f0); vi g;g.resize(m+1); int t=1;f=f0; for(int i=0;i<=m;i++) f0[i]=1ll*fac[i]*f0[i]%mod,g[i]=1ll*t*finv[i]%mod,t=1ll*t*m%mod; reverse(g.begin(),g.end());g=f0*g; for(int i=0;i<=m;i++) g[i]=1ll*g[m+i]*finv[i]%mod; g.resize(m+1);f=f*g; if(n%2==1) { f.push_back(0); for(int i=n-1;i>=0;i--) ad(f[i+1],f[i]),f[i]=1ll*f[i]*(n-1)%mod; } } int main(){ scanf("%d",&n);pre(); A.resize(n+1);gas(A); for(int i=0;i<=n;i++) printf("%d ",A[i]); }第二类斯特林数.列
将递推公式写出来,设
显然有
那么就是要求一个前缀的(1-kx)卷积,多项式平移即可,最后讨论一下系数.
#include<bits/stdc++.h> #define vi vector<int> using namespace std; const int mod=167772161; int ksm(int x,int t){ int tot=1; while(t){ if(t&1) tot=1ll*tot*x%mod; x=1ll*x*x%mod; t/=2; } return tot; } const int N=1<<19;//1<<17 int where[N],limit; int poor[N*5],inv[N],n,m,fac[N],finv[N],h[N]; int*w[2][19],*now=poor; void ad(int&x,int y){x=(x+y>=mod)?(x+y-mod):(x+y);} void pre(){ inv[1]=1;for(int i=2;i<=n;i++) inv[i]=1ll*inv[mod%i]*(mod-mod/i)%mod; fac[0]=1;for(int i=1;i<=n;i++) fac[i]=1ll*fac[i-1]*i%mod; finv[n]=ksm(fac[n],mod-2);for(int i=n-1;i>=0;i--) finv[i]=1ll*finv[i+1]*(i+1)%mod; for(int u=0,l=2;u<19;u++,l<<=1){ w[1][u]=now;now+=l+1; w[0][u]=now;now+=l+1; w[1][u][0]=1;w[1][u][1]=ksm(3,(mod-1)/l); for(int i=2;i<=l;i++) w[1][u][i]=1ll*w[1][u][i-1]*w[1][u][1]%mod; for(int i=0;i<=l;i++) w[0][u][i]=w[1][u][l-i]; } } void prepare(int n){ int l=0;limit=1; while(limit<=n) limit<<=1,l++; for(int i=0;i<limit;i++) where[i]=(where[i>>1]>>1)|((i&1)<<(l-1)); } void DFT(vi&now,int op){ int**w1=w[op]; for(int i=0;i<limit;i++) if(i<where[i]) swap(now[i],now[where[i]]); for(int l=1,u=0,r=2;r<=limit;l=r,r<<=1,u++){ int*w2=w1[u]; for(int i=0;i<limit;i+=r){ for(int v=0;v<l;v++){ int&x=now[i+v],&y=now[i+l+v],b=1ll*y*w2[v]%mod; y=(x>=b)?(x-b):(x+mod-b); x=(x+b>=mod)?(x+b-mod):(x+b); } } } if(!op){ int tmp=ksm(limit,mod-2); for(int i=0;i<limit;i++) now[i]=1ll*now[i]*tmp%mod; } } vi operator*(vi f,vi g){ int tot=f.size()+g.size()-2; prepare(tot);f.resize(limit);g.resize(limit); DFT(f,1);DFT(g,1); for(int i=0;i<limit;i++) f[i]=1ll*f[i]*g[i]%mod; DFT(f,0);f.resize(tot+1); return f; } void INV(vi&f){ int n=f.size(); if(n==1) {f[0]=ksm(f[0],mod-2);return ;} vi f0=f;f0.resize((n+1)/2);INV(f0); prepare(n+2*(n+1)/2-3);f.resize(limit);f0.resize(limit); DFT(f,1);DFT(f0,1); for(int i=0;i<limit;i++) f[i]=1ll*f0[i]*(mod+2-1ll*f[i]*f0[i]%mod)%mod; DFT(f,0);f.resize(n); } vi A,B; void gas(vi&f){ int n=f.size()-1,m=n/2; if(n==1){f[1]=1;return ;} vi f0;f0.resize(m+1);gas(f0); vi g;g.resize(m+1); int t=1;f=f0; for(int i=0;i<=m;i++) f0[i]=1ll*fac[i]*f0[i]%mod,g[i]=1ll*t*finv[i]%mod,t=1ll*t*m%mod; reverse(g.begin(),g.end());g=f0*g; for(int i=0;i<=m;i++) g[i]=1ll*g[m+i]*finv[i]%mod; g.resize(m+1);f=f*g; if(n%2==1) { f.push_back(0); for(int i=n-1;i>=0;i--) ad(f[i+1],f[i]),f[i]=1ll*f[i]*(n-1)%mod; } } int main(){ scanf("%d %d",&n,&m);pre(); if(n<m){for(int i=0;i<=n;i++) printf("0 ");return 0;} A.resize(m+2);gas(A); int t=(m&1)?(mod-1):1; for(int i=0;i<=m;i++) A[i]=1ll*A[i+1]*t%mod,t=(t==1)?(mod-1):(1); A.pop_back();reverse(A.begin(),A.end()); A.resize(n-m+1);INV(A); for(int i=0;i<m;i++) printf("0 "); for(int i=0;i<n-m+1;i++) printf("%d ",A[i]); }第一类斯特林数.列
实际上就是求把若干个个元素拆成k个置换环的方案数,显然可以构造单个置换环的EGF,然后求这个EGF的k次,再除以k的阶乘,然后就可以的得到i个元素拆成k个置换环的EGF
#include<bits/stdc++.h> #define vi vector<int> using namespace std; const int mod=167772161; int ksm(int x,int t){ int tot=1; while(t){ if(t&1) tot=1ll*tot*x%mod; x=1ll*x*x%mod; t/=2; } return tot; } const int N=3000010,I=ksm(3,(mod-1)/4),Iinv=ksm(I,mod-2); int limit,inv[N],where[N],fac[N],finv[N],n,m; int poor[16778000]; int*w[2][22],*now=poor; void pre(){ inv[1]=1;for(int i=2;i<=3000000;i++) inv[i]=1ll*inv[mod%i]*(mod-mod/i)%mod; fac[0]=1;for(int i=1;i<=n;i++) fac[i]=1ll*fac[i-1]*i%mod; finv[n]=ksm(fac[n],mod-2);for(int i=n-1;i>=0;i--) finv[i]=1ll*finv[i+1]*(i+1)%mod; for(int l=2,u=0;u<=21;u++,l<<=1){ w[0][u]=now;now+=l+1; w[1][u]=now;now+=l+1; w[1][u][0]=1;w[1][u][1]=ksm(3,(mod-1)/l); for(int i=2;i<=l;i++) w[1][u][i]=1ll*w[1][u][i-1]*w[1][u][1]%mod; for(int i=0;i<=l;i++) w[0][u][i]=w[1][u][l-i]; } } void prepare(int n){ int l=0;limit=1; while(limit<=n) limit<<=1,l++; for(int i=0;i<limit;i++) where[i]=(where[i>>1]>>1)|((i&1)<<(l-1)); } void DFT(vi&now,int op){ for(int i=0;i<limit;i++) if(i<where[i]) swap(now[i],now[where[i]]); for(int l=1,u=0;l<limit;l<<=1,u++){ for(int i=0;i<limit;i+=l*2) for(int v=0;v<l;v++){ int x=i+v,y=i+l+v,a=now[x],b=1ll*now[y]*w[op][u][v]%mod; now[x]=a+b,now[x]>=mod?now[x]-=mod:0; now[y]=a+mod-b;now[y]>=mod?now[y]-=mod:0; } } if(!op){ int tmp=ksm(limit,mod-2); for(int i=0;i<limit;i++) now[i]=1ll*now[i]*tmp%mod; } } vi operator*(vi f,vi g){ int tot=f.size()+g.size()-2; prepare(tot);f.resize(limit);g.resize(limit); DFT(f,1);DFT(g,1); for(int i=0;i<limit;i++) f[i]=1ll*f[i]*g[i]%mod; DFT(f,0);f.resize(tot+1); return f; } vi operator+(vi f,vi g){ int n=max(f.size(),g.size()); f.resize(n);g.resize(n); for(int i=0;i<n;i++) f[i]=f[i]+g[i],f[i]>=mod?f[i]-=mod:0; return f; } vi operator-(vi f,vi g){ int n=max(f.size(),g.size()); f.resize(n);g.resize(n); for(int i=0;i<n;i++) f[i]=f[i]+mod-g[i],f[i]>=mod?f[i]-=mod:0; return f; } void INV(vi&f){ int n=f.size(); if(n==1) {f[0]=ksm(f[0],mod-2);return ;} vi f0=f;f0.resize((n+1)/2);INV(f0); prepare(n+2*(n+1)/2-3);f.resize(limit);f0.resize(limit); DFT(f,1);DFT(f0,1); for(int i=0;i<limit;i++) f[i]=1ll*f0[i]*(mod+2-1ll*f[i]*f0[i]%mod)%mod; DFT(f,0);f.resize(n); } void DERI(vi&f){for(int i=0;i<f.size()-1;i++) f[i]=1ll*f[i+1]*(i+1)%mod;f.pop_back();} void INTER(vi&f){f.push_back(0);for(int i=f.size()-1;i>=1;i--) f[i]=1ll*f[i-1]*inv[i]%mod;f[0]=0;} void LN(vi&f){int n=f.size();vi g=f;DERI(g);INV(f);f=f*g;INTER(f);f.resize(n);} void EXP(vi&f){ int n=f.size(); if(n==1) {f[0]=1;return ;} vi f0=f;f0.resize((n+1)/2);EXP(f0); vi g=f0;g.resize(n);LN(g); f=f0-(g-f)*f0;f.resize(n); } void KSM(vi&f,int a,int b,bool full){//k mod , k mod-1 , k >= n int n=f.size(),c,d=n; for(int i=0;i<n;i++) if(f[i]) {c=f[i],d=i;break;} if(d && full || 1ll*d*a>=n){for(int i=0;i<n;i++) f[i]=0;return ;} int tmp=ksm(c,mod-2);c=ksm(c,b); for(int i=0;i<n-d;i++) f[i]=1ll*f[i+d]*tmp%mod; f.resize(n-d);LN(f); for(int i=0;i<n-d;i++) f[i]=1ll*f[i]*a%mod; EXP(f);f.resize(n); for(int i=n-1;i>=d*a;i--) f[i]=1ll*f[i-d*a]*c%mod; for(int i=0;i<d*a;i++) f[i]=0; } int main(){ scanf("%d %d",&n,&m);pre(); vi A;A.resize(n+1);A[0]=0; for(int i=1;i<=n;i++) A[i]=1ll*fac[i-1]*finv[i]%mod; KSM(A,m,m,0); for(int i=0;i<=n;i++) printf("%lld ",1ll*A[i]*fac[i]%mod*finv[m]%mod); }
