【笔记】莫比乌斯反演
0 约定
\([n]=[1,n]\cap\mathtt Z\)
1 数论分块
形如 $S(n)=\sum\limits_{i=1}^nf(i)g(\left \lfloor \dfrac{n}{i} \right \rfloor) $。
可以在 \(O(\sqrt n)\) 的时间复杂度内求解。
1.1 解法
-
对于 \(1\le i\le \sqrt n\)，
显然，\(i\) 最多 \(\sqrt n\) 种取值，故而 \(\left \lfloor \dfrac{n}{i} \right \rfloor\) 最多有 \(\sqrt n\) 种取值。
-
对于 \(\sqrt n<i\le n\)，
有 \(\left \lfloor \dfrac{n}{i} \right \rfloor\le\sqrt n\)，最多亦有 \(\sqrt n\) 种取值。
综上，如果可以 \(O(1)\) 得到 \(f\) 的区间和，那么可以通过枚举 \(\left \lfloor \dfrac{n}{i} \right \rfloor\) 的取值，\(O(2\sqrt n)\approx O(\sqrt n)\) 求解。
具体代码如下：
int res = 0;
for (int i = 1, j; i <= n; i = j+1) {
j = n/(n/i);
res += (f(j) - f(i-1)) * g(n/i);
}
其中 \(j\) 为对于当前 \(i\) 满足 \(\left \lfloor \dfrac{n}{i} \right \rfloor=\left \lfloor \dfrac{n}{j} \right \rfloor\) 的最大的 \(j\)。
2 莫比乌斯函数
对于一些函数 \(f(n)\)，如果很难直接求出它的值，而容易求出其倍数和或约数和 \(g(n)\)，那么可以通过莫比乌斯反演，简化运算得到 \(f(n)\) 的值。
2.1 定义
定义莫比乌斯函数：
注：
- 「含有平方因子」是指存在 \(p\) 为质数，\(\alpha\) 为满足 \(p^{\alpha}=n\) 的最大整数，\(\alpha\ge2\)。eg. \(\mu(24)=0\)。
2.2 性质
-
\(\mu(n)\) 是 不完全 积性函数，即 \(\mu(x)\mu(y)=\mu(xy)\) 当且仅当 \((x,y)=1\)。
-
狄利克雷卷积 相关：
-
\((\mu*\operatorname I)=\begin{cases} 1&n=1\\ 0&n\neq 1\\ \end{cases}=\varepsilon\)
证明与下类似。
-
\((\mu*\operatorname{Id})=\varphi\)
证明：
设：\(n=p_1^{\alpha_1}p_2^{\alpha_2}\cdots p_k^{\alpha_k}\)，\(p\) 为质数，\(\alpha\) 为正整数。
\[\begin{aligned} (\mu*\operatorname{Id})(n) &=\sum_{d|n}\mu(d)\dfrac n d\\ &\begin{array}{c} =\left(\mu(1)\frac{p_1^{\alpha_1}}{1}+\mu(p_1)\frac{p_1^{\alpha_1}}{p_1}+\cdots+\mu(p_1^{\alpha_1})\frac{p_1^{\alpha_1}}{p_1^{\alpha_1}}\right) \\ \hspace{1.5em}\vdots \\ \hspace{1.5em}\left(\mu(1)\frac{p_k^{\alpha_k}}{1}+\mu(p_k)\frac{p_k^{\alpha_k}}{p_k}+\cdots+\mu(p_k^{\alpha_k})\frac{p_k^{\alpha_k}}{p_k^{\alpha_k}}\right) \\ \end{array}&\text{【直接把括号展开，就可以发现她是对的】}\\ &=\left(p_1^{\alpha_1}-\frac{p_1^{\alpha_1}}{p_1}\right)\left(p_2^{\alpha_2}-\frac{p_2^{\alpha_2}}{p_2}\right)\cdots\left(p_k^{\alpha_k}-\frac{p_k^{\alpha_k}}{p_k}\right)&\text{【把 } \mu \text{ 的值带进去】}\\ &=\frac{(p_1-1)p_1^{\alpha_1}}{p_1}\times\frac{(p_2-1)p_2^{\alpha_2}}{p_2}\times\cdots\times\frac{(p_k-1)p_k^{\alpha_k}}{p_k}\\ &=(p_1-1)p_1^{\alpha_1-1}(p_2-1)p_2^{\alpha_2-1}\cdots(p_k-1)p_k^{\alpha_k-1}\\ &=\varphi(n)&\text{【 }\varphi\text{ 的一种公式计算法】} \end{aligned}\]
-
3 线性筛求莫比乌斯函数
线性质数筛(欧拉筛)的一大优势在于我们可以知道每个数的 最小质因子，从而通过递推线性求解 积性函数。
void init (int _n) {
for (int i = 2; i <= _n; i++) is_prime[i] = 1;
mu[1] = 1;
for (int i = 2; i <= _n; i++) {
if (is_prime[i])
prime.push_back(i),
mu[i] = -1;
for (int j : prime) {
if (1ll*i*j > _n) break;
is_prime[i*j] = 0;
if (i%j == 0) {
mu[i*j] = 0;
break;
}
mu[i*j] = -mu[i];
}
}
}
4 莫比乌斯反演
设两个个数论函数 \(f(x),F(x)\)，若有 \(F(x)=\sum_{d|x}f(d)\)。
则有：
5 例题
5.1 简单的数学题
5.1.1 Description
已知：\(n\)。
求：
\[\left(\sum_{i=1}^n\sum\limits_{j=1}^n ij \gcd(i,j)\right) \bmod p \]
感性感知一下，\(\gcd(i,j)\) 的取值实际上很少，所以突破口在于枚举她的取值。
设：
\(f(d)=\sum\limits_{i,j\in[n]\wedge\gcd(i,j)=d}ij\)
\(F(x)=\sum_{d|x}f(d)\)
\(G(x)=\sum\limits_{i=1}^ni=\dfrac{n(n+1)}2\)
答案即为 \(Ans=\sum\limits_{d=1}^nf(d)d\)。
下考虑求解 \(f(d)\)：
从而：
推柿子：
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int lim = 1e7;
int mu[lim+5], phi[lim+5];
bool is_prime[lim+5];
vector<int> prime;
ll _S[lim+5], inv2, inv6, ans, p, n;
unordered_map<ll, ll> S;
ll qpow (ll a, ll b) {
ll res = 1;
for (; b; b >>= 1, (a *= a) %= p)
if (b & 1) (res *= a) %= p;
return res;
}
void init (int _n) {
for (int i = 2; i <= _n; i++) is_prime[i] = 1;
phi[1] = 1;
mu[1] = 1;
for (int i = 2; i <= _n; i++) {
if (is_prime[i])
prime.push_back(i),
phi[i] = i-1,
mu[i] = -1;
for (int j : prime) {
if (1ll*i*j > _n) break;
is_prime[i*j] = 0;
if (i%j == 0) {
phi[i*j] = j * phi[i];
mu[i*j] = 0;
break;
}
phi[i*j] = phi[j] * phi[i];
mu[i*j] = -mu[i];
}
}
for (int i = 1; i <= _n; i++)
_S[i] = (_S[i-1] + 1ll*i*i%p * phi[i]%p) % p;
}
ll F1 (ll x) { return x%=p, x*(x+1)%p*inv2%p; }
ll F2 (ll x) { return x%=p, x*(x+1)%p*(x*2+1)%p*inv6%p; }
ll F3 (ll x) { return x%=p, F1(x) * F1(x) % p; }
ll solve (ll x) {
if (x <= lim) return _S[x];
if (S[x]) return S[x];
ll res = F3(x);
for (ll i = 2, j; i <= x; i = j+1) {
j = x/(x/i);
res = (res - solve(x/i) * (F2(j) - F2(i-1) + p) % p + p) % p;
}
return S[x] = res;
}
int main () {
ios::sync_with_stdio(0);
cin.tie(0), cout.tie(0);
cin >> p >> n;
inv2 = qpow(2, p-2), inv6 = qpow(6, p-2);
init(lim);
for (ll i = 1, j; i <= n; i = j+1) {
j = n/(n/i);
(ans += F1(n/i) * F1(n/i)%p * ((solve(j) - solve(i-1)+p)%p) % p) %= p;
}
cout << ans;
return 0;
}
