莫比乌斯反演
简介¶
莫比乌斯反演是数论中的重要内容。对于一些函数 f(n) ,如果很难直接求出它的值,而容易求出其倍数和或约数和 g(n) ,那么可以通过莫比乌斯反演简化运算,求得 f(n) 的值。
开始学习莫比乌斯反演前,我们需要一些前置知识: 积性函数 、 Dirichlet 卷积 、 莫比乌斯函数 。
前置知识¶
引理 1¶
略证:
引理 2¶
|V| 表示集合 V 的元素个数
略证:
对于 d\leq \left\lfloor\sqrt{n}\right\rfloor , \left\lfloor\frac{n}{d}\right\rfloor 有 \left\lfloor\sqrt{n}\right\rfloor 种取值
对于 d> \left\lfloor\sqrt{n}\right\rfloor ,有 \left\lfloor\frac{n}{d}\right\rfloor\leq\left\lfloor\sqrt{n}\right\rfloor ,也只有 \left\lfloor\sqrt{n}\right\rfloor 种取值
综上,得证
数论分块¶
数论分块的过程大概如下:考虑含有 \left\lfloor\frac{n}{i}\right\rfloor 的求和式子( n 为常数)
对于任意一个 i(i\leq n) ,我们需要找到一个最大的 j(i\leq j\leq n) ,使得 \left\lfloor\frac{n}{i}\right\rfloor = \left\lfloor\frac{n}{j}\right\rfloor .
而 j=\left\lfloor\frac{n}{\left\lfloor\frac{n}{i}\right\rfloor}\right\rfloor .
略证:
即 j=\left\lfloor\frac{n}{\left\lfloor\frac{n}{i}\right\rfloor}\right\rfloor .
利用上述结论,我们每次以 [i,j] 为一块,分块求和即可
积性函数¶
定义¶
若 \gcd(x,y)=1 且 f(xy)=f(x)f(y) ,则 f(n) 为积性函数。
性质¶
若 f(x) 和 g(x) 均为积性函数,则以下函数也为积性函数:
例子¶
- 单位函数: \epsilon(n)=[n=1]
- 恒等函数: \operatorname{id}_k(n)=n^k \operatorname{id}_{1}(n) 通常简记作 \operatorname{id}(n) 。
- 常数函数: 1(n)=1
- 除数函数: \sigma_{k}(n)=\sum_{d\mid n}d^{k} \sigma_{0}(n) 通常简记作 \operatorname{d}(n) 或 \tau(n) , \sigma_{1}(n) 通常简记作 \sigma(n) 。
- 欧拉函数: \varphi(n)=\sum_{i=1}^n [\gcd(i,n)=1]
- 莫比乌斯函数: \mu(n) = \begin{cases}1 & n=1 \\ 0 & \exists d>1:d^{2} \mid n \\ (-1)^{\omega(n)} & otherwise\end{cases} 其中 \omega(n) 表示 n 的本质不同质因子个数,是一个积性函数。
Dirichlet 卷积¶
定义¶
定义两个数论函数 f,g 的 Dirichlet 卷积为
性质¶
Dirichlet 卷积满足交换律和结合律。
其中 \varepsilon 为 Dirichlet 卷积的单位元(任何函数卷 \varepsilon 都为其本身)
例子¶
莫比乌斯函数¶
定义¶
\mu 为莫比乌斯函数,定义为
详细解释一下:
令 n=\prod_{i=1}^kp_i^{c_i} ,其中 p_i 为质因子, c_i\ge 1 。上述定义表示:
- n=1 时, \mu(n)=1 ;
- 对于 n\not= 1 时:
- 当存在 i\in [1,k] ,使得 c_i > 1 时, \mu(n)=0 ,也就是说只要某个质因子出现的次数超过一次, \mu(n) 就等于 0 ;
- 当任意 i\in[1,k] ,都有 c_i=1 时, \mu(n)=(-1)^k ,也就是说每个质因子都仅仅只出现过一次时,即 n=\prod_{i=1}^kp_i , \{p_i\}_{i=1}^k 中个元素唯一时, \mu(n) 等于 -1 的 k 次幂,此处 k 指的便是仅仅只出现过一次的质因子的总个数。
性质¶
莫比乌斯函数不但是积性函数,还有如下性质:
即 \sum_{d\mid n}\mu(d)=\varepsilon(n) ,即 \mu * 1 =\varepsilon
证明¶
设 \displaystyle n=\prod_{i=1}^k{p_i}^{c_i},n'=\prod_{i=1}^k p_i
那么 \displaystyle\sum_{d\mid n}\mu(d)=\sum_{d\mid n'}\mu(d)=\sum_{i=0}^k C_k^i\cdot(-1)^i=(1+(-1))^k
根据二项式定理,易知该式子的值在 k=0 即 n=1 时值为 1 否则为 0 ,这也同时证明了 \displaystyle\sum_{d\mid n}\mu(d)=[n=1]=\varepsilon(n) 以及 \mu\ast 1=\varepsilon
补充结论¶
反演结论: \displaystyle [\gcd(i,j)=1] \iff\sum_{d\mid\gcd(i,j)}\mu(d)
直接推导 :如果看懂了上一个结论,这个结论稍加思考便可以推出:如果 \gcd(i,j)=1 的话,那么代表着我们按上个结论中枚举的那个 n 是 1 ,也就是式子的值是 1 ,反之,有一个与 [\gcd(i,j)=1] 相同的值: 0
利用 \varepsilon 函数 :根据上一结论, [\gcd(i,j)=1]\implies \varepsilon(\gcd(i,j)) ,将 \varepsilon 展开即可。
线性筛¶
由于 \mu 函数为积性函数,因此可以线性筛莫比乌斯函数(线性筛基本可以求所有的积性函数,尽管方法不尽相同)。
线性筛实现
void getMu() {
mu[1] = 1;
for (int i = 2; i <= n; ++i) {
if (!flg[i]) p[++tot] = i, mu[i] = -1;
for (int j = 1; j <= tot && i * p[j] <= n; ++j) {
flg[i * p[j]] = 1;
if (i % p[j] == 0) {
mu[i * p[j]] = 0;
break;
}
mu[i * p[j]] = -mu[i];
}
}
}
拓展¶
证明
将 n 分解质因数: \displaystyle n=\prod_{i=1}^k {p_i}^{c_i}
首先,因为 \varphi 是积性函数,故只要证明当 n'=p^c 时 \displaystyle\varphi \ast 1=\sum_{d\mid n'}\varphi(\frac{n'}{d})=\text{ID} 成立即可。
因为 p 是质数,于是 d=p^0,p^1,p^2,\cdots,p^c
易知如下过程:
该式子两侧同时卷 \mu 可得 \displaystyle\varphi(n)=\sum_{d\mid n}d\cdot\mu(\frac{n}{d})
莫比乌斯反演¶
公式¶
设 f(n),g(n) 为两个数论函数。
如果有
那么有
证明¶
方法一:对原式做数论变换。
用 \displaystyle\sum_{d\mid n}g(d) 来替换 f(\dfrac{n}{d}) ,再变换求和顺序。最后一步变换的依据: \displaystyle\sum_{d\mid n}\mu(d)=[n=1] ,因此在 \dfrac{n}{k}=1 时第二个和式的值才为 1 。此时 n=k ,故原式等价于 \displaystyle\sum_{k\mid n}[n=k]\cdot g(k)=g(n)
方法二:运用卷积。
原问题为:已知 f=g\ast1 ,证明 g=f\ast\mu
易知如下转化: f\ast\mu=g*1*\mu\implies f\ast\mu=g (其中 1\ast\mu=\varepsilon )。
问题形式¶
「HAOI 2011」Problem b¶
求值(多组数据)
根据容斥原理,原式可以分成 4 块来处理,每一块的式子都为
考虑化简该式子
因为 \gcd(i,j)=1 时对答案才用贡献,于是我们可以将其替换为 \varepsilon(\gcd(i,j)) ( \varepsilon(n) 当且仅当 n=1 时值为 1 否则为 0 ),故原式化为
将 \varepsilon 函数展开得到
变换求和顺序,先枚举 d\mid \gcd(i,j) 可得
易知 1\sim\lfloor\dfrac{n}{k}\rfloor 中 d 的倍数有 \lfloor\dfrac{n}{kd}\rfloor 个,故原式化为
很显然,式子可以数论分块求解。
时间复杂度 \Theta(N+T\sqrt{n})
代码实现
#include <algorithm>
#include <cstdio>
const int N = 50000;
int mu[N + 5], p[N + 5];
bool flg[N + 5];
void init() {
int tot = 0;
mu[1] = 1;
for (int i = 2; i <= N; ++i) {
if (!flg[i]) {
p[++tot] = i;
mu[i] = -1;
}
for (int j = 1; j <= tot && i * p[j] <= N; ++j) {
flg[i * p[j]] = 1;
if (i % p[j] == 0) {
mu[i * p[j]] = 0;
break;
}
mu[i * p[j]] = -mu[i];
}
}
for (int i = 1; i <= N; ++i) mu[i] += mu[i - 1];
}
int solve(int n, int m) {
int res = 0;
for (int i = 1, j; i <= std::min(n, m); i = j + 1) {
j = std::min(n / (n / i), m / (m / i));
res += (mu[j] - mu[i - 1]) * (n / i) * (m / i);
}
return res;
}
int main() {
int T, a, b, c, d, k;
init();
for (scanf("%d", &T); T; --T) {
scanf("%d%d%d%d%d", &a, &b, &c, &d, &k);
printf("%d\n", solve(b / k, d / k) - solve(b / k, (c - 1) / k) -
solve((a - 1) / k, d / k) +
solve((a - 1) / k, (c - 1) / k));
}
return 0;
}
「SPOJ 5971」LCMSUM¶
求值(多组数据)
易得原式即
将原式复制一份并且颠倒顺序,然后将 n 一项单独提出,可得
根据 \gcd(i,n)=\gcd(n-i,n) ,可将原式化为
两个求和式中分母相同的项可以合并。
即
可以将相同的 \gcd(i,n) 合并在一起计算,故只需要统计 \gcd(i,n)=d 的个数。当 \gcd(i,n)=d 时, \displaystyle\gcd(\frac{i}{d},\frac{n}{d})=1 ,所以 \gcd(i,n)=d 的个数有 \displaystyle\varphi(\frac{n}{d}) 个。
故答案为
变换求和顺序,设 \displaystyle d'=\frac{n}{d} ,合并公因式,式子化为
设 \displaystyle \text{g}(n)=\sum_{d\mid n} d\cdot\varphi(d) ,已知 \text{g} 为积性函数,于是可以 \Theta(n) 筛出。每次询问 \Theta(1) 计算即可。
这个函数筛的时候比较特殊,当 p_j\mid i 的时候,需要根据 p_j\mid\dfrac{i}{p_j} 进行分类讨论。具体可以见代码。
时间复杂度 : \Theta(n+T)
代码实现
#include <cstdio>
const int N = 1000000;
int tot, p[N + 5];
long long g[N + 5];
bool flg[N + 5];
void solve() {
g[1] = 1;
for (int i = 2; i <= N; ++i) {
if (!flg[i]) p[++tot] = i, g[i] = i * (i - 1) + 1;
for (int j = 1; j <= tot && i * p[j] <= N; ++j) {
flg[i * p[j]] = 1;
if (i % p[j] == 0) {
if ((i / p[j]) % p[j] == 0) {
g[i * p[j]] = g[i] + (g[i] - g[i / p[j]]) * p[j] * p[j];
} else {
g[i * p[j]] = g[i] + g[i / p[j]] * (p[j] - 1) * p[j] * p[j] * p[j];
}
break;
}
g[i * p[j]] = g[i] * g[p[j]];
}
}
}
int main() {
int T, n;
solve();
for (scanf("%d", &T); T; --T) {
scanf("%d", &n);
printf("%lld\n", (g[n] + 1) * n / 2);
}
return 0;
}
「BZOJ 2154」Crash 的数字表格¶
求值(对 20101009 取模)
易知原式等价于
枚举最大公因数 d ,显然两个数除以 d 得到的数互质
非常经典的 \gcd 式子的化法
后半段式子中,出现了互质数对之积的和,为了让式子更简洁就把它拿出来单独计算。于是我们记
接下来对 \text{sum}(n,m) 进行化简。首先枚举约数,并将 [\gcd(i,j)=1] 表示为 \varepsilon(\gcd(i,j))
设 i=i'\cdot d , j=j'\cdot d ,显然式子可以变为
观察上式,前半段可以预处理前缀和;后半段又是一个范围内数对之和,记
可以 \Theta(1) 求解
至此
我们可以 \lfloor\frac{n}{\lfloor\frac{n}{d}\rfloor}\rfloor 数论分块求解 \text{sum}(n,m) 函数。
在求出 \text{sum}(n,m) 后,回到定义 \text{sum} 的地方,可得原式为
可见这又是一个可以数论分块求解的式子!
本题除了推式子比较复杂、代码细节较多之外,是一道很好的莫比乌斯反演练习题!(上述过程中,默认 n\leqslant m )
时间复杂度: \Theta(n+m) (瓶颈为线性筛)
代码实现
#include <algorithm>
#include <cstdio>
using std::min;
const int N = 1e7;
const int mod = 20101009;
int n, m, mu[N + 5], p[N / 10 + 5], sum[N + 5];
bool flg[N + 5];
void init() {
mu[1] = 1;
int tot = 0, k = min(n, m);
for (int i = 2; i <= k; ++i) {
if (!flg[i]) p[++tot] = i, mu[i] = -1;
for (int j = 1; j <= tot && i * p[j] <= k; ++j) {
flg[i * p[j]] = 1;
if (i % p[j] == 0) {
mu[i * p[j]] = 0;
break;
}
mu[i * p[j]] = -mu[i];
}
}
for (int i = 1; i <= k; ++i)
sum[i] = (sum[i - 1] + 1LL * i * i % mod * (mu[i] + mod)) % mod;
}
int Sum(int x, int y) {
return (1LL * x * (x + 1) / 2 % mod) * (1LL * y * (y + 1) / 2 % mod) % mod;
}
int func(int x, int y) {
int res = 0;
for (int i = 1, j; i <= min(x, y); i = j + 1) {
j = min(x / (x / i), y / (y / i));
res = (res + 1LL * (sum[j] - sum[i - 1] + mod) * Sum(x / i, y / i) % mod) %
mod;
}
return res;
}
int solve(int x, int y) {
int res = 0;
for (int i = 1, j; i <= min(x, y); i = j + 1) {
j = min(x / (x / i), y / (y / i));
res = (res +
1LL * (j - i + 1) * (i + j) / 2 % mod * func(x / i, y / i) % mod) %
mod;
}
return res;
}
int main() {
scanf("%d%d", &n, &m);
init();
printf("%d\n", solve(n, m));
}
「SDOI2015」约数个数和¶
多组数据,求
其中 d(n) 表示 n 的约数个数
要推这道题首先要了解 d 函数的一个特殊性质
再化一下这个式子
将上述式子代回原式
那么 O(n) 预处理 \mu,d 的前缀和, O(\sqrt{n}) 分块处理询问,总复杂度 O(n+T\sqrt{n}) .
代码实现
#include <algorithm>
#include <cstdio>
#define int long long
using namespace std;
const int N = 5e4 + 5;
int n, m, T, pr[N], mu[N], d[N], t[N], cnt; // t 表示 i 的最小质因子出现的次数
bool bp[N];
void prime_work(int k) {
bp[0] = bp[1] = 1, mu[1] = 1, d[1] = 1;
for (int i = 2; i <= k; i++) {
if (!bp[i]) pr[++cnt] = i, mu[i] = -1, d[i] = 2, t[i] = 1;
for (int j = 1; j <= cnt && i * pr[j] <= k; j++) {
bp[i * pr[j]] = 1;
if (i % pr[j] == 0) {
mu[i * pr[j]] = 0, d[i * pr[j]] = d[i] / (t[i] + 1) * (t[i] + 2),
t[i * pr[j]] = t[i] + 1;
break;
} else
mu[i * pr[j]] = -mu[i], d[i * pr[j]] = d[i] << 1, t[i * pr[j]] = 1;
}
}
for (int i = 2; i <= k; i++) mu[i] += mu[i - 1], d[i] += d[i - 1];
}
int solve() {
int res = 0, mxi = min(n, m);
for (int i = 1, j; i <= mxi; i = j + 1)
j = min(n / (n / i), m / (m / i)),
res += d[n / i] * d[m / i] * (mu[j] - mu[i - 1]);
return res;
}
signed main() {
scanf("%lld", &T);
prime_work(50000);
while (T--) {
scanf("%lld%lld", &n, &m);
printf("%lld\n", solve());
}
return 0;
}
「luogu 3768」简单的数学题¶
求
看似是一道和 \gcd 有关的题,不过由于带有系数,并不容易化简
我们利用 \varphi\ast1=ID 反演
对 \sum_{d=1}^nF\left(\left\lfloor\frac{n}{d}\right\rfloor\right)^2 做数论分块, d^2\varphi(d) 的前缀和用杜教筛处理:
杜教筛(见 杜教筛 - 例 3 )完了是这样的
分块递归求解即可,复杂度 O(n^{\frac{2}{3}}) .
代码实现
#include <cmath>
#include <cstdio>
#include <map>
#define int long long
using namespace std;
const signed N = 5e6, NP = 5e6, SZ = N;
int n, P, inv2, inv6, s[N];
signed phi[N], p[NP], cnt, pn;
bool bp[N];
map<int, int> s_map;
int ksm(int a, int m) { // 求逆元用
int res = 1;
while (m) {
if (m & 1) res = res * a % P;
a = a * a % P, m >>= 1;
}
return res;
}
void prime_work(signed k) { // 线性筛phi,s
bp[0] = bp[1] = 1, phi[1] = 1;
for (signed i = 2; i <= k; i++) {
if (!bp[i]) p[++cnt] = i, phi[i] = i - 1;
for (signed j = 1; j <= cnt && i * p[j] <= k; j++) {
bp[i * p[j]] = 1;
if (i % p[j] == 0) {
phi[i * p[j]] = phi[i] * p[j];
break;
} else
phi[i * p[j]] = phi[i] * phi[p[j]];
}
}
for (signed i = 1; i <= k; i++)
s[i] = (1ll * i * i % P * phi[i] % P + s[i - 1]) % P;
}
int s3(int k) {
return k %= P, (k * (k + 1) / 2) % P * ((k * (k + 1) / 2) % P) % P;
} // 立方和
int s2(int k) {
return k %= P, k * (k + 1) % P * (k * 2 + 1) % P * inv6 % P;
} // 平方和
int calc(int k) { // 计算S(k)
if (k <= pn) return s[k];
if (s_map[k]) return s_map[k]; // 对于超过pn的用map离散存储
int res = s3(k), pre = 1, cur;
for (int i = 2, j; i <= k; i = j + 1)
j = k / (k / i), cur = s2(j),
res = (res - calc(k / i) * (cur - pre) % P) % P, pre = cur;
return s_map[k] = (res + P) % P;
}
int solve() {
int res = 0, pre = 0, cur;
for (int i = 1, j; i <= n; i = j + 1)
j = n / (n / i), cur = calc(j),
res = (res + (s3(n / i) * (cur - pre)) % P) % P, pre = cur;
return (res + P) % P;
}
signed main() {
scanf("%lld%lld", &P, &n);
inv2 = ksm(2, P - 2), inv6 = ksm(6, P - 2);
pn = (int)pow(n, 0.666667); // n^(2/3)
prime_work(pn);
printf("%lld", solve());
return 0;
} // 不要为了省什么内存把数组开小。。。卡了好几次80
解法二
转化一下,可以将式子写成
容易知道
设 sum(n,m)=\sum_{i=1}^{n}\sum_{j=1}^{m}ij ,继续接着前面的往下推
这时我们只要对每个 T 预处理出 T\sum_{d\mid T}d\cdot\mu(d) 的值就行了,考虑如何快速求解
设 f(n)=\sum_{d\mid n}d\cdot\mu(d)
实际上 f 可以用线性筛筛出,具体的是
其中 p 表示 n 的最小质因子,总时间复杂度 O(n+\sqrt n) 。
莫比乌斯反演扩展¶
结尾补充一个莫比乌斯反演的非卷积形式。
对于数论函数 f,g 和完全积性函数 t 且 t(1)=1 :
我们证明一下
参考文献¶
build本页面最近更新:,更新历史
edit发现错误?想一起完善? 在 GitHub 上编辑此页!
people本页面贡献者:hydingsy, hyp1231
copyright本页面的全部内容在 CC BY-SA 4.0 和 SATA 协议之条款下提供,附加条款亦可能应用