常见特殊数的多项式求法

斯特林数和贝尔数的求法

前置知识:

多项式笔记(一)

多项式笔记(二)

第一篇是基础多项式工业。

第二篇是关于各种数性质的推导,如果笔记里有式子这里就不会再证明了。

目录:

第二类斯特林数·行

[egin{Bmatrix}n\mend{Bmatrix}=dfrac{1}{m!}sum_{i=0}^{m}inom{m}{i}(-1)^{m-i}i^n\ =sum_{i=0}^{m}dfrac{i^n}{i!}dfrac{(-1)^{m-i}}{(m-i)!} ]

卷积即可。

注意一下这个毒瘤模数(原根是 (3) ),数组长度是 (2 imes 10^5) 不要开小。

void stirling2_row(int*g,int n){//注意!这个函数左闭右闭
	static int A[M],B[M];
	for(int i=0;i<=n;++i)A[i]=1ll*qpow(i,n)*math::ifc[i]%mod;
	for(int i=0;i<=n;++i)B[i]=i&1?mod-math::ifc[i]:math::ifc[i];
	poly_mul(A,B,A,n+1,n+1);
	for(int i=0;i<=n;++i)g[i]=A[i];
}

第一类斯特林数·行

[x^{overline{n}} =sum_{i=0}^{n}egin{bmatrix}n\iend{bmatrix}x^i ]

发现等式右边就是第一类斯特林数的生成函数。

所以求出 (x^{overline{n}}) 即可。

不知道是谁想出来的,可以倍增搞。

[x^{overline{n+1}}=x^{overline{n}}(x+n) ]

这个可以 (O(n)) ,暴力扫一遍即可。

[x^{overline{2n}}=x^{overline{n}}(x+n)^{overline{n}} ]

有个技巧叫做多项式平移,不会可以看开头的笔记一。

所以我们可以直接由 (x^{overline{n}}) 得到 ((x+n)^{overline{n}}) ,直接平移即可。

注意一下这个毒瘤模数(原根是 (3) ),数组长度是 (262144) 不要开小。

void poly_shift(int*g,int*f,int n,int c){
	static int A[M],B[M];
	for(int i=0;i<n;++i)A[i]=1ll*math::fac[i]*f[i]%mod;
	for(int i=0,j=1;i<n;++i,j=1ll*j*c%mod)B[i]=1ll*j*math::ifc[i]%mod;
	reverse(B,B+n),poly_mul(A,B,A,n,n);
	for(int i=0;i<n;++i)g[i]=1ll*A[n+i-1]*math::ifc[i]%mod;
}
void up_pow(int*g,int n){//注意!这个函数左闭右闭
	static int A[M];
	if(n==1)return g[0]=0,g[1]=1,void();
	if(n&1){
		up_pow(g,n-1);
		for(int i=n;i>0;--i)g[i]=(g[i-1]+1ll*(n-1)*g[i]%mod)%mod;
	}else{
		up_pow(g,n/2);
		clr(A,n/2+1),poly_shift(A,g,n/2+1,n/2),poly_mul(A,g,g,n/2+1,n/2+1);
	}
}
void stirling1_row(int*g,int n){up_pow(g,n);}//注意!这个函数左闭右闭

第一类斯特林数·列

第一类斯特林数列的 ( m EGF)(dfrac{(-ln(1-x))^k}{k!}) ,多项式快速幂即可。

注意 (-ln(1-x)=sumlimits_{i=1}dfrac{x^i}{i}) ,所以系数要左移一位不然没法快速幂。

注意一下这个毒瘤模数(原根是 (3) ),数组长度是 (131072) 不要开小。

void stirling1_column(int*g,int n,int k){//注意!这个函数左闭右闭
	if(n<k){
		for(int i=0;i<=n;++i)g[i]=0;
		return;
	}
	static int A[M];
	for(int i=0;i<=n;++i)A[i]=math::inv[i+1];
	clr(g,n+1),poly_qpow(g,A,k,n+1);
	for(int i=n;i>=k;--i)g[i]=g[i-k];
	for(int i=0;i<k;++i)g[i]=0;
	for(int i=k;i<=n;++i)g[i]=1ll*g[i]*math::fac[i]%mod*math::ifc[k]%mod;
}

第二类斯特林数·列

[egin{Bmatrix}n\mend{Bmatrix}=megin{Bmatrix}n-1\mend{Bmatrix}+egin{Bmatrix}n-1\m-1end{Bmatrix}\ ]

设答案为 (F_{m}(x)=sumlimits_{i}egin{Bmatrix}i\mend{Bmatrix}x^i)

那么

[F_{m}(x)=mxF_{m}(x)+xF_{m-1}(x)\ F_{m}(x)=dfrac{xF_{m-1}(x)}{1-mx} ]

所以现在要求的就是

[dfrac{x^m}{prod_{i=1}^{m}(1-ix)}=dfrac{1}{prod_{i=1}^{m}(frac{1}{x}-i)} ]

把分母算出来求逆再系数平移就好了。

(G(x)=prodlimits_{i=1}^{m}(x-i)) ,那么 (G(dfrac{1}{x})) 就是把 (G(x)) 的系数翻转(reverse)的结果,而 (G(dfrac{1}{x})=prodlimits_{i=1}^{m}(dfrac{1}{x}-i))(G(dfrac{1}{x})) 就是分母,也就是说我们求出 (G(x)) 就行了,显然 (G(x)=(x-1)^{underline{m}}=dfrac{x^{underline{m+1}}}{x}) 可以和上面一样倍增。

[x^{underline{2n}}=x^{underline{n}}(x-n)^{underline{n}}\ x^{underline{n+1}}=x^{underline{n}}(x-n) ]

void down_pow(int*g,int n){//注意!这个函数左闭右闭
	static int A[M];
	if(n==1)return g[0]=0,g[1]=1,void();
	if(n&1){
		down_pow(g,n-1);
		for(int i=n;i>0;--i)g[i]=(g[i-1]+1ll*(mod-n+1)*g[i]%mod)%mod;
	}else{
		down_pow(g,n/2);
		clr(A,n/2+1),poly_shift(A,g,n/2+1,mod-n/2),poly_mul(A,g,g,n/2+1,n/2+1);
	}
}
void stirling2_column(int*g,int n,int k){//注意!这个函数左闭右闭
	static int A[M];
	if(n<k){
		for(int i=0;i<=n;++i)g[i]=0;
		return;
	}
	clr(A,k+2),down_pow(A,k+1);
	for(int i=0;i<k+1;++i)A[i]=A[i+1];A[k+1]=0;
	reverse(A,A+k+1);
	clr(g,k+1),poly_inv(g,A,n-k+1);
	for(int i=n;i>=k;--i)g[i]=g[i-k];
	for(int i=0;i<k;++i)g[i]=0;
}

集合划分计数

贝尔数板子,就是让你求出贝尔数 (B_{1,cdots,100000})

  • 定义:(n) 个元素划分成任意个非空集合的方案数,就是一行第二类斯特林数的和。
  • 递推:

[B_{n+1}=sum_{i=0}^{n}inom{n}{i}B_{i} ]

就是枚举有 (n-i) 个数和第 (n) 个数在同一个集合,有 (dbinom{n}{i}) 种方法选出这么多数。剩下数的划分方案就是 (B_i)

  • 多项式科技

单个集合的 ( m EGF)({0,1,1,1,cdots}=e^x-1)

把它当作元素去生成集合,贝尔数的 ( m EGF) 就是 (exp(e^x-1))

void bell(int*g,int n){
	static int A[M];A[0]=0;
	for(int i=1;i<n;++i)A[i]=math::ifc[i];
	clr(g,n),poly_exp(g,A,n);
	for(int i=0;i<n;++i)g[i]=1ll*g[i]*math::fac[i]%mod;
}

伯努利数

直接用( m EGF)(dfrac{x}{e^x-1}) ,求逆即可。

void bernoulli(int*g,int n){
	static int A[M];
	for(int i=0;i<n;++i)A[i]=ifc[i+1];
	poly_inv(g,A,n);
	for(int i=0;i<n;++i)g[i]=1ll*g[i]*fac[i]%mod;
}

欧拉数 ·行

[left<egin{matrix}n\mend{matrix} ight>=sum_{k=0}^{m}inom{n+1}{k}(m+1-k)^n(-1)^{k} ]

(A_i=dbinom{n+1}{i}(-1)^{i},B_i=(i+1)^n) ,卷积即可。

void eula(int*g,int n){
	static int A[M];
	for(int i=0;i<n;++i)A[i]=qpow(i+1,n);
	for(int i=0;i<n;++i)g[i]=i&1?mod-comb(n+1,i):comb(n+1,i);
	poly_mul(g,A,g,n,n);
	for(int i=n;i<n<<1;++i)g[i]=0;
}

完整代码

附送这部分的完整代码。poly的其他部分在开头《多项式笔记(一)》里面有,加到poly这个namespace最后就可以用了。

void stirling2_row(int*g,int n){//注意!这个函数左闭右闭
	static int A[M],B[M];
	for(int i=0;i<=n;++i)A[i]=1ll*qpow(i,n)*math::ifc[i]%mod;
	for(int i=0;i<=n;++i)B[i]=i&1?mod-math::ifc[i]:math::ifc[i];
	poly_mul(A,B,A,n+1,n+1);
	for(int i=0;i<=n;++i)g[i]=A[i];
}
void poly_shift(int*g,int*f,int n,int c){
	static int A[M],B[M];
	for(int i=0;i<n;++i)A[i]=1ll*math::fac[i]*f[i]%mod;
	for(int i=0,j=1;i<n;++i,j=1ll*j*c%mod)B[i]=1ll*j*math::ifc[i]%mod;
	reverse(B,B+n),poly_mul(A,B,A,n,n);
	for(int i=0;i<n;++i)g[i]=1ll*A[n+i-1]*math::ifc[i]%mod;
}
void up_pow(int*g,int n){//注意!这个函数左闭右闭
	static int A[M];
	if(n==1)return g[0]=0,g[1]=1,void();
	if(n&1){
		up_pow(g,n-1);
		for(int i=n;i>0;--i)g[i]=(g[i-1]+1ll*(n-1)*g[i]%mod)%mod;
	}else{
		up_pow(g,n/2);
		clr(A,n/2+1),poly_shift(A,g,n/2+1,n/2),poly_mul(A,g,g,n/2+1,n/2+1);
	}
}
void stirling1_row(int*g,int n){up_pow(g,n);}//注意!这个函数左闭右闭
void stirling1_column(int*g,int n,int k){//注意!这个函数左闭右闭
	if(n<k){
		for(int i=0;i<=n;++i)g[i]=0;
		return;
	}
	static int A[M];
	for(int i=0;i<=n;++i)A[i]=math::inv[i+1];
	clr(g,n+1),poly_qpow(g,A,k,n+1);
	for(int i=n;i>=k;--i)g[i]=g[i-k];
	for(int i=0;i<k;++i)g[i]=0;
	for(int i=k;i<=n;++i)g[i]=1ll*g[i]*math::fac[i]%mod*math::ifc[k]%mod;
}
void down_pow(int*g,int n){//注意!这个函数左闭右闭
	static int A[M];
	if(n==1)return g[0]=0,g[1]=1,void();
	if(n&1){
		down_pow(g,n-1);
		for(int i=n;i>0;--i)g[i]=(g[i-1]+1ll*(mod-n+1)*g[i]%mod)%mod;
	}else{
		down_pow(g,n/2);
		clr(A,n/2+1),poly_shift(A,g,n/2+1,mod-n/2),poly_mul(A,g,g,n/2+1,n/2+1);
	}
}
void stirling2_column(int*g,int n,int k){//注意!这个函数左闭右闭
	static int A[M];
	if(n<k){
		for(int i=0;i<=n;++i)g[i]=0;
		return;
	}
	clr(A,k+2),down_pow(A,k+1);
	for(int i=0;i<k+1;++i)A[i]=A[i+1];A[k+1]=0;
	reverse(A,A+k+1);
	clr(g,k+1),poly_inv(g,A,n-k+1);
	for(int i=n;i>=k;--i)g[i]=g[i-k];
	for(int i=0;i<k;++i)g[i]=0;
}
void bell(int*g,int n){
	static int A[M];A[0]=0;
	for(int i=1;i<n;++i)A[i]=math::ifc[i];
	clr(g,n),poly_exp(g,A,n);
	for(int i=0;i<n;++i)g[i]=1ll*g[i]*math::fac[i]%mod;
}
void bernoulli(int*g,int n){
	static int A[M];
	for(int i=0;i<n;++i)A[i]=ifc[i+1];
	poly_inv(g,A,n);
	for(int i=0;i<n;++i)g[i]=1ll*g[i]*fac[i]%mod;
}
void eula(int*g,int n){
	static int A[M];
	for(int i=0;i<n;++i)A[i]=qpow(i+1,n);
	for(int i=0;i<n;++i)g[i]=i&1?mod-comb(n+1,i):comb(n+1,i);
	poly_mul(g,A,g,n,n);
	for(int i=n;i<n<<1;++i)g[i]=0;
}
原文地址:https://www.cnblogs.com/zzctommy/p/14274228.html