【9】莫比乌斯反演学习笔记
Eason_cyx
·
2025-07-25 16:38:22
·
算法·理论
好像咕了很久。><
莫比乌斯函数:
\mu(x)=\begin{cases}1&x=1\\0&\exist k > 1 \ \text{s.t.} k^2|x\\-1(w)&x=p_1p_2\dots p_w\end{cases}
莫比乌斯函数是 积性函数 ,可以使用 线性筛 。
狄利克雷卷积
对于两个数论函数 f(x) 和 g(x) ,他们的 狄利克雷 卷积得到的结果 h(x) 为:
h(x)=\displaystyle\sum_{d|x}f(d)g\Big(\frac{x}{d}\Big)
可以简记成 f * g = h 。
常用狄利克雷卷积:
记 I(n)=n ,1(n)=1 ,\varphi(n)=\displaystyle\sum_{d=1}^{n}[\gcd(n,d)=1] ,则有 \varphi * 1 = I 。
记 \epsilon(n)=[n=1] ,则 \mu * 1 = \epsilon 。重点,后面大多依赖这个性质推式子!
莫比乌斯反演
若 f(x)=\displaystyle\sum_{d|x}F(x) ,则 F(x)=\displaystyle\sum_{d|x}\mu(d)f(x) 。这个式子被称作 莫比乌斯反演 。
整除分块
例子:给定 n ,求 \displaystyle\sum_{i=1}^{n}\Big\lfloor\dfrac{n}{i}\Big\rfloor 的值。n \le 10^{12}
首先暴力做是 O(n) 的。
然后不难注意到对于某些连续的 i ,\Big\lfloor\dfrac{n}{i}\Big\rfloor 都是一样的。
可以发现对于一个数 L 使方程 \Big\lfloor\dfrac{n}{L}\Big\rfloor=\Big\lfloor\dfrac{n}{R}\Big\rfloor 成立的最大 R 是 \Big\lfloor\dfrac{n}{\lfloor\frac{n}{L}\rfloor}\Big\rfloor ,那么可以写出如下代码:
for(int l = 1, r = 1;l <= n;l = r + 1) {
r = n / (n / l);
ans += (r - l + 1) * (n / l);
}
时间复杂度 O(\sqrt n) 。
题
P3455 [POI 2007] ZAP-Queries
下面记 T=n,n=a,m=b 。
不妨 n<m 。
&=\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{k}\rfloor}\sum_{d|i,d|j}\mu(d)\\
&= \sum_{d=1}^{\lfloor\frac{n}{k}\rfloor}\mu(d)\sum_{i=1}^{\lfloor\frac{n}{kw}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{kw}\rfloor}1\\
&=\sum_{d=1}^{\lfloor\frac{n}{k}\rfloor}\mu(d)\lfloor\frac{n}{kw}\rfloor\lfloor\frac{m}{kw}\rfloor\end{aligned}
线性筛出莫比乌斯函数,即可做到 O(V+T\sqrt n) ,可以通过。
代码:(高亮处为整除分块)
#include <bits/stdc++.h>
using namespace std;
const int N = 5e4 + 5;
int mu[N]; bool isp[N];
vector<int> pr;
int solve(int n, int m) {
int ans = 0;
for(int i = 1, j = 1;i <= min(n, m);i = j + 1) {
j = min(n / (n / i), m / (m / i)); ans += (mu[j] - mu[i-1]) * (n / i) * (m / i);
} return ans;
} int main() {
mu[1] = 1;
for(int i = 2;i <= 50000;i++) {
if(!isp[i]) pr.push_back(i), mu[i] = -1;
for(auto j : pr) {
if((i * j) > 50000) break;
isp[i * j] = true;
if(i % j == 0) {
mu[i * j] = 0; break;
} else mu[i * j] = -mu[i];
}
} for(int i = 2;i <= 50000;i++) mu[i] += mu[i-1];
int n; cin >> n; while(n--) {
int a, b, k; cin >> a >> b >> k;
cout << solve(a / k, b / k) << endl;
}
return 0;
}
P1829 [国家集训队] Crash的数字表格 / JZPTAB
不妨 n < m 。
&=\sum_{i=1}^{n}\sum_{j=1}^{m}\dfrac{ij}{\gcd(i,j)}\\
&=\sum_{d=1}^{n}\dfrac{1}{d}\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}ijd^2[\gcd(i,j)=1]\\
&=\sum_{d=1}^{n}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}ij[\gcd(i,j)=1]\\
&=\sum_{d=1}^{n}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}ij\sum_{w|i,w|j}\mu(w)\\
&=\sum_{d=1}^{n}d\sum_{w=1}^{\lfloor\frac{n}{d}\rfloor}\mu(w)w^2\sum_{i=1}^{\lfloor\frac{n}{dw}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{dw}\rfloor}ij\\
&=\sum_{d=1}^{n}d\sum_{w=1}^{\lfloor\frac{n}{d}\rfloor}\mu(w)w^2(\dfrac{\lfloor\frac{n}{dw}\rfloor(\lfloor\frac{n}{dw}\rfloor+1)\lfloor\frac{m}{dw}\rfloor(\lfloor\frac{m}{dw}\rfloor+1)}{4})\end{aligned}
预处理 \mu(i)i^2 的前缀和,可以做到 O(n) 。可以通过。
注意这里需要处理一下 2 的逆元。
#include <bits/stdc++.h>
using namespace std;
#define int long long
const int N = 1e7 + 5, mod = 20101009, m2 = 10050505;
int mu[N], s[N]; bool isp[N];
vector<int> pr;
int g(int n, int m) { return 1ll * n * (n + 1) % mod * m2 % mod * 1ll * m % mod * (m + 1) % mod * 10050505 % mod; }
int sum(int n, int m) {
int ans = 0;
for(int l = 1, r = 1;l <= min(n, m);l = r + 1) {
r = min(n / (n / l), m / (m / l));
ans += (s[r] - s[l-1] + mod) * g(n / l, m / l) % mod;
} return ans;
} int solve(int n, int m) {
int ans = 0;
for(int l = 1, r = 1;l <= min(n, m);l = r + 1) {
r = min(n / (n / l), m / (m / l));
ans = (ans + 1ll * (l + r) * (r - l + 1) / 2 % mod * sum(n / l, m / l) % mod) % mod;
} return ans;
} signed main() {
mu[1] = 1;
for(int i = 2;i <= 10000000;i++) {
if(!isp[i]) pr.push_back(i), mu[i] = -1;
for(auto j : pr) {
if((i * j) > 10000000) break;
isp[i * j] = true;
if(i % j == 0) {
mu[i * j] = 0; break;
} else mu[i * j] = -mu[i];
}
} for(int i = 1;i <= 10000000;i++)
s[i] = (s[i-1] + 1ull * i * i % mod * (mu[i] + mod) % mod) % mod;
int n, m; cin >> n >> m;
cout << solve(n, m) << endl;
return 0;
}
P3327 [SDOI2015] 约数个数和
首先有一个神奇式子:
d(ij)=\displaystyle\sum_{x|i}\sum_{y|j}[\gcd(x,y)=1]
直接带入原题:
依然不妨 n<m 。
&=\sum_{x=1}^{n}\sum_{y=1}^{m}[\gcd(x,y)=1]\sum_{i=1}^{\lfloor\frac{n}{x}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{y}\rfloor}1\\
&=\sum_{x=1}^{n}\sum_{y=1}^{m}[\gcd(x,y)=1]\lfloor\frac{n}{x}\rfloor\lfloor\frac{m}{y}\rfloor\\
&=\sum_{x=1}^{n}\sum_{y=1}^{n}\sum_{w|x,w|y}\mu(w)\lfloor\frac{n}{x}\rfloor\lfloor\frac{m}{y}\rfloor\\
&=\sum_{w=1}^{n}\mu(w)\sum_{x=1}^{\lfloor\frac{n}{w}\rfloor}\sum_{y=1}^{\lfloor\frac{m}{w}\rfloor}\lfloor\frac{n}{x}\rfloor\lfloor\frac{m}{y}\rfloor\end{aligned}
套两层整除分块,即可做到 O(n) 。
如果预处理内层整除分块,可以做到 O(n\sqrt n+T\sqrt n) ,可以通过。
高亮部分是预处理。
#include <bits/stdc++.h>
using namespace std;
#define int long long
const int N = 6e4 + 5;
bool isp[N]; int mu[N], s[N]; vector<int> pr;
int h(int n) { int ans = 0;
for(int l = 1, r = 1;l <= n;l = r + 1)
r = n / (n / l), ans += (r-l+1) * (n/l);
return ans;
} int solve(int n, int m) { int ans = 0;
for(int l = 1, r = 1;l <= min(n, m);l = r + 1)
r = min(n / (n / l), m / (m / l)),
ans += (mu[r] - mu[l-1]) * s[n / l] * s[m / l];
return ans;
} signed main() { mu[1] = 1;
for(int i = 2;i <= 60000;i++) {
if(!isp[i]) pr.push_back(i), mu[i] = -1;
for(auto j : pr) {
if(i * j > 60000) break;
isp[i * j] = true;
if(i % j == 0) {
mu[i * j] = 0; break;
} mu[i * j] = -mu[i];
}
} for(int i = 2;i <= 60000;i++) mu[i] += mu[i-1];
for(int i = 1;i <= 60000;i++) s[i] = h(i);
int t; cin >> t; while(t--) {
int n, m; cin >> n >> m; cout << solve(n, m) << endl;
} return 0;
}
P6156 简单题
真 · 简单题
&=\sum_{d=1}^{n}d^{k+1}\mu^2(d)\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d}\rfloor}(i+j)^k[\gcd(i,j)=1]\\
&=\sum_{d=1}^{n}d^{k+1}\mu^2(d)\sum_{w=1}^{\lfloor\frac{n}{d}\rfloor}\mu(w)w^k\sum_{i=1}^{\lfloor\frac{n}{dw}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{dw}\rfloor}(i+j)^k\end{aligned}
记
$$f(x)=\displaystyle\sum_{i=1}^{x}\sum_{j=1}^{x}(i+j)^k$$
考虑如何从 $f(x-1)$ 递推到 $f(x)$。做一下差可以发现:
$$\begin{aligned}f(x)-f(x-1)&=\sum\limits_{i=1}^x\sum\limits_{j=1}^x(i+j)^{k}-\sum\limits_{i=1}^{x-1}\sum\limits_{j=1}^{x-1}(i+j)^{k}\\
&=\sum\limits_{i=1}^{x-1}(i+x)^{k}+\sum\limits_{j=1}^x(x+j)^{k}\end{aligned}$$
哦哈哈那不是又回到预处理 $k$ 次方了那线性筛就完了。
关于线性筛 $k$ 次幂:质数用快速幂,剩下的拿积性函数性质乘。复杂度不会算,但是接近于 $O(n)$。
时间复杂度 $O(n+\sqrt n)$,可以通过。
注意根据欧拉定理 $k$ 和 $k \bmod \varphi(p)$ 时的答案相同,所以我们先模一下 $\varphi(p)=p-1$。
注意精细取模。注意莫比乌斯函数会是负数。啊啊啊啊啊!!
```cpp
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int mod = 998244353;
const int N = 1e7 + 5, _N = 1e7; bool isp[N];
int mu[N], K[N], sk[N], muk[N], mu2k[N];
long long f[N];
vector<int> pr;
ll qpow(ll a, ll b) {
ll c = 1;
while(b) { if(b & 1) c = (c * a) % mod; a = (a * a) % mod; b >>= 1; }
return c;
} ll g(ll x) { ll _ans = 0;
for(int l = 1, r = 1;l <= x;l = r + 1)
r = x / (x / l),
_ans = (_ans + 1ll * (muk[r] - muk[l-1] + 3ll * mod) % mod * f[x / l] % mod) % mod;
return _ans;
} signed main() {
mu[1] = 1; K[1] = 1;
int n; ll k; cin >> n >> k; k %= (mod - 1);
for(int i = 2;i <= _N;i++) {
if(!isp[i]) {
pr.push_back(i); K[i] = qpow(i, k);
mu[i] = -1;
} for(auto j : pr) {
if(i * j > _N) break;
isp[i * j] = 1;
if(i % j == 0) {
mu[i * j] = 0; K[i * j] = (1ll * K[i] * K[j]) % mod; break;
} mu[i * j] = -mu[i]; K[i * j] = (1ll * K[i] * K[j]) % mod;
}
} for(int i = 1;i <= _N;i++) {
muk[i] = (muk[i-1] + 1ll * mu[i] * K[i] % mod) % mod;
mu2k[i] = (mu2k[i-1] + 1ll * mu[i] * mu[i] * K[i] % mod * i % mod) % mod;
sk[i] = (sk[i-1] + K[i]) % mod;
} for(int i = 1;i * 2 <= _N;i++)
f[i] = (f[i-1] + sk[2 * i - 1]) % mod,
f[i] = (f[i] - 2ll * sk[i] % mod + 2ll * mod) % mod,
f[i] = (f[i] + sk[2 * i]) % mod;
ll ans = 0;
for(int l = 1, r = 1;l <= n;l = r + 1)
r = n / (n / l),
ans = (ans + 1ll * (mu2k[r] - mu2k[l-1] + mod) % mod * g(n / l) % mod) % mod;
cout << ans << endl;
return 0;
}
```
### 后记
码字不易,能否给个赞 /wq ><
若对文章有任何问题和建议可以与作者私信交流。