目录
写在前面
- 这是继数论和组合计数类数学相关与多项式类数学相关后的第三篇数学方面内容总结。主要记录自己近期学习的一些数学方法。内容比较杂,同时也起到对前文的一些补充作用。
- 若文章中出现错误,烦请告知。
- 感谢您的造访。
一类反演问题
莫比乌斯反演
这一个内容在数论和组合计数类数学相关里面提到了。
快速莫比乌斯变换(反演)与子集卷积
莫比乌斯变换(反演)
问题提出
若 (h,f,g) 为下标为集合的函数,我们定义
[h=f*g]
表示
[h(S) = sum_{L subseteq S}^{} sum_{R subseteq S}^{} [L cup R = S] f(L) imes g(R)]
容易发现,对于这个问题,我们可以用 (Oleft((2^n)^2 ight)) 的枚举 (L,R) 来计算。
然而这样复杂度较高,我们考虑类比多项式卷积的过程,可以求出 (f,g) 的点值,直接相乘得到 (h) 的点值然后再插回去。
值得注意的是为了便于表述以及规范表达,快速莫比乌斯变换就相当于点值,快速莫比乌斯反演就相当于插值。
算法原理
- 我们定义 (f) 的莫比乌斯变换为 (F) ,其中 (F(S)=sum_{Xsubseteq S}f(X)) ;由这个定义,我们可以推出 (F) 莫比乌斯反演 (f) 为 (f(S) = sum_{X subseteq S} (-1)^{|S| - |X|}F(X)) 。对于莫比乌斯反演的证明,可以带入莫比乌斯变换的式子或容斥来证。
- 我们对于一个函数 (f,g,h) ,记它的点值式为 (F,G,H) 。我们将问题提出中的卷积式两边同时做莫比乌斯变换,得到 [egin{aligned}H(S) &= sum_{L subseteq S} sum_{R subseteq S} [L cup R subseteq S] f(L) imes g(R)\&=sum_{L subseteq S} sum_{R subseteq S}f(L) imes g(R)\&=left(sum_{L subseteq S}f(L) ight) imesleft(sum_{R subseteq S}g(R) ight)\&=F(S) imes G(S)end{aligned}]
至此算法原理及过程已经完全结束。似乎我们可以用 (O(3^n)) 枚举子集来变换和反演,实际上我们可以让复杂度更优。
算法实现
- 设 (hat{f_{S}}^{(i)}) 表示 (sum_{T subseteq S} [(S - T) subseteq {1, 2, ..., i}] f_{T})
- 易得初始状态:(hat{f_{S}}^{(0)} = f_{S})
- 对于每一个不包含 ({i}) 的集合 (S),可知 (hat{f_{S}}^{(i)} = hat{f_{S}}^{(i - 1)})(因为 (S) 并没有 (i) 这位),(hat{f}_{S cup {i}}^{(i)} = hat{f}_{S}^{(i - 1)} + hat{f}_{S cup {i}}^{(i - 1)})(前者的 (T) 没有包含 ({i}),而后者的 (T) 必须包含了 ({i}))。
- 显然,递推了 (n) 轮之后,(hat{f}_{S}^n) 就是所求的变换了。
用高维前缀和可以做到 (O(n imes 2^n)) 的递推,求出点值和插值。
void FMT(int *A, int o) {// o 为识别因子
for (int i = 1; i < ST; i <<= 1)//ST-1 表示全集
for (int j = 0; j < ST; j++)
if (i&j) (A[j] += A[j^i]*o) %= mod;
}
例题
子集卷积
( ext{FWT}) :“你刚才说的那个玩意我也能做啊,要你何用?”
( ext{FMT}) :“……”
问题提出
若 (h,f,g) 为下标为集合的函数,我们定义
[h=f*g]
表示
[h(S) = sum_{X subseteq S} f(X) imes g(S-X)]
算法实现
回顾刚刚的集合并卷积,子集卷积的条件比集合并卷积更苛刻,即 (L) 和 (R) 的集合应该不相交。
我们可以在卷积时多加一维,维护集合的大小,如 (f_{i,S}) 表示集合中有 (i) 个元素,集合表示为 (S) 。显然,当 (i) 和 (S) 的真实元素个数符合时才是对的。记数组 cnt[S]
表示集合 (S) 的模。初始时,我们只把 (f_{cnt[S],S}) 的值赋成原来的 (f(S)) ( (g) 同理),然后每一维做一遍 ( ext{FMT}) ,点值相乘时这么写:(h_{i, S} = sum_{j = 0}^{i} f_{j,S} imes g_{i - j, S}) 。最后扫一遍把不符合实际情况的状态赋成 (0) 即可。
for (int i = 0; i <= n; i++) FMT(g[i], 1);
for (int i = 0; i <= n; i++) FMT(f[i], 1);
for (int i = 1; i <= n; i++) {
for (int j = 0; j <= i; j++)
for (int k = 0; k < ST; k++)
(h[i][k] += 1ll*f[j][k]*g[i-j][k]%mod) %= mod;
FMT(h[i], -1);
for (int k = 0; k < ST; k++) if (cnt[k] != i) h[i][k] = 0;
if (i != n) FMT(h[i], 1);
}
例题
二项式反演
内容
对于函数 (f,g) , (forall pinmathbb{N}) 若 (forall ngeq p) ,满足
[f(n)=sum_{k=p}^{n}inom{n}{k}g(k)]
那么 (forall ngeq p)
[g(n)=sum_{k=p}^{n}(-1)^{n-k}inom{n}{k}f(k)]
证明
为了方便表达,我们取 (p=0) ,实质和取 (pinmathbb{N}) 的证明方法是一样的。
[egin{aligned}g(n)&=sum_{k=0}^{n}(-1)^{n-k}inom{n}{k}f(k)\&=sum_{k=0}^{n}(-1)^{n-k}inom{n}{k}sum_{i=0}^k{kchoose i}g(i)\&=sum_{k=0}^nsum_{i=0}^k(-1)^{n-k}{nchoose k}{kchoose i}g(i)\&=sum_{i=0}^nleft(sum_{k=i}^n(-1)^{n-k}{nchoose k}{kchoose i} ight)g(i)\&=sum_{i=0}^nleft(sum_{k=i}^n(-1)^{n-k}{nchoose i}{n-ichoose k-i} ight)g(i)\&=sum_{i=0}^nleft({nchoose i}sum_{k=i}^n(-1)^{n-k}{n-ichoose n-k} ight)g(i)\&=sum_{i=0}^nleft({nchoose i}(1-1)^{n-i} ight)g(i)\&=g(n)end{aligned}]
故成立。
应用举例
推导错排公式
我们记 (f(n)) 为 (n) 个数字任意放的方案数, (g(n)) 为 (n) 个数没有一个放在自己位置上的方案数。
枚举不在自己位置上的个数,容易得到
[f(n)=sum_{i=0}^n{nchoose i}g(i)]
那么
[egin{aligned}g(n)&=sum_{i=0}^n(-1)^{n-i}{nchoose i}f(i)\&=sum_{i=0}^n(-1)^i{nchoose i}f(n-i)end{aligned}]
注意到 (f(x)=x!) ,那么
[egin{aligned}g(n)&=sum_{i=0}^n(-1)^ifrac{n!}{i!(n-i)!}(n-i)!\&=n!sum_{i=0}^n(-1)^ifrac{1}{i!}end{aligned}]
棋盘染色
有个 (1 imes n) 的格子, (m) 种颜色( (mgeq 2) ),要求相邻格子的颜色不相同且每种颜色都要用到,求染色方案数。
我们记 (f(n)) 为至多用到 (n) 种颜色的方案数, (g(n)) 为 (n) 为恰用到 (n) 种颜色的方案数。
那么
[egin{aligned}f(m)&=sum_{i=2}^m{mchoose i}g(i)\Rightarrow g(m)&=sum_{i=2}^m(-1)^i{mchoose i}f(n-i)end{aligned}]
注意到 (f(x)=x imes(x-1)^{n-1}) 。那么就可以带入直接算了。
另一形式
[a_k=sumlimits_{i=k}^n{ichoose k}b_iRightarrow b_k=sumlimits_{i=k}^n(-1)^{i-k}{ichoose k}a_i]
证明:
[egin{aligned} &sumlimits_{i=k}^n(-1)^{i-k}{ichoose k}a_i\=& sumlimits_{i=k}^n(-1)^{i-k}{ichoose k}sumlimits_{j=k}^n{jchoose i}b_i\=& sumlimits_{i=k}^nsumlimits_{j=k}^n(-1)^{i-k}{ichoose k}{jchoose i}b_i\=& sumlimits_{i=k}^nsumlimits_{j=k}^n(-1)^{i-k}{jchoose i}{ichoose k}b_i\=& sumlimits_{j=k}^nsumlimits_{i=k}^j(-1)^{i-k}{jchoose k}{j-kchoose i-k}b_i\=& sumlimits_{j=k}^n{jchoose k}sumlimits_{i=k}^j(-1)^{i-k}{j-kchoose i-k}b_i\=&sumlimits_{j=k}^n{jchoose k}(1-1)^{j-k}b_i\=&b_kend{aligned}]
例题
斯特林反演
第一类斯特林数
(egin{bmatrix}n\mend{bmatrix}) 表示将 (n) 个元素排成 (m) 个轮换的方法数。
递推公式:(egin{bmatrix}n\mend{bmatrix}=egin{bmatrix}n-1\m-1end{bmatrix}+(n-1)egin{bmatrix}n-1\mend{bmatrix})
含义是考虑第 (n) 个元素的放法:要么新开一个轮换,要么就放在前 (n-1) 个元素的左边。
第二类斯特林数
(egin{Bmatrix}n\mend{Bmatrix}) 表示将 (n) 个元素划分成 (m) 个非空子集的方法数。
递推公式:(egin{Bmatrix}n\mend{Bmatrix}=egin{Bmatrix}n-1\m-1end{Bmatrix}+megin{Bmatrix}n-1\mend{Bmatrix})
含义是考虑第 (n) 个元素的放法:要么新开一个组,要么就放在前 (m) 组内。
通项公式(容斥式): (egin{Bmatrix}n\mend{Bmatrix}=frac{1}{m!}sumlimits_{k=0}^{m}(-1)^kinom{m}{k}(m-k)^n)
有关通项公式的证明及运用可以参考多项式类数学相关这篇文章。
例题
反演公式
[f(x) = sum_{i=0}^x egin{Bmatrix}x\iend{Bmatrix} g(i) Leftrightarrow g(x) = sum_{i=0}^x (-1) ^ {x - i}egin{bmatrix}x\iend{bmatrix} f(i)]
[f(x) = sum_{i=0}^x egin{bmatrix}x\iend{bmatrix} g(i) Leftrightarrow g(x) = sum_{i=0}^x (-1) ^ {x - i}egin{Bmatrix}x\iend{Bmatrix} f(i)]
例题
- 给出 (n) 个点的一张简单图,问有多少个边的子集,满足保留子集中的边后,该图连通。(蒯自Sdchr)
- 大概就是枚举连通块的个数,然后块内随便连,然后容斥就好。
- 考虑如何求容斥系数 (f(i)) 。设实际上是 (x) 个连通块的方案,它应该被计算 ([x=1]) 次,实际上在所有更仔细的分块中被统计,所以
- [[x = 1] = sum_{i=1}^xegin{Bmatrix}x\iend{Bmatrix}f(i)]
- 由斯特林反演
- [egin{aligned}f(x) &= sum_{i=1}^x(-1)^{x-i}egin{bmatrix}x\iend{bmatrix}[i=1]\&=(-1)^{x-1}egin{bmatrix}x\1end{bmatrix}\&=(-1)^{x-1}(x-1)!end{aligned}]
- [BZOJ 4671]异或图
最值反演( ( ext{min-max}) 容斥)
公式
记 (max(S)) 为集合 (S) 中的最大值, (min(S)) 为集合 (S) 中的最小值, (|S|) 为集合 (S) 的元素数量,那么以下两个等式成立
[max(S)=sum_{Tsubseteq S}(-1)^{|T|+1}min(T)]
[min(S)=sum_{Tsubseteq S}(-1)^{|T|+1}max(T)]
证明
这里只证明第一个等式好了,后边的可以自行推出。
其实只需要证明一件事,就是除了 (min(T)=max(S)) 的那个值,其他的 (min) 值都被消掉了就可以了(这里说明一下,我们假定集合中的元素两两相异)
先来说明 (max(S)) 的系数为什么是 (1) ,假设中 (S) 最大的元素是 (a) ,那么我们会发现只有 (min({a})=max(S)) 所以 (max(S)) 的系数必须是 (1) 。
然后再说明为什么别的 (min) 都被消掉了,假设某个元素 (b) 的排名是 (k) ,那么 (min(T)=b) 当且仅当我们选出的集合是后 (n-k) 个的元素构成的集合的子集然后并上 ({b}) 得到的,我们会发现显然这样的集合有 (2^{n-k}) 种,而显然这其中恰有 (2^{n-k-1}) 中是有奇数个元素的,恰有 (2^{n-k-1}) 种是有偶数个元素的,两两相消自然就成 (0) 了,当然上述等式在 (k=n) 的时候不成立,但是此时剩下的刚好是最大值,所以证明完毕。
拉格朗日插值法
简介
给定 (n + 1) 个横坐标不相同的点,可以唯一确定一个 (n) 次的多项式。最直观的求多项式的做法就是列方程求解。但是这样需要 (O(n^3)) 的时间来计算。而拉格朗日插值法则通过构造的方法,得到了一个经过 (n + 1) 个点的 (n) 次多项式。 具体的过程是这样的,假设现在我们得到了 (n + 1) 个点:
[ (x_0, y_0), ;(x_1, y_1),; dots,;(x_n, y_n) ]
设拉格朗日基本多项式为
[ ell_j(x) = prod_{i = 0, i eq j}^n {x - x_i over x_j - x_i} ]
这个基本多项式构造十分巧妙,因为注意到 (ell_j(x_j) = 1) ,并且 (ell_j(x_i) = 0, ;forall i eq j) 。那么,接着构造出这个 (n) 次多项式
[ P(x) = sum_{i = 0}^n y_iell_i(x) ]
根据基本多项式的性质,我们可以知道 (P(x_i) = y_i) ,也就是经过了这 (n + 1) 个点。 通过简单的多项式乘法和多项式除法就可以在 (O(n^2)) 的时间求出这个多项式的系数表达。
求解
已知 (n) 次多项式 (f(n)) 上的 (n+1) 个点 ((x_i,y_i),iin[0,n]) ,求 (f(xi))
int lagrange(int n, int *x, int *y, int xi) {
int ans = 0;
for (int i = 0; i <= n; i++) {
int s1 = 1, s2 = 1;
for (int j = 0; j <= n; j++)
if (i != j) {
s1 = 1ll*s1*(xi-x[j])%mod;
s2 = 1ll*s2*(x[i]-x[j])%mod;
}
ans = (1ll*ans+1ll*y[i]*s1%mod*quick_pow(s2, mod-2)%mod)%mod;
}
return (ans+mod)%mod;
}
如果 (x) 的取值是连续一段的话,我们可以做到 (O(n)) 求解。假设 (forall i<j,x_i<x_j) (具体公式推导的话,如果你有兴趣可以参看之后的内容。因为比较显然,这里不再讲解。)
int lagrange(int n, int *x, int *y, int xi) {
int ans = 0;
s1[0] = (xi-x[0])%mod, s2[n+1] = 1;
for (int i = 1; i <= n; i++) s1[i] = 1ll*s1[i-1]*(xi-x[i])%mod;
for (int i = n; i >= 0; i--) s2[i] = 1ll*s2[i+1]*(xi-x[i])%mod;
ifac[0] = ifac[1] = 1;
for (int i = 2; i <= n; i++) ifac[i] = -1ll*mod/i*ifac[mod%i]%mod;
for (int i = 2; i <= n; i++) ifac[i] = 1ll*ifac[i]*ifac[i-1]%mod;
for (int i = 0; i <= n; i++)
(ans += 1ll*y[i]*(i == 0 ? 1 : s1[i-1])%mod*s2[i+1]%mod
*ifac[i]%mod*(((n-i)&1) ? -1 : 1)*ifac[n-i]%mod) %= mod;
return (ans+mod)%mod;
}
例题
自然数的幂的前缀和
问题提出
给定的 (n) 和 (k) ,求
[sum_{i = 1}^n i^k ]
通常 (n) 比较大,而 (k) 只有几千或者几万。
问题解决
我们可以知道,对于上述式子,推导公式一定是是 (k + 1) 次多项式。对于证明的话,我们可以参考riteme的介绍。
考虑使用拉格朗日插值法来获得答案多项式。
首先如果我们得知了 (n = 0, 1, dots, k + 1) 处的答案 (f(n)) ,那么给定的 (n) 处的答案可以写成
[ egin{aligned} f(n) & = sum_{i = 0}^{k + 1} f(i) {(n - 0)(n - 1)cdots[n - (i - 1)][(n - (i + 1)]cdots[n - (k + 1)] over (i - 0)(i - 1)cdots[i - (i - 1)][(i - (i + 1)]cdots[i - (k + 1)]} \ & = sum_{i = 0}^{k + 1} f(i) {prod_{j = 0}^{i - 1} (n - j) prod_{j = i + 1}^{k + 1} (n - j) over i!(-1)^{k - i + 1}(k + 1 - i)!} \ & = sum_{i = 0}^{k + 1} (-1)^{k - i + 1}f(i) {prod_{j = 0}^{i - 1} (n - j) prod_{j = i + 1}^{k + 1} (n - j) over i!(k + 1 - i)!} end{aligned} ]
注意到后面的分式中,分子是一个前缀积乘以一个后缀积,而分母是两个阶乘。这些都可以在 (O(k)) 的时间内求出。 现在剩下的问题就是如何求出 (f(0), f(1), dots, f(k + 1)) 了。由于 (g(x) = x^k) 是个完全积性函数,所以我们可以通过欧拉筛法求出 (g) 函数前面的一些值。具体的就是对于质数采取直接快速幂,合数则拆出任意一个因子来算,通常是欧拉筛法中可以顺便求得的最小质因子。根据素数定理,素数大约有 (Oleft(frac{k}{ln k} ight)) 个。每次快速幂需要花费 (O(log k)) 的时间,因此总的时间复杂度可以估计为 (O(k)) ,是一个非常优秀的算法。上面的方法具有通用性,只要我们可以快速的求出某个 (k) 次多项式的前 (k + 1) 个值,那么剩下的部分可以使用拉格朗日插值法在 (O(k)) 的时间内完成计算。
代码实现
int lagrange(int k, int *f, int xi) {//k+2 个点对 (i, f[i]), 0 <= i <= k+1
int ans = 0; ++k;
s1[0] = xi, s2[k+1] = 1;
for (int i = 1; i <= k; i++) s1[i] = 1ll*s1[i-1]*(xi-i)%mod;
for (int i = k; i >= 0; i--) s2[i] = 1ll*s2[i+1]*(xi-i)%mod;
for (int i = 0; i <= k; i++)
(ans += 1ll*f[i]*(i == 0 ? 1 : s1[i-1])%mod*s2[i+1]%mod
*ifac[i]%mod*(((k-i)&1) ? -1 : 1)*ifac[k-i]%mod) %= mod;
return (ans+mod)%mod;
}
void pre() {//预处理出阶乘逆元、插值的 k+2 个点
f[1] = ifac[0] = ifac[1] = 1;
for (int i = 2; i <= k+1; i++) ifac[i] = -1ll*mod/i*ifac[mod%i]%mod;
for (int i = 2; i <= k+1; i++) ifac[i] = 1ll*ifac[i-1]*ifac[i]%mod;
memset(isprime, 1, sizeof(isprime));
for (int i = 2; i <= k+1; i++) {
if (isprime[i]) prime[++tot] = i, f[i] = quick_pow(i, k);
for (int j = 1; j <= tot && prime[j]*i <= k+1; j++) {
isprime[i*prime[j]] = 0;
f[i*prime[j]] = 1ll*f[i]*f[prime[j]]%mod;
if (i%prime[j] == 0) break;
}
}
for (int i = 1; i <= k+1; i++) f[i] = (f[i]+f[i-1])%mod;
}
void work() {
scanf("%d", &k); pre();
while (~scanf("%d", &n)) {
if (n <= k+1) printf("%d
", f[n]);
else printf("%d
", lagrange(k, f, n));
}
}
例题