zoukankan      html  css  js  c++  java
  • 【洛谷5408】第一类斯特林数·行

    点此看题面

    • 给定(n),对于所有(i=0sim n),求出(S_1(n,i))
    • (n<262144)

    从递推公式到多项式

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

    我们知道它的递推公式:

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

    考虑这个递推公式的多项式意义,把第二维看成次数,那么一种转移是次数不变系数为(1),一种转移是系数加(1)系数为(n-1)

    构造第(n)行斯特林数的生成函数(F_n(x)=sum_{i=0}^nS_1(n,i)x^i),于是发现(F_n(x))相当于在(F_{n-1}(x))的基础上卷上了(x+n-1)这个二项式。

    得出结论:

    [F_n(x)=prod_{i=0}^{n-1}(x+i) ]

    显然这东西可以直接(O(nlog^2n))的分治(NTT),但是不够优秀,会被卡。

    倍增求第一类斯特林数·行

    我们发现:

    [F_{2n}(x)=F_n(x)*F_n(x+n) ]

    要求(F_n(x+n)),设(F_n(x)=sum_{i=0}^nf_ix^i),这实际上是一个和斯特林数无关的经典多项式问题:

    [egin{align} F_n(x+n)&=sum_{i=0}^{n}f_i(x+n)^i\ &=sum_{i=0}^{n}f_isum_{j=0}^{n}C_i^jx^jn^{i-j}\ &=sum_{i=0}^{n}f_ii!sum_{j=0}^{n}frac{x^j}{j!}frac{n^{i-j}}{(i-j)!}\ &=sum_{i=0}^{n}frac{sum_{j=0}^{n}f_jj!frac{n^{j-i}}{(j-i)!}}{i!}x^i end{align} ]

    于是构造两个生成函数(A(x),B(x))卷积一下:

    [A(x)=sum_{i=0}^nf_ii!x^i\ B(x)=sum_{i=0}^nfrac{n^i}{i!}x^{n-i}\ [x^i]F_n(x+n)=frac{[x^{i+n}](A*B)}{i!} ]

    然后把(F_n(x))和求出的(F_n(x+n))卷起来即可得到(F_{2n}(x))

    复杂度分析:(T(n)=T(n/2)+O(nlogn)=O(nlogn))

    代码:(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 262144
    #define X 167772161
    using namespace std;
    int n,f[N+5],a[N+5],b[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()
    {
    	RI i;for(Fac[0]=i=1;i<=n;++i) Fac[i]=1LL*Fac[i-1]*i%X;
    	for(IFac[i=n]=QP(Fac[n],X-2);i;--i) IFac[i-1]=1LL*IFac[i]*i%X;
    }
    namespace Poly
    {
    	#define PR 3
    	int P,L,R[N<<2],A[N<<2],B[N<<2];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;
    	}
    	I void Mul(CI n,int* a,int* b)//卷积
    	{
    		RI i;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;
    		for(i=0;i<=n;++i) A[i]=a[i],B[i]=b[i];for(NTT(A,1),NTT(B,1),i=0;i^P;++i) A[i]=1LL*A[i]*B[i]%X;
    		RI t=QP(P,X-2);for(NTT(A,X-2),i=0;i<=(n<<1);++i) a[i]=1LL*A[i]*t%X;
    	}
    }
    I void Solve(CI n)//倍增
    {
    	if(n==1) return (void)(f[1]=1);//如果n=1
    	if(n&1) {Solve(n-(n&1));for(RI i=n;~i;--i) f[i]=(1LL*(n-1)*f[i]+(i?f[i-1]:0))%X;return;}//如果n为奇数,暴力卷上x+n-1
    	Solve(n>>1);for(RI i=0,p=1;i<=(n>>1);++i,p=1LL*p*(n>>1)%X) a[i]=1LL*f[i]*Fac[i]%X,b[(n>>1)-i]=1LL*p*IFac[i]%X;//构造A(x),B(x)以求F(x+n/2)
    	Poly::Mul(n>>1,a,b);for(RI i=0;i<=(n>>1);++i) b[i]=1LL*a[i+(n>>1)]*IFac[i]%X;Poly::Mul(n>>1,f,b);//F(x)*F(x+n/2)
    }
    int main()
    {
    	RI i;for(scanf("%d",&n),InitFac(),Solve(n),i=0;i<=n;++i) printf("%d ",f[i]);return 0;
    }
    
    败得义无反顾,弱得一无是处
  • 相关阅读:
    maven打包
    Description Resource Path Location Type Project configuration is not up-to-d
    GoldenGate
    maven打包 把要的依赖也打进去 配置
    mysql如何优化where子句
    根据状态计算操作状态
    kafka direct模式
    Kafka Connect
    Kafka Streams
    如何看源码
  • 原文地址:https://www.cnblogs.com/chenxiaoran666/p/S1_H.html
Copyright © 2011-2022 走看看