莫比乌斯反演：
1. $F(n)=sum_{d|n}f(d) Rightarrow f(n)=sum_{d|n}mu(d)F(frac{n}{d})$
2. $F(n)=sum_{n|d}f(d) Rightarrow f(n)=sum_{n|d}mu(frac{d}{n})F(d)$
1. [POI2007]Zap
题目链接：bzoj - 1101
题意：对于给定的整数a，b和d，有多少正整数对x，y，满足x$leq$a，y$leq$b，并且gcd(x,y)=d
思路：设$F(n)$表示gcd(x,y)为n的倍数(包括n本身)的数量，$f(n)$表示gcd(x,y)为n的数量
显然有$F(n)=sum_{n|d}f(d)$
所以由莫比乌斯反演得$f(n)=sum_{n|d}mu(frac{d}{n})F(d)$，而$F(d)=left lfloor frac{a}{d} ight floor left lfloor frac{b}{d} ight floor$(根据$F(n)$的定义理解)
代入$f(n)=sum_{n|d}mu(frac{d}{n}) left lfloor frac{a}{d} ight floor left lfloor frac{b}{d} ight floor$
令$k=frac{d}{n}$可得$f(n)=sum_{k=1}^{min(left lfloor frac{a}{n} ight floor , left lfloor frac{b}{n} ight floor)}mu(k) left lfloor frac{a}{kn} ight floor left lfloor frac{b}{kn} ight floor$
预处理出$mu(k)$的前缀和，然后整除分块计算即可
#include <iostream> #include <algorithm> #include <cstring> #include <cstdio> using namespace std; const int N = 500010; typedef long long ll; int isp[N], p[N], mu[N], sum[N]; int tot, T, a, b, d; void init(int n) { mu[1] = 1; for (int i = 2; i <= n; i++) { if (!isp[i]) { p[++tot] = i; mu[i] = -1; } for (int j = 1; j <= tot && p[j] <= n / i; j++) { isp[i * p[j]] = 1; if (0 == i % p[j]) { mu[i * p[j]] = 0; break; } else mu[i * p[j]] = -mu[i]; } } for (int i = 1; i <= n; i++) sum[i] = sum[i - 1] + mu[i]; } int main() { init(N - 5); scanf("%d", &T); while (T--) { scanf("%d%d%d", &a, &b, &d); ll res = 0; for (int l = 1, r = 1; l <= min(a / d, b / d); l = r + 1) { r = min(a / (a / l), b / (b / l)); res += 1ll * (sum[r] - sum[l - 1]) * (a / (l * d)) * (b / (l * d)); } printf("%lld ", res); } return 0; }
2. YY的GCD
题目链接：bzoj - 2820
题意：对于给定的整数n，m，有多少正整数对x，y，满足x$leq$n，y$leq$m，并且gcd(x,y)为质数
思路：设$F(n)$表示gcd(x,y)为n的倍数(包括n本身)的数量，$f(n)$表示gcd(x,y)为n的数量
显然有$F(n)=sum_{n|d}f(d)$
所以由莫比乌斯反演得$f(n)=sum_{n|d}mu(frac{d}{n})F(d)$，而$F(d)=left lfloor frac{a}{d} ight floor left lfloor frac{b}{d} ight floor$(根据$F(n)$的定义理解)
代入$f(n)=sum_{n|d}mu(frac{d}{n}) left lfloor frac{a}{d} ight floor left lfloor frac{b}{d} ight floor$
所以$res=sum_{pin prime}f(p)=sum_{pin prime}sum_{p|d}mu(frac{d}{p}) left lfloor frac{a}{d} ight floor left lfloor frac{b}{d} ight floor$
因为对于质数p的每个倍数d求$mu(frac{d}{p}) left lfloor frac{a}{d} ight floor left lfloor frac{b}{d} ight floor$所得到的值和对于d的每个质因子p求$mu(frac{d}{p}) left lfloor frac{a}{d} ight floor left lfloor frac{b}{d} ight floor$所得到的值是一样的
所以$res=sum_{d=1}^{min(n,m)} sum_{pin{prime} p|d} mu(frac{d}{p}) left lfloor frac{a}{d} ight floor left lfloor frac{b}{d} ight floor$
用整数分块化简后得$sum_{d=1}^{min(n,m)}(left lfloor frac{a}{d} ight floor left lfloor frac{b}{d} ight floor sum_{pin{prime} p|d} mu(frac{d}{p}))$
令$c(d)=sum_{pin{prime} p|d} mu(frac{d}{p})$
预处理出$c(d)$的前缀和$sum(d)$，利用整数分块计算即可
#include <iostream> #include <algorithm> #include <cstring> #include <cstdio> using namespace std; typedef long long ll; const int N = 10000010; int isp[N], p[N], mu[N]; int T, tot, n, m; ll sum[N]; void init(int n) { mu[1] = 1; for (int i = 2; i <= n; i++) { if (!isp[i]) { p[++tot] = i; mu[i] = -1; } for (int j = 1; j <= tot && p[j] <= n / i; j++) { isp[i * p[j]] = 1; if (0 == i % p[j]) { mu[i * p[j]] = 0; break; } else mu[i * p[j]] = -mu[i]; } } for (int i = 1; i <= tot; i++) for (int j = 1; p[i] * j <= n; j++) sum[p[i] * j] += mu[j]; for (int i = 1; i <= n; i++) sum[i] += sum[i - 1]; } int main() { init(N - 5); scanf("%d", &T); while (T--) { scanf("%d%d", &n, &m); ll res = 0; for (int l = 1, r = 0; l <= min(n, m); l = r + 1) { r = min(n / (n / l), m / (m / l)); res = res + 1ll * (sum[r] - sum[l - 1]) * (n / l) * (m / l); } printf("%lld ", res); } return 0; }
3. [HAOI2011]Problem b
题目链接：bzoj - 2301
题意：对于给定的整数a，b，c，d和k，有多少正整数对x，y，满足a$leq$x$leq$b，c$leq$y$leq$d，并且gcd(x,y)=k
思路：利用第1题的结果和容斥原理即可得到答案，$res=solve(b,d)-solve(a-1,d)-solve(b,c-1)+solve(a-1,c-1)$
#include <iostream> #include <algorithm> #include <cstring> #include <cstdio> using namespace std; const int N = 50010; int isp[N], p[N], mu[N], sum[N]; int tot, T, a, b, c, d, k; inline int read() { int s = 0, f = 1; char ch = getchar(); while (ch < '0' || ch > '9') { if ('-' == ch) f = -1; ch = getchar(); } while (ch >= '0' && ch <= '9') { s = s * 10 + ch - '0'; ch = getchar(); } return s * f; } inline void write(int x) { if (x < 0) { putchar('-'); x = -x; } if (x > 9) write(x / 10); putchar(x % 10 + '0'); } void init(int n) { mu[1] = 1; for (int i = 2; i <= n; i++) { if (!isp[i]) { p[++tot] = i; mu[i] = -1; } for (int j = 1; j <= tot && p[j] <= n / i; j++) { isp[i * p[j]] = 1; if (0 == i % p[j]) { mu[i * p[j]] = 0; break; } else mu[i * p[j]] = -mu[i]; } } for (int i = 1; i <= n; i++) sum[i] = sum[i - 1] + mu[i]; } int solve(int a, int b, int n) { int res = 0, ta = a / n, tb = b / n; for (int l = 1, r = 1; l <= min(ta, tb); l = r + 1) { r = min(a / (a / l), b / (b / l)); res += (sum[r] - sum[l - 1]) * (a / (n * l)) * (b / (n * l)); } return res; } int main() { init(N - 5); T = read(); while (T--) { a = read(), b = read(), c = read(), d = read(), k = read(); int res = solve(b, d, k) - solve(a - 1, d, k) - solve(b, c - 1, k) + solve(a - 1, c - 1, k); write(res); puts(""); } return 0; }
4. [Sdoi2014]数表
题目链接：bzoj - 3529
题意：有一张n*m的数表，其第i行第j列的数值为能同时整除i和j的所有自然数之和，给定a，计算表中不大于a的数之和
思路：假设n<m，先不考虑不大于a的这个约束条件，那么$res=sum_{i=1}^{n}sum_{j=1}^{m}sigma (gcd(i,j))$，其中$sigma $表示约数之和
令$d=gcd(i,j)$，枚举d可得$res=sum_{d=1}^{n}sum_{i=1}^{n}sum_{j=1}^{m}sigma(d)[gcd(i,j)=d]$
将$sigma(d)$提前得$res=sum_{d=1}^{n}[sigma(d)sum_{i=1}^{n}sum_{j=1}^{m}[gcd(i,j)=d]]$
后半部分$sum_{i=1}^{n}sum_{j=1}^{m}[gcd(i,j)=d]$相当于求有多少对i，j满足i$leq n$，$jleq m$且$gcd(i,j)=d$
后半部分用莫比乌斯化简一下得$res=sum_{d=1}^{n}[sigma(d)sum_{k=1}^{left lfloor frac{n}{d} ight floor} mu(k) left lfloor frac{n}{kd} ight floor left lfloor frac{m}{kd} ight floor]$
令$T=kd$得$res=sum_{d=1}^{n}[sigma(d)sum_{T=1 d|T}^{n} mu(frac{T}{d}) left lfloor frac{n}{T} ight floor left lfloor frac{m}{T} ight floor]$
交换一下顺序得$res=sum_{T=1}^{n}[left lfloor frac{n}{T} ight floor left lfloor frac{m}{T} ight floor sum_{d=1 d|T}^{n}sigma(d) mu(frac{T}{d})]$
现在考虑限制条件，只有当$sigma(d)leq a$时，$sigma(d) mu(frac{T}{d})$才对答案有贡献，所以我们先预处理出$sigma(d)$，对询问按a排序，对$sigma(d)$排序，当a变大时，某些$sigma(d) mu(frac{T}{d})$开始对答案产生贡献，此时枚举d得倍数T，然后用树状数组维护$sigma(d) mu(frac{T}{d})$的前缀和就可以计算出$sum_{d=1 d|T}^{n}sigma(d) mu(frac{T}{d})$，这个d对别的询问也会产生贡献(询问是按a排序的)，对于每个询问用整数分块计算
#include <iostream> #include <algorithm> #include <cstring> #include <cstdio> using namespace std; const int N = 100010; struct node { int n, m, a, id; }; int T, c[N], mu[N], res[N]; int isp[N], p[N], tot; node q[N], d[N]; bool cmp(node a, node b) { return a.a < b.a; } void init(int n) { mu[1] = 1; for (int i = 2; i <= n; i++) { if (!isp[i]) { p[++tot] = i; mu[i] = -1; } for (int j = 1; j <= tot && p[j] <= n / i; j++) { isp[i * p[j]] = 1; if (0 == i % p[j]) { mu[i * p[j]] = 0; break; } mu[i * p[j]] = -mu[i]; } } for (int i = 1; i <= n; i++) { d[i].id = i; for (int j = i; j <= n; j += i) { d[j].a += i; } } sort(d + 1, d + n + 1, cmp); } inline int lowbit(int x) { return x & (-x); } void add(int x, int v) { while (x <= N - 5) { c[x] += v; x += lowbit(x); } } int ask(int x) { int res = 0; while (x > 0) { res += c[x]; x -= lowbit(x); } return res; } int solve(int n, int m) { int res = 0; for (int l = 1, r = 1; l <= n; l = r + 1) { r = min(n / (n / l), m / (m / l)); res += (n / l) * (m / l) * (ask(r) - ask(l - 1)); } return res; } int main() { init(N - 5); scanf("%d", &T); for (int i = 1; i <= T; i++) { scanf("%d%d%d", &q[i].n, &q[i].m, &q[i].a); if (q[i].n > q[i].m) swap(q[i].n, q[i].m); q[i].id = i; } sort(q + 1, q + T + 1, cmp); for (int i = 1, j = 1; i <= T; i++) { while (j <= N - 5 && d[j].a <= q[i].a) { for (int k = d[j].id; k <= N - 5; k += d[j].id) { add(k, d[j].a * mu[k / d[j].id]); } j++; } res[q[i].id] = solve(q[i].n, q[i].m); } for (int i = 1; i <= T; i++) printf("%d ", res[i] & (~(1 << 31))); return 0; }
5. Crash的数字表格
题目链接：bzoj - 2154
题意：有一张n*m的数表，其第i行第j列的数值为lcm(i,j)，求数表内所有数的和
思路：假设n<m，$res=sum_{i=1}^{n} sum_{j=1}^{m} lcm(i,j)=sum_{i=1}^{n}sum_{j=1}^{m}frac{ij}{gcd(i,j)}$
令$d=gcd(i,j)$，枚举d得$res=sum_{d=1}^{n}sum_{i=1}^{n} sum_{j=1}^{m} frac{ij}{d} [gcd(i,j)=d]=sum_{d=1}^{n}sum_{i=1}^{left lfloor frac{n}{d} ight floor} sum_{j=1}^{left lfloor frac{m}{d} ight floor} ijd [gcd(i,j)=1]$
把d移到前面有$res=sum_{d=1}^{n}[dsum_{i=1}^{left lfloor frac{n}{d} ight floor} sum_{j=1}^{left lfloor frac{m}{d} ight floor} ij [gcd(i,j)=1]]$
根据莫比乌斯得性质得$res=sum_{d=1}^{n}[dsum_{i=1}^{left lfloor frac{n}{d} ight floor} sum_{j=1}^{left lfloor frac{m}{d} ight floor} ij sum_{k|gcd(i,j)}mu(k)]$
枚举k得$res=sum_{d=1}^{n}[dsum_{k=1}^{left lfloor frac{n}{d} ight floor} sum_{i=1}^{left lfloor frac{n}{dk} ight floor} sum_{j=1}^{left lfloor frac{m}{dk} ight floor} ijk^{2}mu(k)]$
移一下项得$res=sum_{d=1}^{n}[dsum_{k=1}^{left lfloor frac{n}{d} ight floor} [ k^{2}mu(k) sum_{i=1}^{left lfloor frac{n}{dk} ight floor} i sum_{j=1}^{left lfloor frac{m}{dk} ight floor} j]]$
可以发现$sum_{i=1}^{left lfloor frac{n}{dk} ight floor} i$和$sum_{j=1}^{left lfloor frac{m}{dk} ight floor} j$都是等差数列求和，$sum_{k=1}^{left lfloor frac{n}{d} ight floor} k^{2}mu(k)$可以通过预处理前缀和来求，两次整数分块即可
#include <iostream> #include <algorithm> #include <cstring> #include <cstdio> using namespace std; const int N = 10000010; const int mod = 20101009; int isp[N], p[N], sum[N], mu[N]; int tot, n, m; void init(int n) { mu[1] = 1; for (int i = 2; i <= n; i++) { if (!isp[i]) { p[++tot] = i; mu[i] = -1; } for (int j = 1; j <= tot && p[j] <= n / i; j++) { isp[i * p[j]] = 1; if (0 == i % p[j]) { mu[i * p[j]] = 0; break; } mu[i * p[j]] = -mu[i]; } } for (int i = 1; i <= n; i++) sum[i] = (sum[i - 1] + 1ll * i * i % mod * (mu[i] + mod)) % mod; } int qs(int x, int y) { return (1ll * x * (x + 1) / 2 % mod) * (1ll * y * (y + 1) / 2 % mod) % mod; } int calc(int n, int m) { int res = 0; for (int l = 1, r = 1; l <= n; l = r + 1) { r = min(n / (n / l), m / (m / l)); res = (res + 1ll * (sum[r] - sum[l - 1] + mod) * qs(n / l, m / l) % mod) % mod; } return res; } int solve(int n, int m) { int res = 0; for (int l = 1, r = 1; l <= n; l = r + 1) { r = min(n / (n / l), m / (m / l)); res = (res + 1ll * (r - l + 1) * (l + r) / 2 % mod * calc(n / l, m / l) % mod) % mod; } return res; } int main() { init(N - 5); scanf("%d%d", &n, &m); if (n > m) swap(n, m); printf("%d ", solve(n, m)); return 0; }
6. AT5200 [AGC038C] LCMs
题目链接：AtCoder - agc038_c
题意：给定一个长度为n得序列$A_{1}，A_{2}，...，A_{n}$，求$sum_{i=1}^{n} sum_{j=i+1}^{n} lcm(A_{i},A_{j})$
思路：令M=max($A_{i}$)，$sum_{i=1}^{n} sum_{j=i+1}^{n} lcm(A_{i},A_{j})=frac{sum_{i=1}^{n} sum_{j=1}^{n} lcm(A_{i},A_{j})-sum_{i=1}^{n}A_{i}}{2}$
考虑如何求$res=sum_{i=1}^{n} sum_{j=1}^{n} lcm(A_{i},A_{j})$
$sum_{i=1}^{n} sum_{j=1}^{n} lcm(A_{i},A_{j})=sum_{i=1}^{n} sum_{j=1}^{n} frac{A_{i}A_{j}}{gcd(A_{i},A_{j})}$
令$d=gcd(A_{i},A_{j})$得$sum_{i=1}^{n} sum_{j=1}^{n} lcm(A_{i},A_{j})=sum_{d=1}^{M} frac{1}{d} sum_{i=1}^{n} sum_{j=1}^{n} A_{i}A_{j}[gcd(A_{i},A_{j})=d]=sum_{d=1}^{M} frac{1}{d} sum_{i=1}^{n} sum_{j=1}^{n} A_{i}A_{j}[gcd(frac{A_{i}}{d},frac{A_{j}}{d})=1][d|A_{i}][d|A_{j}]$
根据莫比乌斯反演的性质得$res=sum_{d=1}^{M} frac{1}{d} sum_{i=1}^{n} sum_{j=1}^{n} A_{i}A_{j} sum_{k|gcd(frac{A_{i}}{d},frac{A_{j}}{d})} mu(k) [d|A_{i}][d|A_{j}]$
将k提到前面得$res=sum_{d=1}^{M} frac{1}{d} sum_{k=1}^{left lfloor frac{M}{d} ight floor} mu(k) sum_{i=1}^{n} sum_{j=1}^{n} A_{i}A_{j} [kd|A_{i}][kd|A_{j}]$
所以最后$res=sum_{d=1}^{M} frac{1}{d} sum_{k=1}^{left lfloor frac{M}{d} ight floor} mu(k) (sum_{i=1}^{n}[kd|A_{i}]A_{i})^{2}$
设$f(t)=(sum_{i=1}^{n}[t|A_{i}]A_{i})^{2}$，显然$f(t)$可以在$O(nlnn)$的时间内预处理出来，再通过枚举d和k暴力求解即可
#include <iostream> #include <algorithm> #include <cstring> #include <cstdio> using namespace std; typedef long long ll; const int N = 1000010; const ll mod = 998244353; int isp[N], p[N], mu[N], n, tot, a[N], imax; ll inv[N], f[N], c[N], sum; void init(int n) { mu[1] = 1; for (int i = 2; i <= n; i++) { if (!isp[i]) { p[++tot] = i; mu[i] = -1; } for (int j = 1; j <= tot && p[j] <= n / i; j++) { isp[i * p[j]] = 1; if (0 == i % p[j]) { mu[i * p[j]] = 0; break; } else mu[i * p[j]] = -mu[i]; } } inv[0] = inv[1] = 1; for (int i = 2; i <= n; i++) inv[i] = (mod - mod / i) * inv[mod % i] % mod; } int main() { scanf("%d", &n); for (int i = 1; i <= n; i++) { scanf("%lld", &a[i]); imax = max(imax, a[i]); sum = (sum + a[i]) % mod; c[a[i]] += 1; } init(imax); for (int i = 1; i <= imax; i++) { for (int j = i; j <= imax; j += i) f[i] = (f[i] + c[j] * j) % mod; f[i] = (f[i] * f[i]) % mod; } ll res = 0; for (int d = 1; d <= imax; d++) { for (int k = 1, t = d; t <= imax; t += d, k++) { res = (res + inv[d] * mu[k] % mod * f[t] % mod) % mod; } } res = ((res - sum) % mod + mod) % mod * inv[2] % mod; printf("%lld ", res); return 0; }