前段时间教练讲了莫比乌斯反演，决定发一篇博客来加深自己的理解。（众所周知我是个很懒的人，不愿意去打非常复杂的(latex)，所以我的题解里很少有需要很复杂的推式子的题。。。）
【hdu1695】GCD
题面
题解
题意：求
[sum_{i=1}^a sum_{j=1}^b [gcd(i,j)==d]
]
不难发现和这个式子是等价的
[sum_{i=1}^{a/d} sum_{j=1}^{b/d} [gcd(i,j)==1]
]
为了方便之后的表达，后面的(a)和(b)分别为(a/d)和(b/d)
我们设
[f(d)=sum_{i=1}^{a} sum_{j=1}^{b} [gcd(i,j) == d]
]
[F(n) = sum_{n|d}f(d)
]
此时
[F(x)=lfloor frac{a}{x}
floor lfloor frac{b}{x}
floor
]
根据莫比乌斯反演公式
[F(n)=sum_{n|d}f(d)
]
[f(n)=sum_{n|d} mu(frac{d}{n})F(d)
]
那么我们的答案其实就是(f(1))
这里还有(1)个问题。
那就是(gcd(x,y))和(gcd(y,x))是一样的。
那我们要将总答案减去这种情况(/2)。
（当时写这份代码的时候我还不会整除分块，所以代码里并没有...）
代码
#include <bits/stdc++.h>
const int maxn = 1e5 + 10;
typedef long long ll;
bool vis[maxn];
int pri[maxn], mu[maxn];
int tot, n, m, i, j, k, T, a, b, ctot;
inline void prep(int n) {
mu[1] = 1;
for(int i = 2;i <= n;i++) {
if(!vis[i]) { pri[ ++pri[0] ] = i; mu[i] = -1; }
for(int j = 1;j <= pri[0] && i * pri[j] <= n;j++) {
vis[ i * pri[j] ] = 1;
if(i % pri[j])
mu[ i * pri[j] ] = -mu[i];
else {
mu[ i * pri[j] ] = 0;
break;
}
}
}
}
int main() {
prep(maxn - 10);
scanf("%d",&T);
while(T--) {
scanf("%d %d %d",&a,&b,&k);
printf("Case %d: ",++ctot);
if(!k) { puts("0"); continue; }
a /= k; b /= k;
if(a > b)
a = b;
ll r1 = 0, r2 = 0;
for(int i = 1;i <= std::min(a,b);i++)
r1 += 1ll * mu[i] * (a / i) * (b / i);
for(int i = 1;i <= std::min(a,b);i++)
r2 += 1ll * mu[i] * (std::min(a,b) / i) * (std::min(a,b) / i);
printf("%lld
",r1 - r2 / 2);
}
return 0;
}
[HAOI2011] Problem b
题面
题解
好像就是上面那道题加个差分...
具体看代码吧。
代码
// luogu-judger-enable-o2
#include <bits/stdc++.h>
const int maxn = 1e5 + 10;
typedef long long ll;
bool vis[maxn];
int pri[maxn], mu[maxn], s[maxn];
int tot, n, m, i, j, k, T, a, b, ctot, c, d;
inline void prep(int n) {
mu[1] = 1;
for(int i = 2;i <= n;i++) {
if(!vis[i]) { pri[ ++pri[0] ] = i; mu[i] = -1; }
for(int j = 1;j <= pri[0] && i * pri[j] <= n;j++) {
vis[ i * pri[j] ] = 1;
if(i % pri[j])
mu[ i * pri[j] ] = -mu[i];
else {
mu[ i * pri[j] ] = 0;
break;
}
}
}
for(int i = 1;i <= n;i++)
s[i] = s[i - 1] + mu[i];
}
inline ll get(int a,int b) {
int c = std::min(a,b); int res = 0;
for(int l = 1, r;l <= c;l = r + 1) {
r = std::min(a / (a / l),b / (b / l));
res += (1ll * a / (1ll * l * k)) * (1ll * b / (1ll * r * k)) * (s[r] - s[l - 1]);
}
return res;
}
int main() {
prep(maxn - 10);
scanf("%d",&T);
while(T--) {
scanf("%d %d %d %d %d",&a,&b,&c,&d,&k);
printf("%lld
",get(b,d) - get(b,c - 1) - get(a - 1,d) + get(a - 1,c - 1));
}
return 0;
}
YY的GCD
题面
题解
题意：
[sum_{i=1}^n sum_{j=1}^m sum_{k} [gcd(i,j)==k &&k in prime ]
]
第一道题设同样的东西：
[f(d)=sum_{i=1}^n sum_{j=1}^m [gcd(i,j)==d]
]
[F(n) = sum_{n|d}f(d)
]
[F(d) = lfloor frac{n}{d}
floor lfloor frac{m}{d}
floor
]
此时
[f(n)=sum_{n|d} mu(lfloor frac{d}{n}
floor) F(d)
]
开始解这道题，求的就是
[sum_{p in prime}f(p)
]
代入公式
[sum_{p in prime} sum_{p|d} mu(lfloor frac{d}{p}
floor) F(d)
]
枚举(lfloor frac{d}{p}
floor)
则
[sum_{p in prime} sum_{d=1}^{min(lfloor frac{n}{p}
floor, lfloor frac{m}{p}
floor)}mu(d)F(dp)
]
[sum_{p in prime} sum_{d=1}^{min(lfloor frac{n}{p}
floor, lfloor frac{m}{p}
floor)}mu(d) lfloor frac{n}{dp}
floor lfloor frac{m}{dp}
floor
]
设(dp)为(T)，然后枚举(T)
[sum_{T=1}^{min(n,m)} sum_{t|T,t in prime} mu (lfloor frac{T}{t}
floor) lfloor frac{n}{T}
floor lfloor frac{m}{T}
floor
]
[sum_{T=1}^{min(n,m)} lfloor frac{n}{T}
floor lfloor frac{m}{T}
floor (sum_{t|T} mu(lfloor frac{T}{t}
floor))
]
然后就可以做了。
代码
#include <bits/stdc++.h>
const int maxn = 1e7 + 10;
typedef long long ll;
inline void _swap(int& a,int& b) {
a ^= b ^= a ^= b;
}
int n, m, i, j, k, T;
int mu[maxn], f[maxn], s[maxn], pri[maxn];
bool vis[maxn];
inline void sieve(int n) {
mu[1] = 1;
for(int i = 2;i <= n;i++) {
if(!vis[i]) { pri[ ++pri[0] ] = i; mu[i] = -1; }
for(int j = 1;j <= pri[0] && i * pri[j] <= n;j++) {
vis[ i * pri[j] ] = 1;
if(i % pri[j])
mu[ i * pri[j] ] = -mu[i];
else {
mu[ i * pri[j] ] = 0;
break;
}
}
}
for(int i = 1;i <= pri[0];i++)
for(int j = 1;j * pri[i] <= n;j++)
f[ j * pri[i] ] += mu[j];
for(int i = 1;i <= n;i++)
s[i] = s[i - 1] + f[i];
}
int main() {
sieve(maxn - 10);
scanf("%d",&T);
while(T--) {
scanf("%d %d",&n,&m);
if(n > m)
_swap(n,m);
ll ans = 0;
for(int l = 1, r = 0;l <= n;l = r + 1) {
r = std::min(n / (n / l),m / (m / l));
ans += 1ll * (s[r] - s[l - 1]) * (ll)(n / l) * (ll)(m / l);
}
printf("%lld
",ans);
}
return 0;
}