莫比乌斯反演学习笔记
莫比乌斯反演学习笔记
不建议初学者看此垃圾博客，因为我是蒟蒻（
前置知识
数论函数
其为定义域为整数的函数。
常见的数论函数有：
\(\varphi(x)\) 欧拉函数，\([1,x]\) 中与 \(x\) 互质的数的个数。
\(\sigma_k(x)\) \(x\) 的因数的\(k\)次方和。
\(\omega(x)\) \(x\) 的本质不同质因子。
\(Id_k(x)\) 表示 \(x^k\)。
\(I_k(x)\) 在任意位置都是 \(1\)。
\(\epsilon(x)\) 仅在 \(x=1\) 时为 \(1\)，其余为 \(0\)，即为 \([x=1]\)。
\(\mu(x)\) 后面介绍。
基本运算：
\((f+g)(n)=f(n)+g(n)\)
\((f\cdot g)(n)=f(n)\times g(n)\)
积性函数
定义：
如果 \(gcd(x,y)=1\) 时 \(f(xy)=f(x)f(y)\)，则 \(f\) 为积性函数。
如果任意 \((x,y)\) 均满足，则 \(f\) 为完全积性函数。
常见积性函数：\(\varphi,\sigma,Id,I,\epsilon,\mu\)。
常见完全积性函数：\(Id,I,\epsilon\)。
Part.1 狄利克雷卷积
这是莫比乌斯反演中最重要的东西，就像多项式卷积在多项式中一般。
定义两个数论函数的卷积为：
计算可以枚举倍数做到 \(O(n\log n)\)
容易发现 \(\epsilon\) 在狄利克雷卷积中起到单位元的作用。
因为只有 \(d=1,n/d=n\) 时有值。
而定义 \((f*f^{-1})=\epsilon\)，为 \(f\) 的逆函数。
一些性质：
如果 \(f(n),g(n)\) 均为积性函数，则 \((f*g)(n)\) 也为积性函数。
如果 \(f(n)\) 为积性函数，且 \(f(1) \neq 0\)，则 \(f^{-1}(n)\) 也为积性函数。
证明：读者自证不难。
狄利克雷卷积满足结合律、交换律以及加法上的分配律。
一些常见的卷积结果：
证明：考虑展开得 \(\sum_{d|n}\varphi (n/d)\)，而与 \(n/d\) 互质，代表 \(gcd(x,n)=d\)，而把所有这样的数加起来，就得到了 \(n\)。
证明：展开一下使用 \((\varphi *I)(n)=n\) 即可。
Part.2 莫比乌斯函数/反演
定义莫比乌斯函数为：
\( mu(x)= \begin{cases} 1 & x=1\\ 0 & x 含有平方因子\\ (-1)^{w(x)} & others \end{cases}\)
这个定义一看就很奇怪，能有什么用呢？
一件十分有趣的事情，请看 vcr
是不是很神奇？\(\mu\) 竟然是常函数的逆函数。
咋证？
考虑 \(\sum_{d|n}\mu(d)\) 是啥。
由于含有平方因子就没有贡献了，所以其实就变成了 \(k\) 个元素的集合了，选出若干个数，假设选出 \(t\) 个，则贡献为 \((-1)^t\)。
所以上一发二项式反演或者随便推一下就发现是对的。
由此发现一个卷积就是 \(\mu*id_1=\varphi\)。
莫比乌斯反演又是啥呢？
如果有：
则有：
这其实是啥，就是卷积啊。
所以常说的莫反其实就是 \(\mu*I=\epsilon\)。
Part.3 例题
会了这些其实就可以做题了。
luogu P2257 YY的GCD
式子：
一个常见的思路（套路）就是枚举 \(\gcd\)。
考虑什么时候 \(\gcd(i,j)=k\)，就是在除去 \(k\) 后互质的数，所以可以只枚举 \(k\) 的倍数（不是倍数的一定不合法）。
发现了吗？其实后面的就是 \(\epsilon(\gcd(i,j))\)，可以直接莫反。
发现由于 \(d\) 是 \(\gcd(i,j)\) 的因数，所以一定是 \(i,j\) 的因数，就是上面的 \([d|i][d|j]\)。
所以可以只枚举 \(d\) 的倍数。
发现 \(kd\) 的乘积固定上界的，所以可以枚举 \(kd\)，再枚举 \(k\)。
更改一下得到：
发现最后的 \(\sum_{k|T,k\in prime}{\mu(T/k)}\)，本质就是卷积，可以枚举质数倍数做到 \(O(n\log\log n)\)。
现在钦定读者会数论分块（雾，其实挺简单的）。
设 \(f(n)=\sum_{d|n}\mu(d)\)。
直接处理前缀和后数论分块就做完了。
点击查看代码
const int M = 1e6 + 100, mod = 19930726, B = 300, MAXK = 22, N = 1e7 + 10; // !!!!!!!!
// #define int long long
int f[N], sum[N], mu[N], prime[N], cnt;
bool notPrime[N];
void Euler()
{
mu[1] = 1;
for (int i = 2; i < N; i++)
{
if (!notPrime[i]) prime[++cnt] = i, mu[i] = -1;
for (int j = 1; j <= cnt; j++)
{
if (i * prime[j] >= N) break;
notPrime[i * prime[j]] = 1;
if (i % prime[j]) mu[i * prime[j]] = -mu[i];
else break;
}
}
for (int i = 1; i <= cnt; i++)
{
for (int j = 1; prime[i] * j < N; j++)
{
f[prime[i] * j] += mu[j];
}
}
for (int i = 1; i < N; i++) sum[i] = sum[i - 1] + f[i];
}
ll calc(ll n, ll m)
{
ll ans = 0;
for (ll l = 1, r; l <= n; l = r + 1)
{
r = min(n / (n / l), m / (m / l));
// cerr << l << ' ' << r << '\n';
// cerr << sum[r] << '\n';
ans += (n / l) * (m / l) * (sum[r] - sum[l - 1]);
}
return ans;
}
signed main()
{
// freopen("data.in", "r", stdin); freopen("data.out", "w", stdout);
ios :: sync_with_stdio(false), cin.tie(0), cout.tie(0);
Euler();
int T; cin >> T;
while (T--)
{
int n, m; cin >> n >> m;
if (n > m) swap(n, m);
cout << calc(n, m) << '\n';
}
QED;
}
// 😡😡😡
luogu P3312 [SDOI2014] 数表
要求：
依旧枚举 \(\gcd(i,j)\)。
化简出来是：
这东西你发现没有上界十分简单，但是这个上界咋办？
十分简单，只需按 \(a\) 离线即可。
发现每个点会对 \(n/i\) 个点产生影响，这刚好对应着单点修改，而数论分块要求前缀和，这不就是树状数组板子吗？
具体就是将 \(\sigma_1(x)\) 升序排序，然后每次插入就可以了。
复杂度 \(O(n\log^2 n+T\sqrt n)\)。
点击查看代码
const int M = 1e6 + 100, B = 300, MAXK = 22, N = 1e5 + 10; // !!!!!!!!
#define int long long
const int mod = 2147483648ll;
int mu[N], sum[N];
struct Node
{
int n, m, a, id;
}q[N];
int ans[N], prime[N], cnt, d[N], n, m, c[N];
pii num[N];
bool notPrime[N];
void Update(int x, int k)
{
for (int i = x; i < N; i += lowbit(i)) (c[i] += k) %= mod;
}
int Query(int x)
{
int ans = 0;
for (int i = x; i; i -= lowbit(i)) (ans += c[i]) %= mod;
return ans;
}
int calc(int n, int m)
{
if (n > m) swap(n, m);
int ans = 0;
for (int l = 1, r; l <= n; l = r + 1)
{
r = min(n / (n / l), m / (m / l));
(ans += (n / l) * (m / l) % mod * (Query(r) - Query(l - 1)) % mod) %= mod;
}
return (ans + mod) % mod;
}
void Euler()
{
mu[1] = 1;
for (int i = 2; i < N; i++)
{
if (!notPrime[i]) prime[++cnt] = i, mu[i] = -1;
for (int j = 1; j <= cnt; j++)
{
if (i * prime[j] >= N) break;
notPrime[i * prime[j]] = 1;
if (i % prime[j]) mu[i * prime[j]] = -mu[i];
else break;
}
}
for (int i = 1; i < N; i++)
{
for (int j = 1; i * j < N; j++)
{
(d[i * j] += i) %= mod;
}
num[i] = {d[i], i};
}
sort(num + 1, num + N), sort(q + 1, q + 1 + m, [](Node x, Node y)
{
return x.a < y.a;
});
for (int i = 1, l = 1; i < N; i++)
{
int t = num[i].se;
if (num[i].fi != num[i - 1].fi)
{
while (l <= m && q[l].a < num[i].fi) ans[q[l].id] = calc(q[l].n, q[l].m), ++l;
}
for (int j = 1; t * j < N; j++)
{
Update(t * j, d[t] * mu[j]);
}
}
}
signed main()
{
// freopen("data.in", "r", stdin); freopen("data.out", "w", stdout);
ios :: sync_with_stdio(false), cin.tie(0), cout.tie(0);
cin >> m;
for (int i = 1; i <= m; i++) cin >> q[i].n >> q[i].m >> q[i].a, q[i].id = i;
Euler();
for (int i = 1; i <= m; i++) cout << ans[i] << '\n';
QED;
}
// 😡😡😡
luogu P3704 [SDOI2017] 数字表格
就是将乘法变换为加法，使用 exp 和 log 即可。
最终式子就是：
点击查看代码
const int M = 1e6 + 100, B = 300, MAXK = 22, mod = 1e9 + 7, N = 1e6 + 10; // !!!!!!!!
#define int long long
int f[N], h[N], mu[N], prime[N], cnt;
bool notPrime[N];
int qpow(int x, int b)
{
if (b == 1) return x;
if (b == -1) return qpow(x, mod - 2);
int res = 1;
while (b)
{
if (b & 1) res = res * x % mod;
x = x * x % mod;
b >>= 1;
}
return res;
}
void Euler()
{
f[1] = f[2] = 1;
for (int i = 3; i < N; i++) f[i] = (f[i - 1] + f[i - 2]) % mod;
mu[1] = 1;
for (int i = 2; i < N; i++)
{
if (!notPrime[i]) prime[++cnt] = i, mu[i] = -1;
for (int j = 1; j <= cnt; j++)
{
if (i * prime[j] >= N) break;
notPrime[i * prime[j]] = 1;
if (i % prime[j]) mu[i * prime[j]] = -mu[i];
else break;
}
}
for (int i = 0; i < N; i++) h[i] = 1;
for (int i = 1; i < N; i++)
{
for (int j = 1; i * j < N; j++)
{
h[i * j] = h[i * j] * qpow(f[i], mu[j]) % mod;
}
}
for (int i = 1; i < N; i++) h[i] = h[i] * h[i - 1] % mod;
}
int calc(int n, int m)
{
if (n > m) swap(n, m);
int ans = 1;
for (int l = 1, r; l <= n; l = r + 1)
{
r = min(n / (n / l), m / (m / l));
ans = ans * qpow(h[r] * qpow(h[l - 1], mod - 2) % mod, n / l * (m / l)) % mod;
}
return ans;
}
signed main()
{
// freopen("data.in", "r", stdin); freopen("data.out", "w", stdout);
ios :: sync_with_stdio(false), cin.tie(0), cout.tie(0);
Euler();
int T; cin >> T;
for (int i = 1; i <= T; i++)
{
int n, m; cin >> n >> m;
cout << calc(n, m) << '\n';
}
QED;
}
// 😡😡😡
luogu P4240 毒瘤之神的考验
十分好的题，先咕了。
可以看 bjt 写的 link
\( \begin{aligned} awa \end{aligned} \)
浙公网安备 33010602011771号