ICode9

精准搜索请尝试: 精确搜索
首页 > 其他分享> 文章详细

斯特林数、欧拉数学习笔记

2021-03-10 19:04:34  阅读:277  来源: 互联网

标签:ch end 斯特林 sum 笔记 int rg now1 欧拉


第一类斯特林数

定义

\(\begin{bmatrix}n\\m\end{bmatrix}\) 表示将 \(n\) 个元素划分为 \(m\) 个圆排列的方案数,也可以写作 \(s(n,m)\)。

性质

\(x^{\overline n}=\sum_{i=0}^nx^is(n,i)\)

\(x^{\underline n}=\sum_{i=0}^n(-1)^{n-i}x^is(n,i)\)

\(\begin{aligned}n!=\sum_{i=0}^n s(n,i)\end{aligned}\)

一个排列唯一对应一个置换,而一个置换唯一对应一组轮换,

所以直接枚举划分为 \(i\) 个轮换的方案数,求和后就是排列总数。

求法

递推式

\(s(n,m)=s(n-1,m-1)+(n-1)s(n-1,m)\)

含义分别是把当前的元素作为一个新的圆排列和把当前元素插入到其它元素的后面。

第一类斯特林数·行

\(x^{\overline n}=\prod\limits_{i=0}^{n-1}(x+i)=\sum\limits_{i=0}^n\begin{bmatrix} n\\i\end{bmatrix}x^i\)

分治去乘可以做到 \(n log^2 n\),更优秀的做法是倍增。

因为 \(x^{\overline {2n}}=x^{\overline n} (x+n)^{\overline n}\),

所以我们需要由 \(f(x)\) 的值快速推出 \(f(x+d)\) 的值。

\(\begin{aligned} f(x+d) &=\sum_{i=0}^d a_i (x+d)^i \\ &=\sum_{i=0}^d a_i \sum_{j=0}^i C_i^j x^j d^{i-j} \\ &=\sum_{j=0}^d \sum_{i=j}^d a_iC_i^jx^j d^{i-j} \\ &=\sum_{j=0}^d \frac{x^j}{j!} \sum_{i=j}^d a_ii! \frac{d^{i-j}}{(i-j)!} \end{aligned}\)

后面的部分就是一个减法卷积的形式。

最终的时间复杂度为 \(O(nlogn)\)。

代码

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<iostream>
#define rg register
inline int read(){
	rg int x=0,fh=1;
	rg char ch=getchar();
	while(ch<'0' || ch>'9'){
		if(ch=='-') fh=-1;
		ch=getchar();
	}
	while(ch>='0' && ch<='9'){
		x=(x<<1)+(x<<3)+(ch^48);
		ch=getchar();
	}
	return x*fh;
}
const int maxn=524288,G=3,mod=167772161;
inline int addmod(rg int now1,rg int now2){
	return now1+=now2,now1>=mod?now1-mod:now1;
}
inline int delmod(rg int now1,rg int now2){
	return now1-=now2,now1<0?now1+mod:now1;
}
inline int mulmod(rg long long now1,rg int now2){
	return now1*=now2,now1>=mod?now1%mod:now1;
}
inline int ksm(rg int ds,rg int zs){
	rg int nans=1;
	while(zs){
		if(zs&1) nans=mulmod(nans,ds);
		ds=mulmod(ds,ds);
		zs>>=1;
	}
	return nans;
}
int w[23][maxn],wz[maxn];
void ntt(rg int A[],rg int lim,rg int typ){
	for(rg int i=0;i<lim;i++) if(i<wz[i]) std::swap(A[i],A[wz[i]]);
	for(rg int len=1,t0=0;len<lim;len<<=1,t0++){
		for(rg int j=0,now=len<<1;j<lim;j+=now){
			for(rg int k=0;k<len;k++){
				rg int x=A[j+k],y=mulmod(A[j+k+len],w[t0][k]);
				A[j+k]=addmod(x,y),A[j+k+len]=delmod(x,y);
			}
		}
	}
	if(typ==-1){
		std::reverse(A+1,A+lim);
		rg int ny=ksm(lim,mod-2);
		for(rg int i=0;i<lim;i++) A[i]=mulmod(A[i],ny);
	}
}
int n,jc[maxn],jcc[maxn],ny[maxn];
void pre(){
	ny[1]=1;
	for(rg int i=2;i<=n;i++) ny[i]=mulmod(mod-mod/i,ny[mod%i]);
	jc[0]=jcc[0]=1;
	for(rg int i=1;i<=n;i++) jc[i]=mulmod(jc[i-1],i),jcc[i]=mulmod(jcc[i-1],ny[i]);
}
int getC(rg int nn,rg int mm){
	return mulmod(jc[nn],mulmod(jcc[mm],jcc[nn-mm]));
}
int a[maxn],b[maxn],c[maxn],d[maxn];
void solve(rg int l,rg int r){
	if(l==r){
		a[0]=0,a[1]=1;
		return;
	}
	rg int mids=(l+r)>>1,len=r-l+1;
	if(len&1) mids--;
	solve(l,mids);
	len=mids+1;
	rg int lim=1,bit=0;
	for(;lim<=len+len;lim<<=1) bit++;
	for(rg int i=0;i<lim;i++) wz[i]=(wz[i>>1]>>1)|((i&1)<<(bit-1));
	for(rg int i=0,tmp=1;i<=len;i++,tmp=mulmod(tmp,len)) b[i]=mulmod(jcc[i],tmp);
	for(rg int i=0;i<=len;i++) c[i]=mulmod(jc[i],a[i]);
	for(rg int i=len+1;i<lim;i++) b[i]=c[i]=0;
	ntt(c,lim,-1),ntt(b,lim,1);
	for(rg int i=0;i<lim;i++) c[i]=mulmod(c[i],b[i]);
	ntt(c,lim,1);
	for(rg int i=0;i<=len;i++) c[i]=mulmod(c[i],jcc[i]);
	for(rg int i=len+1;i<lim;i++) c[i]=0;
	ntt(a,lim,1),ntt(c,lim,1);
	for(rg int i=0;i<lim;i++) a[i]=mulmod(a[i],c[i]);
	ntt(a,lim,-1);
	for(rg int i=r+2;i<lim;i++) a[i]=0;
	if((r-l+1)&1){
		for(rg int i=0;i<=r;i++) d[i+1]=a[i];
		for(rg int i=0;i<=r;i++) d[i]=addmod(d[i],mulmod(a[i],r));
		for(rg int i=0;i<=r+1;i++) a[i]=d[i];
	}
}
int main(){
	n=read();
	rg int lim=1;
	for(;lim<=n+n;lim<<=1);
	for(rg int len=1,t0=0;len<lim;len<<=1,t0++){
		w[t0][0]=1,w[t0][1]=ksm(G,(mod-1)/(len<<1));
		for(rg int i=2;i<len;i++) w[t0][i]=mulmod(w[t0][i-1],w[t0][1]);
	}
	pre();
	solve(0,n-1);
	for(rg int i=0;i<=n;i++) printf("%d ",a[i]);
	printf("\n");
	return 0;
}

第二类斯特林数

定义

\(\left\{\begin{matrix} n\\k\end{matrix}\right\}\)表示的是把 \(n\) 个不同的小球放在 \(m\) 个相同的盒子里的方案数,也可以写作 \(S(n,m)\)。

性质

\(n^k=\sum_ { i=0}^k S(k,i)×i!×C(n,i)\)

左边的含义就是 \(k\) 个球放在 \(n\) 个盒子里,

右边的含义是枚举有多少个盒子放。

也可以写成下降幂的形式:

\(n^k=\sum_{i=0}^k S(k,i) n^{\underline i}\)

求法

递推式

$S(n,m)=S(n-1,m-1)+mS(n-1,m) $

含义分别是新开了一个集合和在已有的集合中选择一个放。

容斥

\(S(n,m)=\frac{1}{m!}\sum_{i=0}^m(-1)^iC_m^i(m-i)^n\)

考虑用二项式反演。

假设一共有 \(n\) 个小球,

设 \(f[i]\) 为恰好有 \(i\) 个盒子放了小球的方案数,

\(g[i]\) 为至多有 \(i\) 个盒子放了小球的方案数。

那么 \(g[k]=\sum_{i=0}^k C_k^if[i]\)。

则 \(f[k]=\sum_{i=0}^k(-1)^{k-i}C_k^ig[i]=\sum_{i=0}^k(-1)^{k-i}C_k^i i^n\),

化得好看一点就成了 \(f[k]=\sum_{i=0}^k(-1)^iC_k^i (k-i)^n\) ,

那么 \(f[m]=\sum_{i=0}^m(-1)^iC_m^i(m-i)^n\)。

因为我们算组合数的时候是按照盒子不同算的,所以还要除以一个 \(m!\)。

第二类斯特林数·行

容斥的式子可以化成卷积的形式。

\(S(n,m)=\frac{1}{m!}\sum_{i=0}^m(-1)^iC_m^i(m-i)^n\)

\(S(n,m)=\sum_{i=0}^m\frac{(-1)^i}{i!}\frac{(m-i)^n}{(m-i)!}\)

就可以用 \(ntt\) 在 \(nlogn\) 的时间复杂度内求出一行的第二类斯特林数

代码

#include<cstdio>
#include<algorithm>
#include<cstring>
#define rg register
inline int read(){
	rg int x=0,fh=1;
	rg char ch=getchar();
	while(ch<'0' || ch>'9'){
		if(ch=='-') fh=-1;
		ch=getchar();
	}
	while(ch>='0' && ch<='9'){
		x=(x<<1)+(x<<3)+(ch^48);
		ch=getchar();
	}
	return x*fh;
}
const int maxn=524289,mod=167772161,G=3;
inline int addmod(rg int now1,rg int now2){
	return now1+=now2,now1>=mod?now1-mod:now1;
}
inline int delmod(rg int now1,rg int now2){
	return now1-=now2,now1<0?now1+mod:now1;
}
inline int mulmod(rg long long now1,rg int now2){
	return now1*=now2,now1>=mod?now1%mod:now1;
}
int ksm(rg int ds,rg int zs){
	rg int nans=1;
	while(zs){
		if(zs&1) nans=mulmod(nans,ds);
		ds=mulmod(ds,ds);
		zs>>=1;
	}
	return nans;
}
int w[25][maxn],wz[maxn];
void ntt(rg int A[],rg int lim,rg int typ){
	for(rg int i=0;i<lim;i++) if(i<wz[i]) std::swap(A[i],A[wz[i]]);
	for(rg int len=1,t0=0;len<lim;len<<=1,t0++){
		for(rg int j=0,now=len<<1;j<lim;j+=now){
			for(rg int k=0;k<len;k++){
				rg int x=A[j+k],y=mulmod(A[j+k+len],w[t0][k]);
				A[j+k]=addmod(x,y),A[j+k+len]=delmod(x,y);
			}
		}
	}
	if(typ==-1){
		std::reverse(A+1,A+lim);
		rg int ny=ksm(lim,mod-2);
		for(rg int i=0;i<lim;i++) A[i]=mulmod(A[i],ny);
	}
}
int a[maxn],b[maxn],jc[maxn],jcc[maxn],ny[maxn],n;
void pre(){
	ny[1]=1;
	for(rg int i=2;i<=n;i++) ny[i]=mulmod(mod-mod/i,ny[mod%i]);
	jc[0]=jcc[0]=1;
	for(rg int i=1;i<=n;i++) jc[i]=mulmod(jc[i-1],i),jcc[i]=mulmod(jcc[i-1],ny[i]);
}
int main(){
	n=read();
	pre();
	for(rg int i=0;i<=n;i++){
		if(i&1) a[i]=mod-jcc[i];
		else a[i]=jcc[i];
		b[i]=mulmod(ksm(i,n),jcc[i]);
	}
	rg int lim=1,bit=0;
	for(;lim<=n+n;lim<<=1) bit++;
	for(rg int i=0;i<lim;i++) wz[i]=(wz[i>>1]>>1)|((i&1)<<(bit-1));
	for(rg int len=1,t0=0;len<lim;len<<=1,t0++){
		w[t0][0]=1,w[t0][1]=ksm(G,(mod-1)/(len<<1));
		for(rg int i=2;i<len;i++) w[t0][i]=mulmod(w[t0][1],w[t0][i-1]);
	}
	ntt(a,lim,1),ntt(b,lim,1);
	for(rg int i=0;i<lim;i++) a[i]=mulmod(a[i],b[i]);
	ntt(a,lim,-1);
	for(rg int i=0;i<=n;i++) printf("%d ",a[i]);
	printf("\n");
	return 0;
}

第二类斯特林数·列

设 \(F(n,x)=\sum_{i} \begin{Bmatrix} i \\ n\end{Bmatrix} x^i\)。

\(\begin{aligned} F(n,x)&=\sum_{i} \begin{Bmatrix} i \\ n\end{Bmatrix}x^i \\ &=\sum_{i} (n\begin{Bmatrix} i-1 \\ n\end{Bmatrix}+\begin{Bmatrix} i-1 \\ n-1\end{Bmatrix})x^i \\ &=n\sum_{i}\begin{Bmatrix} i-1 \\ n\end{Bmatrix}x^i+\sum_{i} \begin{Bmatrix} i-1 \\ n-1\end{Bmatrix}x^i \\ &=nxF(n,x)+xF(n-1,x)\end{aligned}\)

所以 \(F(n,x)=\frac{x}{1-nx}F(n-1,x)\)。

又因为 \(F(0,x)=1\),

所以 \(F(n,x)=\frac{x^n}{\prod_{i=1}^n (1-ix)}\)。

把下面的部分分治乘起来再求逆即可。

代码

#include<cstdio>
#include<cmath>
#include<algorithm>
#include<iostream>
#include<vector>
#define rg register
inline int read(){
	rg int x=0,fh=1;
	rg char ch=getchar();
	while(ch<'0' || ch>'9'){
		if(ch=='-') fh=-1;
		ch=getchar();
	}
	while(ch>='0' && ch<='9'){
		x=(x<<1)+(x<<3)+(ch^48);
		ch=getchar();
	}
	return x*fh;
}
const int maxn=1<<19,mod=167772161,G=3;
inline int addmod(rg int now1,rg int now2){
	return now1+=now2,now1>=mod?now1-mod:now1;
}
inline int delmod(rg int now1,rg int now2){
	return now1-=now2,now1<0?now1+mod:now1;
}
inline int mulmod(rg long long now1,rg int now2){
	return now1*=now2,now1>=mod?now1%mod:now1;
}
int ksm(rg int ds,rg int zs){
	rg int nans=1;
	while(zs){
		if(zs&1) nans=mulmod(nans,ds);
		ds=mulmod(ds,ds);
		zs>>=1;
	}
	return nans;
}
int w[25][maxn],wz[maxn];
void ntt(std::vector<int>& A,rg int lim,rg int typ){
	for(rg int i=0;i<lim;i++) if(i<wz[i]) std::swap(A[i],A[wz[i]]);
	for(rg int len=1,t0=0;len<lim;len<<=1,t0++){
		for(rg int j=0,now=len<<1;j<lim;j+=now){
			for(rg int k=0;k<len;k++){
				rg int x=A[j+k],y=mulmod(w[t0][k],A[j+k+len]);
				A[j+k]=addmod(x,y),A[j+k+len]=delmod(x,y);
			}
		}
	}
	if(typ==-1){
		std::reverse(A.begin()+1,A.end());
		rg int ny=ksm(lim,mod-2);
		for(rg int i=0;i<lim;i++) A[i]=mulmod(A[i],ny);
	}
}
int n,k;
std::vector<int> g[maxn];
void dfs(rg int da,rg int l,rg int r){
	if(l==r){
		g[da].resize(2);
		g[da][0]=1,g[da][1]=mod-l;
		return;
	}
	rg int mids=(l+r)>>1;
	dfs(da<<1,l,mids);
	dfs(da<<1|1,mids+1,r);
	rg int lim=1,bit=0,len=r-l+1;
	for(;lim<=len+len;lim<<=1) bit++;
	for(rg int i=0;i<lim;i++) wz[i]=(wz[i>>1]>>1)|((i&1)<<(bit-1));
	g[da].resize(lim),g[da<<1].resize(lim),g[da<<1|1].resize(lim);
	ntt(g[da<<1],lim,1),ntt(g[da<<1|1],lim,1);
	for(rg int i=0;i<lim;i++) g[da][i]=mulmod(g[da<<1][i],g[da<<1|1][i]);
	ntt(g[da],lim,-1);
	for(rg int i=len+1;i<lim;i++) g[da][i]=0;
}
std::vector<int> finv,ans;
void getinv(std::vector<int>&A,std::vector<int>&B,rg int deg){
	if(deg==1){
		B.resize(1);
		B[0]=ksm(A[0],mod-2);
		return;
	}
	getinv(A,B,(deg+1)>>1);
	rg int lim=1,bit=0;
	for(;lim<deg+deg;lim<<=1) bit++;
	for(rg int i=0;i<lim;i++) wz[i]=(wz[i>>1]>>1)|((i&1)<<(bit-1));
	finv.resize(lim);
	if(A.size()<lim) A.resize(lim);
	for(rg int i=0;i<deg;i++) finv[i]=A[i];
	for(rg int i=deg;i<lim;i++) finv[i]=0;
	B.resize(lim);
	ntt(finv,lim,1),ntt(B,lim,1);
	for(rg int i=0;i<lim;i++) B[i]=delmod(mulmod(2,B[i]),mulmod(B[i],mulmod(B[i],finv[i])));
	ntt(B,lim,-1);
	for(rg int i=deg;i<lim;i++) B[i]=0;
}
int main(){
	n=read(),k=read();
	rg int lim=1;
	for(;lim<=n+n;lim<<=1);
	for(rg int len=1,t0=0;len<lim;len<<=1,t0++){
		w[t0][0]=1,w[t0][1]=ksm(G,(mod-1)/(len<<1));
		for(rg int i=2;i<len;i++) w[t0][i]=mulmod(w[t0][i-1],w[t0][1]);
	}
	dfs(1,1,k);
	if(n<k){
		for(rg int i=0;i<=n;i++) printf("0 ");
		printf("\n");
		return 0;
	}
	getinv(g[1],ans,n-k+2);
	for(rg int i=0;i<k && i<=n;i++) printf("0 ");
	for(rg int i=k;i<=n;i++) printf("%d ",ans[i-k]);
	printf("\n");
	return 0;
}

斯特林反演

\(F_i=\sum_{j=0}^{i}\begin{Bmatrix}i\\ j\end{Bmatrix}G_j\Leftrightarrow G_i=\sum_{j=0}^{i}(-1)^{i-j}\begin{bmatrix}i\\ j\end{bmatrix}F_j\)

\(F_i=\sum_{j=0}^{i}\begin{bmatrix}i\\ j\end{bmatrix}G_j\Leftrightarrow G_i=\sum_{j=0}^{i}(-1)^{i-j}\begin{Bmatrix}i\\ j\end{Bmatrix}F_j\)

欧拉数

定义

欧拉数 \(\langle\begin{smallmatrix}n\\ k\end{smallmatrix}\rangle\) 是 \(\{1,2,...,n\}\) 的有 \(k\) 个升高的排列 \(\pi_1 \pi_2 \cdots \pi_n\) 的个数,也就是说,在其中 \(k\) 个地方有 \(\pi_j<\pi_{j+1}\)

性质

\(\left\langle\begin{matrix}n\\m\end{matrix}\right\rangle=\sum_{k=0}^{m}\binom{n+1}{k}(-1)^k(m+1-k)^n\)

求法

递推式

\(\left\langle\begin{matrix}n\\m\end{matrix}\right\rangle=(n-m) \left\langle\begin{matrix}n-1 \\ m-1 \end{matrix}\right\rangle+(m+1)\left\langle\begin{matrix}n-1\\m\end{matrix}\right\rangle\)

强制从小到大放数。

如果放在最前面,升高的排列数不会增加,如果放在之前某个升高的排列中间,那么会减少一个排列增加一个排列,总数不变,否则放在其它的位置升高的排列数都会增加一个。

求一行欧拉数

构造两个多项式,令

\(f(x)=\sum_k \binom{n+1}{k}(-1)^kx^k\)

\(g(x)=\sum_k (k+1)^nx^k\)

这两个多项式卷积得到的结果就是第 \(n\) 行的欧拉数。

代码

#include<cstdio>
#include<cmath>
#include<algorithm>
#include<iostream>
#define rg register
inline int read(){
	rg int x=0,fh=1;
	rg char ch=getchar();
	while(ch<'0' || ch>'9'){
		if(ch=='-') fh=-1;
		ch=getchar();
	}
	while(ch>='0' && ch<='9'){
		x=(x<<1)+(x<<3)+(ch^48);
		ch=getchar();
	}
	return x*fh;
}
const int maxn=1<<19,mod=998244353,G=3;
inline int addmod(rg int now1,rg int now2){
	return now1+=now2,now1>=mod?now1-mod:now1;
}
inline int delmod(rg int now1,rg int now2){
	return now1-=now2,now1<0?now1+mod:now1;
}
inline int mulmod(rg long long now1,rg int now2){
	return now1*=now2,now1>=mod?now1%mod:now1;
}
int ksm(rg int ds,rg int zs){
	rg int nans=1;
	while(zs){
		if(zs&1) nans=mulmod(nans,ds);
		ds=mulmod(ds,ds);
		zs>>=1;
	}
	return nans;
}
int w[21][maxn],wz[maxn];
void ntt(rg int A[],rg int lim,rg int typ){
	for(rg int i=0;i<lim;i++) if(i<wz[i]) std::swap(A[i],A[wz[i]]);
	for(rg int len=1,t0=0;len<lim;len<<=1,t0++){
		for(rg int j=0,now=len<<1;j<lim;j+=now){
			for(rg int k=0;k<len;k++){
				rg int x=A[j+k],y=mulmod(w[t0][k],A[j+k+len]);
				A[j+k]=addmod(x,y),A[j+k+len]=delmod(x,y);
			}
		}
	}
	if(typ==-1){
		std::reverse(A+1,A+lim);
		rg int ny=ksm(lim,mod-2);
		for(rg int i=0;i<lim;i++) A[i]=mulmod(A[i],ny);
	}
}
int ny[maxn],jc[maxn],jcc[maxn];
int getC(rg int nn,rg int mm){
	return mulmod(jc[nn],mulmod(jcc[nn-mm],jcc[mm]));
}
int n,a[maxn],b[maxn];
int main(){
	n=read();
	rg int lim=1,bit=0;
	for(;lim<=n+n;lim<<=1) bit++;
	for(rg int i=0;i<lim;i++) wz[i]=(wz[i>>1]>>1)|((i&1)<<(bit-1));
	for(rg int len=1,t0=0;len<lim;len<<=1,t0++){
		w[t0][0]=1,w[t0][1]=ksm(G,(mod-1)/(len<<1));
		for(rg int i=2;i<len;i++) w[t0][i]=mulmod(w[t0][i-1],w[t0][1]);
	}
	ny[1]=1;
	for(rg int i=2;i<=n+1;i++) ny[i]=mulmod(mod-mod/i,ny[mod%i]);
	jc[0]=jcc[0]=1;
	for(rg int i=1;i<=n+1;i++) jc[i]=mulmod(jc[i-1],i),jcc[i]=mulmod(jcc[i-1],ny[i]);
	for(rg int i=0;i<=n;i++) a[i]=ksm(i+1,n),b[i]=i&1?mod-getC(n+1,i):getC(n+1,i);
	ntt(a,lim,1),ntt(b,lim,1);
	for(rg int i=0;i<lim;i++) a[i]=mulmod(a[i],b[i]);
	ntt(a,lim,-1);
	for(rg int i=0;i<=n;i++) printf("%d ",a[i]);
	printf("\n");
	return 0;
}

标签:ch,end,斯特林,sum,笔记,int,rg,now1,欧拉
来源: https://www.cnblogs.com/liuchanglc/p/14364020.html

本站声明: 1. iCode9 技术分享网(下文简称本站)提供的所有内容,仅供技术学习、探讨和分享;
2. 关于本站的所有留言、评论、转载及引用,纯属内容发起人的个人观点,与本站观点和立场无关;
3. 关于本站的所有言论和文字,纯属内容发起人的个人观点,与本站观点和立场无关;
4. 本站文章均是网友提供,不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属;如您发现该文章侵犯了您的权益,可联系我们第一时间进行删除;
5. 本站为非盈利性的个人网站,所有内容不会用来进行牟利,也不会利用任何形式的广告来间接获益,纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。

专注分享技术,共同学习,共同进步。侵权联系[81616952@qq.com]

Copyright (C)ICode9.com, All Rights Reserved.

ICode9版权所有