笔记 组合数学魔法: 斯特林数
笔记 组合数学魔法: 斯特林数
最基本的
斯特林数是什么?
和组合数类似(写法也类似),表示一种计数
有两类斯特林数
第一类:\({n}\brack{m}\),也读作“n 轮换 m”,表示把 \(n\) 划分成 \(m\) 个环排列数的方案数。
第二类:\({n}\brace{m}\),也读作“n 子集 m”,表示把 \(n\) 划分成 \(m\) 个子集的方案数。
如何求?
存在递推。
第一类:考虑前 \(n-1\) 个,如果已经放进了 \(m\) 个圆排列,考虑最后一个插进来的方案数即可。每个排列里的每个位置都可以插入,并且都是不一样的方案(显然)。那么这样的方案数就是 \({n-1\brack m}\times (n-1)\)。另一种,前面只放了 \(m-1\) 个圆排列,最后一个单独出来:\({n-1\brack m-1}\)。综合两种,有:
第二类:同样这么考虑,易得:
魔法
第二类 & 下降幂
下降幂:\(x^{\underline{n}}=x(x-1)(x-2)...(x-n+1)\),即 \(x\) 往下 \(n\) 个数字的连乘积。
(后面的上升幂类似,就是往上 \(n\) 个,记作 \(x^{\overline{n}}\))
证:考虑组合意义
\(n^m\) 即选 \(n\) 个数,\(m\) 次,每次选出来后不会拿走这个数。于是每次都有 \(n\) 种选法,方法数当然是 \(n^m\)。
这个转化在我的 这篇文章 也有
然后这样显然可能会有重复的。枚举 \(m\) 次选出来的数中,去重过后剩下 \(i\) 个。
那 \(m\) 次选择一定可以划分出 \(i\) 个等价类 (即相同的分到一类,分出来的结果)。
然后先选出等价类的分布,相当于 \(m\) 分成 \(i\) 个子集,就是 \({m\brace i}\)
然后考虑等价类中的选择:第一个类能选 \(n\) 个,但是第二类就只能选 \(n-1\) 个,第三类就只能选 \(n-2\) 个...这个方案显然是 \(n^{\underline{i}}\)
然后就有了这个式子。
第一类 & 上升幂
形象记忆:上一个式子每个幂都往上变一下,然后斯特林数第二类变第一类
这个证明就不好考虑组合意义了 (因为第一类斯特林数很逊),但是可以归纳 —— 证明略 (懒得打了qaq)
反转公式
俩证明差不多(因为第一类第二类斯特林数有很优美的对称性,一个证明稍微改一下就可以证另一个)
证:
首先,\(a^{\underline{b}}=(-1)^b(-a)^{\overline{b}}\)
显然吧 (每一项取负即可)
考虑 \(n^m\)
\[=\sum\limits_{i=0}^{m} {m\brace i}n^{\underline{i}}\\ =\sum\limits_{i=0}^{m} {m\brace i}(-n)^{\overline{i}}\times (-1)^{i}\\ =\sum\limits_{i=0}^{m} {m\brace i} \sum\limits_{j=0}^{i} (-n)^{j}{i\brack j}\times (-1)^i\\ =\sum\limits_{i=0}^{m} \sum\limits_{j=0}^{i} {m\brace i}{i\brack j}n^j (-1)^{i+j}\\ =\sum\limits_{j=0}^{m} n^j \sum\limits_{i=j}^{m}{m\brace i}{i\brack j}(-1)^{i+j} \]然后发现这玩意恰好是 \(j\) 从 \(0\) 到 \(m\) 求和,\(n^j\) 后面跟了一堆系数。
对于任意的 \(n,m\) ,它应该都成立。所以后面的系数,当且仅当 \(j=m\) 的时候能取到 \(1\) ,否则都是 \(0\)。
后面的系数就是斯特林反转公式的式子。
证毕。
备注:\((-1)^{i+j}\) 有时也被写作 \((-1)^{i-j}\)。它俩显然等价。
斯特林反演
证:很 simple,把右式暴力代入左式,瞎 jb 换求和符号,你会发现一个斯特林反转公式 —— 然后就证完了。
例题
arc096C —— 基本应用
设 \(S=\{1,2...n\}\)。
求:从 \(S\) 集合中选出若干个子集,使得 \([1,n]\) 中所有元素在所有子集中出现的次数都 \(\ge 2\) 的方案数。
\(n\le 3000\)。模大质数(给定)。
考虑容斥。设 \(f(x)\) 表示选若干子集,使得有 \(x\) 个元素出现次数 \(<2\) 的方案数。然后原答案为
然后考虑 \(f(x)\) 如何求。显然划分子集数不知道,我们要枚举子集数 \(m\)
首先钦定 \(x\) 个位置 \(<2\)。然后剩下 \(n-x\) 个随便选。
钦定 \(x\) 个位置 \(<2\),相当于出现 \(1\) 次或者不出现。然后还要将出现了的划分到 \(m\) 个子集中。
考虑那些不出现的,我们新建一个“垃圾桶集合”:独立于 \(m\) 个集合,如果有元素被划分到了这个集合,说明它被丢掉了。
但是垃圾桶集合可能为空,与斯特林数的定义相违背。一个巧妙的处理方法:加入一个不会出现的元素 \(114514\),谁和 \(114514\) 分在一组,说明它被丢掉了。化粪池集合
这样显然就不会有空的集合了(因为就算没有元素被丢掉,\(114514\) 的存在也保证了化粪池集合非空)。
这样就有 \(x+1\) 个元素,和 \(m+1\) 个集合。这一部分答案是 \({x+1\brace m+1}\)
所以这个处理的精髓在于加入了一个辅助元素,而它是多少都没必要,只要没有歧义就行
比如你可以说它是 \(-1\),\(0\),或者 \(1919810\)
然后“随便选”的方案显然是 \(2^{2^{n-x}}\)
有 \(2^{n-x}\) 个子集,每个自己可以决定是否选,所以再来一层指数
然后这样直接算就行了,是 \(O(n^2)\) 的。
枚举 \(i\) ,枚举 \(m\) ,两层枚举。
代码:
#include <bits/stdc++.h>
using namespace std;
namespace Flandre_Scarlet
{
#define N 3003
#define int long long
#define F(i,l,r) for(int i=l;i<=r;++i)
#define D(i,r,l) for(int i=r;i>=l;--i)
#define Fs(i,l,r,c) for(int i=l;i<=r;c)
#define Ds(i,r,l,c) for(int i=r;i>=l;c)
#define MEM(x,a) memset(x,a,sizeof(x))
#define FK(x) MEM(x,0)
#define Tra(i,u) for(int i=G.st(u),v=G.to(i);~i;i=G.nx(i),v=G.to(i))
#define p_b push_back
#define sz(a) ((int)a.size())
#define all(a) a.begin(),a.end()
#define iter(a,p) (a.begin()+p)
#define PUT(a,n) F(i,1,n) printf("%d ",a[i]); puts("");
int I() {char c=getchar(); int x=0; int f=1; while(c<'0' or c>'9') f=(c=='-')?-1:1,c=getchar(); while(c>='0' and c<='9') x=(x<<1)+(x<<3)+(c^48),c=getchar(); return ((f==1)?x:-x);}
template <typename T> void Rd(T& arg){arg=I();}
template <typename T,typename...Types> void Rd(T& arg,Types&...args){arg=I(); Rd(args...);}
void RA(int *p,int n) {F(i,1,n) *p=I(),++p;}
int n,mod;
void Input()
{
Rd(n,mod);
}
int cc[N][N],s2[N][N];
void init()
{
int n=3001;
cc[0][0]=1;
F(i,1,n)
{
cc[i][0]=1;
F(j,1,i) cc[i][j]=(cc[i-1][j-1]+cc[i-1][j])%mod;
}
s2[0][0]=1;
F(i,1,n)
{
F(j,1,i)
{
s2[i][j]=(s2[i-1][j-1]+s2[i-1][j]*j%mod)%mod;
}
}
}
int qpow(int a,int b,int m=mod) {int r=1; while(b) {if (b&1) r=r*a%m; a=a*a%m,b>>=1;} return r;}
void Sakuya()
{
init();
int ans=0;
F(i,0,n)
{
int sum=0;
F(k,0,i)
{
sum=(sum+s2[i+1][k+1]*qpow(qpow(2,n-i),k)%mod)%mod;
}
sum%=mod;
sum=sum*cc[n][i]%mod*qpow(2,qpow(2,n-i,mod-1))%mod;
if (i&1) ans=(ans-sum+mod)%mod;
else ans=(ans+sum)%mod;
}
printf("%lld\n",ans);
}
void IsMyWife()
{
Input();
Sakuya();
}
}
#undef int //long long
int main()
{
Flandre_Scarlet::IsMyWife();
getchar();
return 0;
}
hdu4625 —— 斯特林变形幂
给一颗树,和常数 \(k\) ,对于每个点 \(u\),求:
\(\sum\limits_{v\neq u} dis(u,v)^k\)
复杂度 \(O(nk)\),模 \(10007\)
用第二类斯特林数把幂变成下降幂
\(n^{\underline{m}}=\dfrac{\binom{n}{m}}{m!}\)。
然后转化成:维护 \(\sum\limits_{v\neq u} \binom{dis(u,v)}{i}\),对于每个 \(i\in [0,k]\)
换根 \(dp\) 一遍,接下来考虑维护子树中这玩意。
设子树中这个和为 \(f(u,i)=\sum\limits_{v \in \operatorname{subtree} u} \binom{dis(u,v)}{i}\)
从儿子到自己,其实就是每个距离都 \(+1\)。
\(\binom{dis}{i}=\binom{dis-1}{i}+\binom{dis-1}{i-1}\)(组合数递推)
每个都这么变,得: \(f(u,i)=\sum\limits_{v\in \operatorname{son}u}f(v,i)+f(v,i-1)\)
它显然可以去除贡献,然后就可以换根 \(dp\) 了。
预处理斯特林数和阶乘,两遍 \(DFS\)。复杂度显然是 \(O(nk)\) 的。
注意边界 (转移的时候)
代码:
#include <bits/stdc++.h>
using namespace std;
namespace Flandre_Scarlet
{
#define N 50004
#define K 503
#define mod 10007
#define F(i,l,r) for(int i=l;i<=r;++i)
#define D(i,r,l) for(int i=r;i>=l;--i)
#define Fs(i,l,r,c) for(int i=l;i<=r;c)
#define Ds(i,r,l,c) for(int i=r;i>=l;c)
#define MEM(x,a) memset(x,a,sizeof(x))
#define FK(x) MEM(x,0)
#define Tra(i,u) for(int i=G.st(u),v=G.to(i);~i;i=G.nx(i),v=G.to(i))
#define p_b push_back
#define sz(a) ((int)a.size())
#define all(a) a.begin(),a.end()
#define iter(a,p) (a.begin()+p)
#define PUT(a,n) F(i,1,n) printf("%d ",a[i]); puts("");
int I() {char c=getchar(); int x=0; int f=1; while(c<'0' or c>'9') f=(c=='-')?-1:1,c=getchar(); while(c>='0' and c<='9') x=(x<<1)+(x<<3)+(c^48),c=getchar(); return ((f==1)?x:-x);}
template <typename T> void Rd(T& arg){arg=I();}
template <typename T,typename...Types> void Rd(T& arg,Types&...args){arg=I(); Rd(args...);}
void RA(int *p,int n) {F(i,1,n) *p=I(),++p;}
class Graph
{
public:
struct edge{int v,nx;} pe[N<<1]; edge *e;
int ph[N]; int* head;
int ecnt;
void clear()
{
MEM(pe,-1); e=pe+2;
MEM(ph,-1); head=ph+2;
ecnt=-1;
}
int& st(int u) {return head[u];}
int& to(int i) {return e[i].v;}
int& nx(int i) {return e[i].nx;}
void add(int u,int v)
{
e[++ecnt]=(edge){v,head[u]}; head[u]=ecnt;
}
void add2(int u,int v) {add(u,v); add(v,u);}
}G;
int n,k;
void Input()
{
Rd(n,k); G.clear();
F(i,1,n-1)
{
int u,v; Rd(u,v);
G.add2(u,v);
}
}
int S2[K][K],fc[K];
int dp[N][K];
void DFS(int u,int f)
{
dp[u][0]=1;
Tra(i,u) if (v!=f)
{
DFS(v,u);
F(c,0,k)
{
dp[u][c]=(dp[u][c]+dp[v][c]+(c?dp[v][c-1]:0))%mod;
}
// C(dis,c)=C(dis-1,c)+C(dis,c-1)
}
}
int tmp[K];
void DFS2(int u,int f) // 换根dp
{
if (u!=f) // 非根
{
F(i,0,k) tmp[i]=dp[f][i];
F(c,0,k) tmp[c]=(tmp[c]-(dp[u][c]+(c?dp[u][c-1]:0))%mod+mod)%mod; // 去掉贡献
F(c,0,k) dp[u][c]=(dp[u][c]+tmp[c]+(c?tmp[c-1]:0))%mod; // 把贡献加到自己上来
}
Tra(i,u) if (v!=f)
{
DFS2(v,u);
}
}
void Sakuya()
{
fc[0]=1;
F(i,1,k)
{
S2[i][1]=1;
fc[i]=fc[i-1]*i%mod;
F(j,2,i)
{
S2[i][j]=(S2[i-1][j-1]+S2[i-1][j]*j%mod)%mod;
}
}
// 预处理
FK(dp);
DFS(1,1);
DFS2(1,1);
F(i,1,n)
{
int ans=0;
F(j,0,k) ans=(ans+S2[k][j]*dp[i][j]%mod*fc[j]%mod)%mod;
printf("%d\n",ans);
}
}
void IsMyWife()
{
int t=I();
while(t-->0)
{
Input();
Sakuya();
}
}
}
#undef int //long long
int main()
{
Flandre_Scarlet::IsMyWife();
getchar();
return 0;
}
bzoj 4671 —— 斯特林反演魔法
有 \(s\) 个图,每个图都有 \(n\) 个点。然后用一个压缩过的长度为 \(n(n-1)/2\) 的 01 串,来表示每条边是否存在
压缩方式:
for i = 1 to n do
for j = i + 1 to n do
if G contains edge (i, j) then
print 1
else
print 0
end if
end for
end for
然后两张图的“异或”定义为,两个图的压缩01串异或起来再解压回去
求有多少个选择图的方案,使得选出来所有图异或出来是一个连通图。
\(n\le 10,s\le 60\)
草这个题真的太 NB 了
设 \(g(m)\) 表示图有 \(m\)个联通块的方案数。
设 \(f(m)\) 表示图能划分成 \(m\) 个集合,集合之间两两没有边的方案数。
因为集合内部也可能不连通,所以 \(f\) 的限制比 \(g\) 弱,也就是 \(f(m)\ge g(m)\)
然后考虑它俩的关系。\(f\) 中应该有 \(\ge m\) 个连通块,枚举连通块数量为 \(i\),加起来得
反演
我们要求的是 \(g(1)=\)
换句话说只要能求 \(f\) 就行了
\(f\) 如何求呢?
首先爆搜划分集合方案
总划分方案是贝尔数,最大是第 \(10\) 项,1e5 级别的
对于一个方案,\(O(n^2)\) 枚举一条边,使得它两个点不在一个集合里,也就是这条边一定不存在。
然后去找给定的图中,哪些图包含了这些边。设 \(x_1,x_2...x_s\) 表示每个图是否选择。那么包含这条边的所有图的 \(x\) 值异或起来显然为 \(0\) ,因为最后异或的结果不能包含这条边。
然后就得到了若干异或方程,形如:\(x_{a_1}\oplus x_{a_2}...\oplus x_{a_k}=0\),\(a_1,a_2...a_k\) 表示含这条边的图的编号。
用线性基求出它有多少解,就是有多少方案
线性基相当于 %2 下高斯消元
消完之后可以得到有多少主元,也就是 \(p[i]\neq 0\)
数出来假设有 \(c\) 个,剩下的 \(x\) 都可以随便取 (从而确定这 \(c\) 个,一一对应一个解)。
于是解数为 \(2^{s-c}\)
然后把它累加到 \(f\) 里就行了。当然也可以一边搜一边累加,不显示记录 \(f\)。
(我就这么写的)
代码:
#include <bits/stdc++.h>
using namespace std;
namespace Flandre_Scarlet
{
#define int long long
#define F(i,l,r) for(int i=l;i<=r;++i)
#define D(i,r,l) for(int i=r;i>=l;--i)
#define Fs(i,l,r,c) for(int i=l;i<=r;c)
#define Ds(i,r,l,c) for(int i=r;i>=l;c)
#define MEM(x,a) memset(x,a,sizeof(x))
#define FK(x) MEM(x,0)
#define Tra(i,u) for(int i=G.st(u),v=G.to(i);~i;i=G.nx(i),v=G.to(i))
#define p_b push_back
#define sz(a) ((int)a.size())
#define all(a) a.begin(),a.end()
#define iter(a,p) (a.begin()+p)
#define PUT(a,n) F(i,1,n) printf("%d ",a[i]); puts("");
int I() {char c=getchar(); int x=0; int f=1; while(c<'0' or c>'9') f=(c=='-')?-1:1,c=getchar(); while(c>='0' and c<='9') x=(x<<1)+(x<<3)+(c^48),c=getchar(); return ((f==1)?x:-x);}
template <typename T> void Rd(T& arg){arg=I();}
template <typename T,typename...Types> void Rd(T& arg,Types&...args){arg=I(); Rd(args...);}
void RA(int *p,int n) {F(i,1,n) *p=I(),++p;}
int vn,gn; // 点数, 图数
char g[100][500];
void Input()
{
Rd(gn);
int len;
F(i,1,gn)
{
scanf("%s",g[i]+1);
len=strlen(g[i]+1);
}
for(vn=1;vn<=10;++vn)
{
if (vn*(vn-1)/2==len) {break;}
}
}
int eid[20][20];
int m,col[20];
class Linear_basis
{
public:
int p[64];
void clear()
{
FK(p);
}
void ins(int x)
{
D(i,61,0) if ((x>>i)&1)
{
if (!p[i]) {p[i]=x; break;}
else x^=p[i];
}
}
int sakuya() // 求解数
{
int cnt=0;
F(i,0,61) if (p[i]) ++cnt;
return 1ll<<(gn-cnt);
}
}LB; // 线性基
int fc[20];
int ans=0;
void DFS(int cur,int m)
{
if (cur==vn+1)
{
LB.clear();
F(i,1,vn) F(j,i+1,vn) if (col[i]!=col[j]) // 钦定(i,j)边不能出现
{
int e=eid[i][j];
int cur=0; // cur 状压记录方程
F(k,1,gn) if (g[k][e]=='1') // 包含这条边
{
cur|=(1ll<<k); // 方程中算上 x_k
}
LB.ins(cur); // 新建一条异或方程
}
int cans=LB.sakuya()*fc[m-1];
if (m&1) ans+=cans;
else ans-=cans;
return;
}
F(i,1,min(m+1,vn)) // 小技巧:一边枚举颜色一边更新集合数,常数小,并且保证了每个划分都非空
{
col[cur]=i;
DFS(cur+1,max(m,i));
}
}
void Sakuya()
{
int tot=0;
F(i,1,vn) F(j,i+1,vn) eid[i][j]=++tot;
fc[0]=1; F(i,1,vn) fc[i]=fc[i-1]*i;
ans=0;
DFS(1,0);
printf("%lld\n",ans);
}
void IsMyWife()
{
Input();
Sakuya();
}
}
#undef int //long long
int main()
{
Flandre_Scarlet::IsMyWife();
getchar();
return 0;
}