zoukankan      html  css  js  c++  java
  • 【洛谷5396】第二类斯特林数·列

    点此看题面

    • 给定(n,m),对于所有(i=0sim n),求出(S_2(i,m))
    • (n<131072)

    法一:倍增(NTT)

    关于第二类斯特林数的基础性质可见:斯特林数的基础性质与斯特林反演的初步入门

    我们知道它的递推公式:

    [S_2(n,m)=S_2(n-1,m-1)+m imes S_2(n-1,m) ]

    构造第(m)列第二类斯特林数的生成函数(F_m(x)=sum_{i=0}^{+infty}S_2(i,m)x^i),于是发现:

    [F_m(x)=F_{m-1}(x) imes x+F_m(x) imes mx\ F_m(x)=frac{x}{1-mx}F_{m-1}(x) ]

    得出结论:

    [F_m(x)=prod_{i=1}^mfrac x{1-ix} ]

    进一步推导,发现(frac x{1-ix}=frac{1}{x^{-1}-i})

    我们可以先求出(prod_{i=1}^m(x^{-1}-i))再求逆,而这家伙其实就是(prod_{i=1}^m(x-i))的转置。

    接下来的操作类似于【洛谷5408】第一类斯特林数·行,因为比较麻烦就不提供代码了(可参考那题代码)。

    法二:多项式快速幂

    考虑单独一个盒子,可以构造出第一列第二类斯特林数的指数型生成函数

    [F(x)=sum_{i=1}^nfrac{x^i}{i!} ]

    现在我们需要(m)个盒子,相当于是选出(m)个盒子拼在一起,然后除以(m!)去掉盒子的顺序,也就是说答案的指数型生成函数就应该是:

    [frac{F(x)^m}{m!} ]

    这里涉及到多项式快速幂,就类似于【洛谷5409】第一类斯特林数·列

    由于这种写法在拷了个多项式板子后显得非常好写,于是我选用了这种写法。

    代码:(O(nlogn))

    #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 131072
    #define X 167772161
    using namespace std;
    int n,m,f[N+5];
    I int QP(RI x,RI y) {RI t=1;W(y) y&1&&(t=1LL*t*x%X),x=1LL*x*x%X,y>>=1;return t;}
    int Fac[N+5],IFac[N+5];I void InitFac(CI S)
    {
    	RI i;for(Fac[0]=i=1;i<=S;++i) Fac[i]=1LL*Fac[i-1]*i%X;
    	for(IFac[i=S]=QP(Fac[S],X-2);i;--i) IFac[i-1]=1LL*IFac[i]*i%X;
    }
    namespace Poly//多项式模板
    {
    	#define PR 3
    	#define Init(n) P=1,L=0;W(P<=((n)<<1)) P<<=1,++L;for(i=0;i^P;++i) R[i]=R[i>>1]>>1|((i&1)<<L-1),A[i]=B[i]=0; 
    	int P,L,R[N<<2],A[N<<2],B[N<<2],t1[N+5],t2[N+5],t3[N+5];
    	I void NTT(int* s,CI op)
    	{
    		RI i,j,k,x,y,U,S;for(i=0;i^P;++i) i<R[i]&&(swap(s[i],s[R[i]]),0);
    		for(i=1;i^P;i<<=1) for(U=QP(QP(PR,op),(X-1)/(i<<1)),j=0;j^P;j+=i<<1) for(S=1,
    			k=0;k^i;++k,S=1LL*S*U%X) s[j+k]=((x=s[j+k])+(y=1LL*S*s[i+j+k]%X))%X,s[i+j+k]=(x-y+X)%X;
    		if(op==X-2) {RI t=QP(P,X-2);for(i=0;i^P;++i) s[i]=1LL*s[i]*t%X;}
    	}
    	I void Inv(CI n,int* a,int* b)
    	{
    		if(!n) return (void)(b[0]=QP(a[0],X-2));RI i;Inv(n>>1,a,b);
    		Init(n);for(i=0;i<=n;++i) A[i]=a[i],B[i]=b[i];NTT(A,1),NTT(B,1);
    		for(i=0;i^P;++i) A[i]=(2LL*B[i]-1LL*A[i]*B[i]%X*B[i]%X+X)%X;for(NTT(A,X-2),i=0;i<=n;++i) b[i]=A[i];
    	}
    	I void Ln(CI n,int* a,int* b)
    	{
    		RI i,j;for(i=0;i<=n;++i) b[i]=0;for(Inv(n,a,b),i=0;i^n;++i) t1[i]=1LL*a[i+1]*(i+1)%X;t1[n]=0;
    		Init(n);for(i=0;i<=n;++i) A[i]=b[i],B[i]=t1[i];NTT(A,1),NTT(B,1);
    		for(i=0;i^P;++i) A[i]=1LL*A[i]*B[i]%X;for(NTT(A,X-2),b[0]=0,i=1;i<=n;++i) b[i]=1LL*A[i-1]*QP(i,X-2)%X;
    	}
    	I void Exp(CI n,int* a,int* b)
    	{
    		if(!n) return (void)(b[0]=1);RI i;Exp(n>>1,a,b);Ln(n,b,t2);
    		Init(n);for(i=0;i<=n;++i) A[i]=b[i],B[i]=(!i-t2[i]+a[i]+X)%X;NTT(A,1),NTT(B,1);
    		for(i=0;i^P;++i) A[i]=1LL*A[i]*B[i]%X;for(NTT(A,X-2),i=0;i<=n;++i) b[i]=A[i];
    	}
    	I void Pow(CI n,int* a,CI k)
    	{
    		Ln(n,a,t3);for(RI i=0;i<=n;++i) t3[i]=1LL*t3[i]*k%X,a[i]=0;Exp(n,t3,a);
    	}
    }
    int main()
    {
    	RI i;for(scanf("%d%d",&n,&m),InitFac(n),i=0;i^n;++i) f[i]=IFac[i+1];//求出第一列除去原第0项的指数型生成函数
    	for(Poly::Pow(n-1,f,m),i=0;i<=min(n,m-1);++i) printf("0 ");//求出F(x)^m,在开头添m个0
    	for(i=0;i<=n-m;++i) printf("%d ",1LL*f[i]*IFac[m]%X*Fac[i+m]%X);return 0;//除以m!去掉盒子的顺序,由于是EGF答案要乘(i+m)!
    }
    
    败得义无反顾,弱得一无是处
  • 相关阅读:
    IDEA最常用快捷键汇总+快速写出Main函数
    设计模式之代理模式
    Java多线程中join、yield、sleep方法详解
    git基础命令详解
    用友网络科技Java高级开发面试题(2019)
    Java内部类超详细总结(含代码示例)
    构造器中绝对绝对不能调用可被覆盖的方法
    写给小白看的Git的安装配置与使用
    Vue学习笔记5--前端工程化
    Vue学习笔记4--vue-router路由
  • 原文地址:https://www.cnblogs.com/chenxiaoran666/p/S2_L.html
Copyright © 2011-2022 走看看