莫比乌斯反演题目结（上）
首先是喜闻乐见的几个入门基础题，连题面基本都是一样的，流程是预处理mu函数，得到输入数据后整除分块（如果时间复杂度需要的话），思路主要是套用下图第二个公式：
（截图来自电科bilibili上的算法讲堂大家可以去看哦）
洛谷2257
题意：全部gcd（1～n，1～m）为质数的个数。
刚开始不用管他要求的数有什么特性，就列一般式子：
（截图来自pengym的博客题解）
但这个题用第二条推出的式子还得变一变，就直观感受去变！d的范围？n是质数？把F(d)提出去？
下图中T即为上图中的d，t即为n。
式子写到这里就可以写代码了，先把mu全加上，这样写：
1 rep(j, 1, cnt) { 2 for (int i = 1; i * primes[j] <= n; i++) { 3 g[primes[j] * i] += mu[i]; 4 } 5 }
然后整除分块是某一块乘上的（n/T）*（m/T）是不变的，所以这一块一起算，代码是这样的：
1 rep(i, 1, n) sum[i] = sum[i-1] + g[I];//这个是预处理
1 ll ans = 0ll; 2 for (int l = 1, r; l <= n; l = r + 1) { 3 r = min(n / (n / l), m / (m / l)); 4 ans += (sum[r] - sum[l - 1]) * (n / l) * (m / l); 5 }
总代码（部分）：
1 void pre(int n) { 2 mu[1] = 1; 3 rep(i, 2, n) { 4 if (!vis[i]) { 5 mu[i] = -1; 6 primes[++cnt] = i; 7 } 8 for (int j = 1; j <= cnt && primes[j] * i <= n; j++) { 9 vis[primes[j] * i] = true; 10 if (i % primes[j] == 0) break; 11 mu[primes[j] * i] = -mu[i]; 12 } 13 } 14 15 rep(j, 1, cnt) { 16 for (int i = 1; i * primes[j] <= n; i++) { 17 g[primes[j] * i] += mu[i]; 18 } 19 } 20 rep(i, 1, n) sum[i] = sum[i-1] + g[i]; 21 } 22 23 int main() { 24 pre(maxn); 25 for (T = ri; T; T--) { 26 n = ri, m = ri; 27 if (n > m) swap(n, m); 28 29 ll ans = 0ll; 30 for (int l = 1, r; l <= n; l = r + 1) { 31 r = min(n / (n / l), m / (m / l)); 32 ans += (sum[r] - sum[l - 1]) * (n / l) * (m / l); 33 } 34 printf("%lld\n", ans); 35 } 36 return 0; 37 }
HDU1695
（复习一下式子）
题意：gcd（1～b，1～d）== k的个数。
跟上一题很像，甚至这道题根本不用整除分块，也不用进一步推式子……就是求f(k)，直接拿第二个式子扫一遍加和就行了。另外去个重，本题（x，y）和（y，x）算一种。
1 void pre(int n) { 2 mu[1] = 1; 3 rep(i, 2, n) { 4 if (!vis[i]) { 5 primes[++cnt] = i; 6 mu[i] = -1; 7 } 8 for (int j = 1; j <= cnt && primes[j] * i <= n; j++) { 9 vis[primes[j] * i] = true; 10 if (i % primes[j] == 0) break; 11 mu[primes[j] * i] = -mu[i]; 12 } 13 } 14 } 15 16 ll cal(int b, int d) { 17 ll ret = 0ll; 18 for (int i = 1; i * k <= b; i++) { 19 int t = i * k; 20 ret += (ll)mu[i] * (b / t) * (d / t); 21 } 22 return ret; 23 } 24 25 int main() { 26 pre(maxn - 5); 27 for (int T = ri, i = 1; i <= T; i++) { 28 a = ri, b = ri, c = ri, d = ri, k = ri; 29 if (b > d) swap(b, d); 30 printf("Case %d: %lld\n", i, k ? cal(b, d) - cal(b, b) / 2 : 0); 31 } 32 return 0; 33 }
BZOJ1101
n次询问gcd(1~a, 1~b) == d。同一道题。
1 #pragma comment(linker, "/STACK:1024000000,1024000000") 2 #include <cstdio> 3 #include <cstring> 4 #include <cstdlib> 5 #include <cmath> 6 #include <ctime> 7 #include <cctype> 8 #include <climits> 9 #include <iostream> 10 #include <iomanip> 11 #include <algorithm> 12 #include <string> 13 #include <sstream> 14 #include <stack> 15 #include <queue> 16 #include <set> 17 #include <map> 18 #include <vector> 19 #include <list> 20 #include <fstream> 21 #include <bitset> 22 #define init(a, b) memset(a, b, sizeof(a)) 23 #define rep(i, a, b) for (int i = a; i <= b; i++) 24 #define irep(i, a, b) for (int i = a; i >= b; i--) 25 using namespace std; 26 27 typedef double db; 28 typedef long long ll; 29 typedef unsigned long long ull; 30 typedef pair<int, int> P; 31 const int inf = 0x3f3f3f3f; 32 const ll INF = 1e18; 33 34 template <typename T> void read(T &x) { 35 x = 0; 36 int s = 1, c = getchar(); 37 for (; !isdigit(c); c = getchar()) 38 if (c == '-') s = -1; 39 for (; isdigit(c); c = getchar()) 40 x = x * 10 + c - 48; 41 x *= s; 42 } 43 44 template <typename T> void write(T x) { 45 if (x < 0) x = -x, putchar('-'); 46 if (x > 9) write(x / 10); 47 putchar(x % 10 + '0'); 48 } 49 50 template <typename T> void writeln(T x) { 51 write(x); 52 puts(""); 53 } 54 55 const int maxn = 5e4 + 5; 56 int n, a, b, d; 57 int mu[maxn], primes[maxn], tot, sum[maxn]; 58 bool vis[maxn]; 59 60 void pre(int n) { 61 mu[1] = 1; 62 rep(i, 2, n) { 63 if (!vis[i]) { 64 primes[++tot] = i; 65 mu[i] = -1; 66 } 67 for (int j = 1; j <= tot && primes[j] * i <= n; j++) { 68 vis[primes[j] * i] = true; 69 if (i % primes[j] == 0) break; 70 mu[primes[j] * i] = - mu[i]; 71 } 72 } 73 rep(i, 1, n) sum[i] = sum[i - 1] + mu[i]; 74 } 75 76 ll solve(int n, int m) { 77 if (n > m) swap(n, m); 78 ll ret = 0ll; 79 for (int l = 1, r; l <= n; l = r + 1) { 80 r = min(n / (n / l), m / (m / l)); 81 ret += (ll)(sum[r] - sum[l - 1]) * (n / l) * (m / l); 82 } 83 return ret; 84 } 85 86 int main() { 87 pre(maxn - 5); 88 for (read(n); n; n--) { 89 read(a), read(b), read(d); 90 writeln(solve(a / d, b / d)); 91 } 92 return 0; 93 }
BZOJ2301
跟HDU1695的区别是变成了gcd（a～b，c～d）== k的个数。
还是那么做，但是ans容斥一下：
1 printf("%lld\n", cal(b, d) - cal(b, c) - cal(a, d) + cal(a, c));
唔，实际上应该是cal(b,c-1)、cal(a-1,d)……主要是上一行代码直接操作掉了：
1 int main() { 2 pre(maxn - 5); 3 for (T = ri; T; T--) { 4 a = ri, b = ri, c = ri, d = ri, k = ri; 5 a = (a - 1) / k, b = b / k, c = (c - 1) / k, d = d / k; 6 printf("%lld\n", cal(b, d) - cal(b, c) - cal(a, d) + cal(a, c)); 7 } 8 return 0; 9 }
除以k有降低复杂度的功效，gcd和范围同时除以k使得问题变成求f(1)！也就是所有的mu[d]*F(d)加和，不用预处理k的倍数了。上一题也可以这么写的。
先对照一下公式，旧书不厌百回读：
1 void pre(int n) { 2 mu[1] = 1; 3 rep(i, 2, n) { 4 if (!vis[i]) { 5 primes[++cnt] = i; 6 mu[i] = -1; 7 } 8 for (int j = 1; j <= cnt && primes[j] * i <= n; j++) { 9 vis[primes[j] * i] = true; 10 if (i % primes[j] == 0) break; 11 mu[primes[j] * i] = -mu[i]; 12 } 13 } 14 rep(i, 1, n) sum[i] = sum[i - 1] + mu[i]; 15 } 16 17 ll cal(int n, int m) { 18 if (n > m) swap(n, m); 19 ll ans = 0; 20 21 for (int l = 1, r; l <= n; l = r + 1) { 22 r = min(n / (n / l), m / (m / l)); 23 ans += (sum[r] - sum[l - 1]) * (n / l) * (m / l); 24 } 25 return ans; 26 }
POJ3904
题意：首先，这次范围为给定数列，不是从1～n、从1～m所有的数这种大范围了；然后抽取四个数使他们gcd为1。也就是共有多少四元组的gcd是1，还是计数～
通过这个题可以得知套路公式中的F(n)一般是根据莫比乌斯而跟f(d)有关的一个式子，另外又会等于一个与题目相关的、直观感受的、可以直接求出的式子。
比如本题还是使用套路公式，但是F(d)不再是(n/d)*(m/d)，而是数列中d的倍数的这些数中，随便抽4个，有多少种取法：
1 ll ans = 0ll; 2 rep(i, 1, maxn - 5) { 3 ans += mu[i] * C(c[i], 4);//gcd是1的莫比乌斯公式 4 }
总代码：
1 void mobi(int n) { 2 mu[1] = 1; 3 rep(i, 2, n) { 4 if (!vis[i]) { 5 primes[++cnt] = i; 6 mu[i] = -1; 7 } 8 for (int j = 1; j <= cnt && primes[j] * i <= n; j++) { 9 vis[primes[j] * i] = true; 10 if (i % primes[j] == 0) break; 11 mu[primes[j] * i] = -mu[i]; 12 } 13 } 14 } 15 16 void get_cnt() { 17 init(c, 0); 18 rep(i, 1, n) { 19 int x = a[i]; 20 for (int j = 1; j * j <= x; j++) { 21 if (x % j == 0) { 22 c[j]++; 23 if (j * j != x) c[x / j]++; 24 } 25 } 26 } 27 } 28 29 ll C(int n, int k) { 30 if (n < k) return 0; 31 32 ll ret = 1ll; 33 rep(i, 1, k) { 34 ret = ret * (n - i + 1) / i; 35 } 36 return ret; 37 } 38 39 int main() { 40 ios_base::sync_with_stdio(false); 41 cin.tie(0); 42 43 mobi(maxn - 5); 44 while (cin >> n) { 45 rep(i, 1, n) cin >> a[i]; 46 get_cnt(); 47 48 ll ans = 0ll; 49 rep(i, 1, maxn - 5) { 50 ans += mu[i] * C(c[i], 4); 51 } 52 cout << ans << endl; 53 } 54 return 0; 55 }
BZOJ2005
做完以上几题以后此题可独立做出。还是个gcd（1～n，1～m），但这次不是个数加和而是值加和。参考代码：
1 #pragma comment(linker, "/STACK:1024000000,1024000000") 2 #include <cstdio> 3 #include <cstring> 4 #include <cstdlib> 5 #include <cmath> 6 #include <ctime> 7 #include <cctype> 8 #include <climits> 9 #include <iostream> 10 #include <iomanip> 11 #include <algorithm> 12 #include <string> 13 #include <sstream> 14 #include <stack> 15 #include <queue> 16 #include <set> 17 #include <map> 18 #include <vector> 19 #include <list> 20 #include <fstream> 21 #include <bitset> 22 #define init(a, b) memset(a, b, sizeof(a)) 23 #define rep(i, a, b) for (int i = a; i <= b; i++) 24 #define irep(i, a, b) for (int i = a; i >= b; i--) 25 using namespace std; 26 27 typedef double db; 28 typedef long long ll; 29 typedef unsigned long long ull; 30 typedef pair<int, int> P; 31 const int inf = 0x3f3f3f3f; 32 const ll INF = 1e18; 33 34 template <typename T> void read(T &x) { 35 x = 0; 36 int s = 1, c = getchar(); 37 for (; !isdigit(c); c = getchar()) 38 if (c == '-') s = -1; 39 for (; isdigit(c); c = getchar()) 40 x = x * 10 + c - 48; 41 x *= s; 42 } 43 44 template <typename T> void write(T x) { 45 if (x < 0) x = -x, putchar('-'); 46 if (x > 9) write(x / 10); 47 putchar(x % 10 + '0'); 48 } 49 50 template <typename T> void writeln(T x) { 51 write(x); 52 puts(""); 53 } 54 55 const int maxn = 1e5 + 5; 56 int n, m; 57 int mu[maxn], primes[maxn], tot; 58 bool vis[maxn]; 59 ll ans; 60 61 void pre(int n) { 62 mu[1] = 1; 63 rep(i, 2, n) { 64 if (!vis[i]) { 65 primes[++tot] = i; 66 mu[i] = -1; 67 } 68 for (int j = 1; j <= tot && i * primes[j] <= n; j++) { 69 vis[i * primes[j]] = true; 70 if (i % primes[j] == 0) break; 71 mu[i * primes[j]] = -mu[i]; 72 } 73 } 74 } 75 76 ll cal(int n, int m) { 77 ll ret = 0ll; 78 rep(i, 1, n) { 79 ret += (ll)mu[i] * (n / i) * (m / i); 80 } 81 return ret; 82 } 83 84 int main() { 85 read(n), read(m); 86 if (n > m) swap(n, m); 87 pre(maxn - 5); 88 for (int l = 1, r; l <= n; l = r + 1) { 89 r = min(n / (n / l), m / (m / l)); 90 ans += (ll)(l + r) * (r - l + 1) / 2 * cal(n / l, m / l); 91 } 92 writeln(2 * ans - (ll)n * m); 93 return 0; 94 }
BZOJ2154
这个我就推不出来了Orz。这式子还蛮复杂的，推荐这个博客，证明得蛮清楚的。
核心：
#pragma comment(linker, "/STACK:1024000000,1024000000") #include <cstdio> #include <cstring> #include <cstdlib> #include <cmath> #include <ctime> #include <cctype> #include <climits> #include <iostream> #include <iomanip> #include <algorithm> #include <string> #include <sstream> #include <stack> #include <queue> #include <set> #include <map> #include <vector> #include <list> #include <fstream> #include <bitset> #define init(a, b) memset(a, b, sizeof(a)) #define rep(i, a, b) for (int i = a; i <= b; i++) #define irep(i, a, b) for (int i = a; i >= b; i--) using namespace std; typedef double db; typedef long long ll; typedef unsigned long long ull; typedef pair<int, int> P; const int inf = 0x3f3f3f3f; const ll INF = 1e18; template <typename T> void read(T &x) { x = 0; int s = 1, c = getchar(); for (; !isdigit(c); c = getchar()) if (c == '-') s = -1; for (; isdigit(c); c = getchar()) x = x * 10 + c - 48; x *= s; } template <typename T> void write(T x) { if (x < 0) x = -x, putchar('-'); if (x > 9) write(x / 10); putchar(x % 10 + '0'); } template <typename T> void writeln(T x) { write(x); puts(""); } const int maxn = 1e7 + 5; const int mod = 20101009; int n, m; ll ans; int mu[maxn], primes[maxn], tot; ll sum[maxn]; bool vis[maxn]; void pre(int n) { mu[1] = 1; rep(i, 2, n) { if (!vis[i]) { primes[++tot] = i; mu[i] = -1; } for (int j = 1; j <= tot && primes[j] <= n / i; j++) { vis[primes[j] * i] = true; if (i % primes[j] == 0) break; mu[primes[j] * i] = -mu[i]; } } rep(i, 1, n) { if (mu[i] < 0) mu[i] = mod - 1; sum[i] = sum[i - 1] + (ll)i * i % mod * mu[i] % mod; sum[i] %= mod; } } ll S(ll x, ll y) { return (x * (x + 1) / 2 % mod) * (y * (y + 1) / 2 % mod) % mod; } ll F(int x, int y) { ll ret = 0ll; for (int l = 1, r; l <= x; l = r + 1) { r = min(x / (x / l), y / (y / l)); ret += (sum[r] - sum[l - 1] + mod) % mod * S(x / l, y / l); ret %= mod; } return ret; } int main() { read(n), read(m); if (n > m) swap(n, m); pre(m); for (int l = 1, r; l <= n; l = r + 1) { r = min(n / (n / l), m / (m / l)); ans += (ll)(r - l + 1) * (l + r) / 2 % mod * F(n / l, m / l); ans %= mod; } writeln(ans); return 0; }
