莫比乌斯反演

一句话总结：题面极其简洁并且大多能一眼出算法，代码短但是式子又臭又长。

前置知识：数论分块：数论全家桶 - lfyszy（别问为啥没有狄利克雷卷积)

莫比乌斯函数

μ\mu 为莫比乌斯函数，定义为：

μ(n)={1n=10n 含有平方因子(−1)kk 为 n 的本质不同质因子个数\mu(n)= \begin{cases} 1&n=1\\ 0&n\text{ 含有平方因子}\\ (-1)^k&k\text{ 为 }n\text{ 的本质不同质因子个数}\\ \end{cases}

同时莫比乌斯函数有以下性质：

∑d∣nμ(d)={1n=10n≠1\sum_{d\mid n}\mu(d)= \begin{cases} 1&n=1\\ 0&n\neq 1\\ \end{cases}

证明：

设 n=∏i=1kpici,n′=∏i=1kpi\displaystyle n=\prod_{i=1}^k{p_i}^{c_i},n'=\prod_{i=1}^k p_i

那么 ∑d∣nμ(d)=∑d∣n′μ(d)=∑i=0k(ki)⋅(−1)i=(1+(−1))k\displaystyle\sum_{d\mid n}\mu(d)=\sum_{d\mid n'}\mu(d)=\sum_{i=0}^k \dbinom{k}{i}\cdot(-1)^i=(1+(-1))^k

这样由二项式定理推导出，所以易得只有在 n=1n=1 时此式等于 1{1}。

观察到，莫比乌斯函数是一个积性函数，所以我们可以在欧拉筛的过程中预处理：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
bool st[N];
int mu[N];
vector<int> p;
void getprime(int n)
{
    mu[1] = 1;
    for(int i = 2; i <= n; i ++)
    {
        if(!st[i]) p.push_back(i), mu[i] = -1;
        for(auto j : p)
        {
            if(i * j > n) break;
            st[i * j] = 1;
            if(i % j == 0) {mu[i * j] = 0; break;}
            mu[i * j] = -mu[i];
        }
        mu[i] += mu[i - 1];
    }
}

莫比乌斯变换

公式

形式一：如果有 f(n)=∑d∣ng(d)f(n)=\sum_{d\mid n}g(d)，那么有 g(n)=∑d∣nμ(d)f(nd)g(n)=\sum_{d\mid n}\mu(d)f(\frac{n}{d})。

形式二：如果有 f(n)=∑n∣dg(d)f(n)=\sum_{n|d}g(d)，那么有 g(n)=∑n∣dμ(dn)f(d)g(n)=\sum_{n|d}\mu(\frac{d}{n})f(d)。

证明

形式一：

∑d∣nμ(d)f(nd)=∑d∣nμ(d)∑k∣ndg(k)=∑k∣ng(k)∑d∣nkμ(d)=g(n)\sum_{d\mid n}\mu(d)f(\frac{n}{d})=\sum_{d\mid n}\mu(d)\sum_{k\mid \frac{n}{d}}g(k)=\sum_{k\mid n}g(k)\sum_{d\mid \frac{n}{k}}\mu(d)=g(n)

考虑倒推，首先将 f(nd)f(\frac{n}{d}) 替换为 ∑k∣ndg(k)\sum_{k\mid \frac{n}{d}}g(k)，然后改变枚举顺序，而 ∑d∣nkμ(d)\sum_{d\mid \frac{n}{k}}\mu(d) 根据上文中提到的性质显然在 n=kn=k 时才取 1{1}，所以化简为 g(n)g(n)。

形式二

∑n∣dμ(dn)f(d)=∑k=1+∞μ(k)f(kn)=∑k=1+∞μ(k)∑kn∣dg(d)=∑n∣dg(d)∑k∣dnμ(k)=g(n)\sum_{n|d}{\mu(\frac{d}{n})f(d)}=\sum_{k=1}^{+\infty}{\mu(k)f(kn)}=\sum_{k=1}^{+\infty}{\mu(k)\sum_{kn|d}{g(d)}}=\sum_{n|d}{g(d)\sum_{k|\frac{d}{n}}{\mu(k)}}=g(n)

其实本质是一样的，考虑枚举 k=ndk=\frac{n}{d}，接着替换掉 ff，便可以改变枚举顺序，化简为 g(n)g(n)。

应用

注：本文中默认 n<mn < m

[HAOI2011] Problem b

可以转化成二维平面上思考，则可以利用容斥分成四个部分解决。

对于每个部分：

∑i=1n∑j=1m[gcd⁡(i,j)=k]=∑i=1⌊nk⌋∑j=1⌊mk⌋∑d∣gcd⁡(i,j)μ(d)=∑d=1μ(d)∑i=1⌊nk⌋[d∣i]∑j=1⌊mk⌋[d∣j]=∑d=1⌊nk⌋μ(d)⌊nkd⌋⌊mkd⌋\begin{aligned} &\sum_{i=1}^{n}\sum_{j=1}^{m}[\gcd(i,j)=k]\\ &=\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{k}\rfloor}\sum_{d\mid \gcd(i,j)}\mu(d)\\ &=\sum_{d=1}\mu(d)\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}[d\mid i]\sum_{j=1}^{\lfloor\frac{m}{k}\rfloor}[d\mid j]\\ &=\sum_{d=1}^{\lfloor \frac{n}{k}\rfloor}\mu(d)\lfloor\frac{n}{kd}\rfloor\lfloor\frac{m}{kd}\rfloor\\ \end{aligned}

直接用数论分块解决就可以了。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
// problem: P3455 [POI2007] ZAP-Queries
// date: 2024.2.19-16:35
#include <bits/stdc++.h>
#define int long long
#define FH signed
using namespace std;
const int N = 5e4 + 10;
bool st[N];
int mu[N];
vector<int> p;
void getprime(int n)
{
    mu[1] = 1;
    for(int i = 2; i <= n; i ++)
    {
        if(!st[i]) p.push_back(i), mu[i] = -1;
        for(auto j : p)
        {
            if(i * j > n) break;
            st[i * j] = 1;
            if(i % j == 0) {mu[i * j] = 0; break;}
            mu[i * j] = -mu[i];
        }
        mu[i] += mu[i - 1];
    }
}
int solve(int n, int m, int k)
{
    int res = 0; n /= k, m /= k;
    if(n > m) swap(n, m);
    for(int i = 1; i <= n; i ++)
    {
        int up = min(n / (n / i), m / (m / i));
        res += (mu[up] - mu[i - 1]) * (n / i) * (m / i);
        i = up;
    }
    return res;
}
FH main()
{
    // freopen(".in", "r", stdin);
    // freopen(".out", "w", stdout);
    ios::sync_with_stdio(0);
    cin.tie(0); cout.tie(0);
    getprime(N - 1);
    int t; cin >> t;
    while(t --)
    {
        int a1, b1, a2, b2, k;
        cin >> a1 >> b1 >> a2 >> b2 >> k;
        cout << solve(b1, b2, k) - solve(b1, a2 - 1, k) - solve(b2, a1 - 1, k) + solve(a1 - 1, a2 - 1, k) << "\n";
    }
    return 0;
} //FH

[国家集训队] Crash的数字表格 / JZPTAB

∑i=1n∑j=1mlcm(i,j)=∑i=1n∑j=1mijgcd(i,j)=∑i=1n∑j=1m∑d∣i,d∣j,gcd(id,jd)=1ijd=∑d=1n∑d∣in∑d∣jm[gcd(i,j)=d]ijd=∑d=1nd∑i=1⌊nd⌋∑j=1⌊md⌋[gcd(i,j)=1]ij=∑d=1nd∑i=1⌊nd⌋∑j=1⌊md⌋∑k∣gcd(i,j)μ(k)ij=∑d=1nd∑k=1⌊nd⌋μ(k)∑k∣i⌊nd⌋∑k∣j⌊md⌋ij=∑d=1nd∑k=1⌊nd⌋μ(k)k2∑i=1⌊ndk⌋∑j=1⌊mdk⌋ij=∑d=1nd∑k=1⌊nd⌋μ(k)k2ndk(ndk−1)2ndk(ndk−1)2\begin{aligned} &\sum_{i=1}^n\sum_{j=1}^mlcm(i,j)\\ &=\sum_{i=1}^n\sum_{j=1}^m\dfrac{ij}{gcd(i,j)}\\ &=\sum_{i=1}^n\sum_{j=1}^m\sum_{d|i,d|j,gcd(\frac{i}{d}, \frac{j}{d})=1}\dfrac{ij}{d}\\ &=\sum_{d=1}^{n}\sum_{d|i}^n\sum_{d|j}^m[gcd(i,j)=d]\dfrac{ij}{d}\\ &=\sum_{d=1}^{n}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}[gcd(i,j)=1]ij\\ &=\sum_{d=1}^{n}d\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}\sum_{k|gcd(i,j)}\mu(k)ij\\ &=\sum_{d=1}^{n}d\sum_{k=1}^{\lfloor\frac{n}{d}\rfloor}\mu(k)\sum_{k|i}^{\lfloor\frac{n}{d}\rfloor}\sum_{k|j}^{\lfloor\frac{m}{d}\rfloor}ij\\ &=\sum_{d=1}^{n}d\sum_{k=1}^{\lfloor\frac{n}{d}\rfloor}\mu(k)k^2\sum_{i=1}^{\lfloor\frac{n}{dk}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{dk}\rfloor}ij\\ &=\sum_{d=1}^{n}d\sum_{k=1}^{\lfloor\frac{n}{d}\rfloor}\mu(k)k^2\frac{\frac{n}{dk}(\frac{n}{dk}-1)}{2}\frac{\frac{n}{dk}(\frac{n}{dk}-1)}{2}\\ \end{aligned}

令 mul(n,m)=n(n−1)2m(m−1)2mul(n,m)=\frac{n(n-1)}{2}\frac{m(m-1)}{2}，sum(n,m)=∑k=1nμ(k)k2mul(nk,mk)sum(n,m)=\sum_{k=1}^{n}\mu(k)k^2mul(\frac{n}{k},\frac{m}{k})

原式=∑d=1nd×sum(nd,md)\text{原式} =\sum_{d=1}^{n}d\times sum(\frac{n}{d},\frac{m}{d})

这就比较好说了，用数论分块计算 sum()sum() 与答案，解决~~~

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
// problem: P1829 [国家集训队] Crash的数字表格 / JZPTAB
// date: 2024.2.20-16:18
#include <bits/stdc++.h>
#define int long long
#define FH signed
using namespace std;
const int N = 1e7 + 10, mod = 20101009;
bool st[N];
int mu[N];
vector<int> p;
void getprime(int n)
{
    mu[1] = 1;
    for(int i = 2; i <= n; i ++)
    {
        if(!st[i]) p.push_back(i), mu[i] = -1;
        for(auto j : p)
        {
            if(i * j > n) break;
            st[i * j] = 1;
            if(i % j == 0) {mu[i * j] = 0; break;}
            mu[i * j] = -mu[i];
        }
        mu[i] = (mu[i - 1] + (mu[i] + mod) % mod * i % mod * i % mod) % mod;
    }
}
int mul(int n, int m) {return (n * (n + 1) / 2 % mod) * (m * (m + 1) / 2 % mod) % mod;}
int sum(int n , int m)
{
    int res = 0;
    for(int k = 1; k <= n; k ++)
    {
        int up = min(n / (n / k), m / (m / k));
        res = (res + (mu[up] - mu[k - 1] + mod) % mod * mul(n / k, m / k) % mod) % mod;
        k = up;
    }
    return res;
}
FH main()
{
    // freopen(".in", "r", stdin);
    // freopen(".out", "w", stdout);
    // ios::sync_with_stdio(0);
    // cin.tie(0); cout.tie(0);
    int n, m, res = 0; cin >> n >> m;
    if(n > m) swap(n, m);
    getprime(n);
    for(int d = 1; d <= n; d ++)
    {
        int up = min(n / (n / d), m / (m / d));
        res = (res + (d + up) * (up - d + 1) / 2 % mod * sum(n / d, m / d) % mod) % mod;
        d = up;
    }
    cout << res << "\n";
    return 0;
} //FH

[SDOI2017] 数字表格

这道题稍有难度（

∏i=1n∏j=1mf(gcd(i,j))=∏i=1n∏j=1m∏gcd(i,j)=df(d)=∏d=1n∏i=1⌊nd⌋∏j=1⌊md⌋[gcd(i,j)=1]f(d)=∏d=1nf(d)∑i=1⌊nd⌋∑j=1⌊md⌋[gcd(i,j)=1]∑i=1⌊nd⌋∑j=1⌊md⌋[gcd(i,j)=1]=∑k=1⌊nd⌋∑i=1⌊nkd⌋∑j=1⌊mkd⌋μ(k)=∑k=1⌊nd⌋μ(k)⌊nkd⌋⌊mkd⌋原式=∏d=1nf(d)∑k=1⌊nd⌋μ(k)⌊nkd⌋⌊mkd⌋=∏d=1n(∏k=1⌊nd⌋f(d)μ(k))⌊nkd⌋⌊mkd⌋=∏T=1n(∏d∣Tf(d)μ(Td))⌊nT⌋⌊mT⌋\begin{aligned} &\prod_{i=1}^{n}\prod_{j=1}^{m}f(gcd(i,j))\\ &=\prod_{i=1}^{n}\prod_{j=1}^{m}\prod_{gcd(i,j)=d}f(d)\\ &=\prod_{d=1}^{n}\prod_{i=1}^{\lfloor\frac{n}{d}\rfloor}\prod_{j=1}^{\lfloor\frac{m}{d}\rfloor}[gcd(i,j)=1]f(d)\\ &=\prod_{d=1}^{n}f(d)^{\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}[gcd(i,j)=1]}\\ &\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}[gcd(i,j)=1]\\ &=\sum_{k=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{i=1}^{\lfloor\frac{n}{kd}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{kd}\rfloor}\mu(k)\\ &=\sum_{k=1}^{\lfloor\frac{n}{d}\rfloor}\mu(k)\lfloor\frac{n}{kd}\rfloor\lfloor\frac{m}{kd}\rfloor\\ &\text{原式}=\prod_{d=1}^{n}f(d)^{\sum_{k=1}^{\lfloor\frac{n}{d}\rfloor}\mu(k)\lfloor\frac{n}{kd}\rfloor\lfloor\frac{m}{kd}\rfloor}\\ &=\prod_{d=1}^{n}(\prod_{k=1}^{\lfloor\frac{n}{d}\rfloor}f(d)^{\mu(k)} )^{\lfloor\frac{n}{kd}\rfloor\lfloor\frac{m}{kd}\rfloor}\\ &=\prod_{T=1}^{n}(\prod_{d\mid T}f(d)^{\mu(\frac{T}{d})} )^{\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor}\\ \end{aligned}

这样就可以直接枚举 dd 预处理括号中部分与其逆元的前缀积，使用数论分块解决剩余部分。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
// problem: P3704 [SDOI2017] 数字表格
// date: 2024.2.21
#include<bits/stdc++.h>
#define int long long
#define FH signed
using namespace std;
const int N = 1e6 + 10, mod = 1e9 + 7;
int qpow(int a, int b, int p) {int res = 1; while(b){if(b & 1) res = res * a % p; a = a * a % p; b >>= 1;} return res;}
bool st[N];
int mu[N], f[N], g[N];
vector<int> p;
void getprime(int n)
{
    mu[1] = 1;
    for(int i = 2; i <= n; i ++)
    {
        if(!st[i]) p.push_back(i), mu[i] = -1;
        for(auto j : p)
        {
            if(i * j > n) break;
            st[i * j] = 1;
            if(i % j == 0) {mu[i * j] = 0; break;}
            mu[i * j] = -mu[i];
        }
    }
}
void init()
{
	getprime(N - 1);
	for(int i = 0; i < N; i ++) f[i] = 1, g[i] = 1;
	for(int i = 1, l1 = 0, l2 = 1; i < N; i ++)
	{
		l1 = (l1 + l2) % mod, l2 = (l1 - l2 + mod) % mod;
		int tp[] = {qpow(l1, mod - 2, mod), 1, l1};
		for(int j = 1; j * i < N; j ++)
			f[i * j] = f[i * j] * tp[mu[j] + 1] % mod,
			g[i * j] = g[i * j] * tp[1 - mu[j]] % mod;
	}
	for(int i = 2; i < N; i ++)
		f[i] = f[i] * f[i - 1] % mod,
		g[i] = g[i] * g[i - 1] % mod;
}

void solve()
{
	int n, m, res = 1;
	cin >> n >> m;
	if(n > m) swap(n, m);
	for(int i = 1; i <= n; i ++)
	{
		int up = min(n / (n / i), m / (m / i));
		res = res * qpow(f[up] * g[i - 1] % mod, (n / i) * (m / i) % (mod - 1), mod) % mod;
		i = up;
	}
	cout << res << "\n";
}

FH main()
{
	ios::sync_with_stdio(0);
	cin.tie(0); cout.tie(0);
	init();
	int t; cin >> t;
	while(t --) solve();
    return 0;
}//FH

讲个笑话：

image.png

YY的GCD

接着 problem B 的证明，因为推到到这里的部分相同。

∑k=1n∑d=1⌊nk⌋μ(d)⌊nkd⌋⌊mkd⌋(k∈prime)令T=kd原式=∑T=1n∑k∣Tμ(k)⌊nT⌋⌊mT⌋(k∈prime)=∑T=1n⌊nT⌋⌊mT⌋∑k∣Tμ(Tk)(k∈prime)\begin{aligned} &\sum_{k=1}^{n}\sum_{d=1}^{\lfloor \frac{n}{k}\rfloor}\mu(d)\lfloor\frac{n}{kd}\rfloor\lfloor\frac{m}{kd}\rfloor(k\in prime)\\ &\text{令}T=kd\\ &\text{原式}=\sum_{T=1}^{n}\sum_{k|T}\mu(k)\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor(k\in prime)\\ &=\sum_{T=1}^{n}\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor\sum_{k|T}\mu(\frac{T}{k})(k\in prime)\\ \end{aligned}

这竟然是可以预处理后一个 ∑\sum 的内容，就可以直接数论分块。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
// problem: P2257 YY的GCD
// date: 2024.2.23-14:41
#include <bits/stdc++.h>
#define int long long
#define FH signed
using namespace std;
const int N = 1e7 + 10;
bool st[N];
int mu[N], f[N];
vector<int> p;
void getprime(int n)
{
    mu[1] = 1;
    for(int i = 2; i <= n; i ++)
    {
        if(!st[i]) p.push_back(i), mu[i] = -1;
        for(auto j : p)
        {
            if(i * j > n) break;
            st[i * j] = 1;
            if(i % j == 0) {mu[i * j] = 0; break;}
            mu[i * j] = -mu[i];
        }
    }
    for(auto i : p)
        for(int j = 1; j * i <= n; j ++) f[i * j] += mu[j];
    for(int i = 1; i <= n; i ++) f[i] += f[i - 1];
}
void solve()
{
    int n, m, res = 0; cin >> n >> m; if(n > m) swap(n, m);
    for(int i = 1; i <= n; i ++)
    {
        int up = min(n / (n / i), m / (m / i));
        res = res + (n / i) * (m / i) * (f[up] - f[i - 1]);
        i = up;
    }
    cout << res << "\n";
}
FH main()
{
    // freopen(".in", "r", stdin);
    // freopen(".out", "w", stdout);
    // ios::sync_with_stdio(0);
    // cin.tie(0); cout.tie(0);
    getprime(N - 1);
    int t; cin >> t;
    while(t --) solve();
    return 0;
} //FH