常见特殊数的多项式求法

斯特林数和贝尔数的求法

前置知识：

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

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

第二类斯特林数·行

\[\begin{Bmatrix}n\\m\end{Bmatrix}=\dfrac{1}{m!}\sum_{i=0}^{m}\binom{m}{i}(-1)^{m-i}i^n\\ =\sum_{i=0}^{m}\dfrac{i^n}{i!}\dfrac{(-1)^{m-i}}{(m-i)!} \]

卷积即可。

注意一下这个毒瘤模数（原根是 \(3\) ），数组长度是 \(2\times 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}\begin{bmatrix}n\\i\end{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);}//注意！这个函数左闭右闭

第一类斯特林数·列

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

注意 \(-\ln(1-x)=\sum\limits_{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;
}

第二类斯特林数·列

\[\begin{Bmatrix}n\\m\end{Bmatrix}=m\begin{Bmatrix}n-1\\m\end{Bmatrix}+\begin{Bmatrix}n-1\\m-1\end{Bmatrix}\\ \]

设答案为 \(F_{m}(x)=\sum\limits_{i}\begin{Bmatrix}i\\m\end{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)=\prod\limits_{i=1}^{m}(x-i)\) ，那么 \(G(\dfrac{1}{x})\) 就是把 \(G(x)\) 的系数翻转（reverse）的结果，而 \(G(\dfrac{1}{x})=\prod\limits_{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}\binom{n}{i}B_{i} \]

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

多项式科技

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

把它当作元素去生成集合，贝尔数的 \(\rm 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;
}

伯努利数

直接用\(\rm 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<\begin{matrix}n\\m\end{matrix}\right>=\sum_{k=0}^{m}\binom{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;
}

posted @ 2021-01-21 19:05 zzctommy 阅读(508) 评论(0) 编辑收藏举报

会员力量，点亮园子希望

刷新页面返回顶部

Loading

zzctommy