【题目描述】
给定 $n. m, k$,求: $$ \sum_{i = 1}^n \sum_{j = 1}^m [\gcd(i, j) = k] $$ $n, m, k \leq 5 \times 10^4$
有 $T$ 组询问,$T \leq 5 \times 10^4$。
【题目链接】
BZOJ 1011 Zap
【解题思路】
考虑到 $\gcd(i, j) = k$ 等价于 $\gcd(\frac i k, \frac j k) = 1$,那么我们先将 $n, m$ 分别除以 $k$,为方便,将除后的结果仍记作 $n, m$ ,那么我们要求的就是:
$$ \text{ans}(n, m) = \sum_{i = 1}^n\sum_{j = 1}^m[\gcd(i, j) = 1] $$ 通用思路,将$[\gcd(i, j) = 1]$ 通过 $\sum\limits_{d \mid x}^x \mu(d) = [x = 1]$ 即$\mu \times u = e$替换: $$ \text{ans}(n, m) =\sum_{i = 1}^n\sum_{j = 1}^m\sum_{d | \gcd(i, j)} \mu(d) $$ 接下来更改求和指标,改为外层枚举 $d$,并只有在 $d \mid \gcd(i, j)$ 时将$\mu(d)$计入答案,即: $$ ans(n, m) = \sum_{d = 1}^{min(n, m)}\mu(d) \sum_{i = 1}^n \sum_{j = 1}^m [d \mid \gcd(i, j)] $$ 考虑后面那两个 $\sum$,其实它们计算的是 $d$ 成为 $\gcd(i, j)$ 的次数,这个其实是可以快速计算的。
我们知道 $d \mid \gcd(i, j)$,当且仅当 $d \mid i \land d \mid j$ ,而这两个条件是独立的,我们可以分别计算它们各自成立的次数,然后将其相乘,即: $$ \begin{align*} \sum_{i = 1}^n \sum_{j = 1}^m[d \mid \gcd(i, j)] &= \sum_{i = 1}^n \sum_{j = 1}^m[d \mid i \land d \mid j] \\ &= \sum_{i = 1}^n [d \mid i] \cdot \sum_{j = 1}^m [d \mid j] \end{align*} $$ 而这样一个式子 $\sum\limits_{i = 1}^n[d \mid i]$ ,其含义是 $1 \sim n$ 中,$d$ 的倍数的个数,那么其值也就很明了了: $$ \sum_{i = 1}^n[d \mid i] = \left \lfloor \frac n d \right \rfloor $$ 至此,我们的式子成了这个样子: $$ \text {ans}(n, m) = \sum_{d}^{\min(n, m)}\mu(d) \cdot \left \lfloor \frac n d \right \rfloor \cdot \left \lfloor \frac m d \right \rfloor $$ 由于 $\left \lfloor \frac n d \right \rfloor$ 至多有 $2\sqrt n$ 种取值,$\left \lfloor \frac m d \right \rfloor$同理,所以我们可以枚举使其各自取值相同的 $d$ 的区间,然后并在一起计算,一起计算过程的进行需要 $\mu$ 的区间和,可以通过预处理 $\mu$ 的前缀和来计算 。
总结:
利用 $\mu \times u = e$ 对形如 $[\gcd(i, j) = 1]$ 的式子进行代换。
更改求和指标,枚举约数。
考虑每个约数满足条件的次数。
最后一般是利用下取整的性质,成块统计答案。
【AC 代码】
#include <cstdio>
#include <algorithm>
typedef long long int64;
const int MAX_N = 50000;
bool isNotPrime[MAX_N + 1];
int primes[MAX_N], m;
int64 mu[MAX_N + 1];
inline void seive(int n = MAX_N){
isNotPrime[1] = true, mu[1] = 1;
for(int i = 2; i <= n; i++){
if(!isNotPrime[i]){
primes[m++] = i;
mu[i] = -1;
}
for(int j = 0, p, x; j < m; j++){
if((x = i * (p = primes[j])) > n) break;
isNotPrime[x] = true;
if(i % p == 0){
mu[x] = 0;
break;
} else{
mu[x] = -mu[i];
}
}
}
for(int i = 1; i <= n; i++) mu[i] += mu[i - 1];
}
inline int64 solve(int n, int m, int k){
n /= k, m /= k;
if(n > m) std::swap(n, m);
int64 ans = 0;
for(int l = 1, r; l <= n; l = r + 1){
r = std::min(n / (n / l), m / (m / l));
ans += (mu[r] - mu[l - 1]) * (int64)(n / l) * (m / l);
}
return ans;
}
int main(){
seive();
int T_T;
scanf("%d", &T_T);
while(T_T--){
int n, m, k;
scanf("%d%d%d", &n, &m, &k);
printf("%lld\n", solve(n, m, k));
}
return 0;
}