P3455 [POI2007]ZAP-Queries

$题意$
题目链接

$求{\sum }_{i=1}^a{\sum }_{j=1}^b[gcd(i,j)=d]$


$题解$


$莫比乌斯函数有个性质$
$$\sum _{d|n}\mu (d)=[n=1]$$
$我们使用gcd(i,j)代替n得到如下的式子$
$$\sum _{d|gcd(i,j)}\mu (d)=[gcd(i,j)=1]$$
$接下来我们来推式子,我们假设a<b$

$\begin{array}{rl}原式& =\sum _{i=1}^a\sum _{j=1}^b[gcd(i,j)=d]\\ & =\sum _{i=1}^n\sum _{j=1}^m[gcd(i,j)=1](令n=a/d,m=b/d)\\ & =\sum _{i=1}^n\sum _{j=1}^m\cdot \sum _{k|gcd(i,j)}\mu (k)\\ & =\sum _{k=1}^n\mu (k)\cdot \sum _{i=1}^n\sum _{j=1}^m[k|gcd(i,j)]\\ & =\sum _{k=1}^n\mu (k)\cdot \lfloor \frac{n}{k}\rfloor \lfloor \frac{m}{k}\rfloor \\ & =\boxed{\sum _{k=1}^{a/d}\mu (k)\cdot \lfloor \frac{a}{dk}\rfloor \lfloor \frac{b}{dk}\rfloor }\end{array}$

$$\begin{gathered} 那么我们可以在O(N)的时间复杂度内来计算出这个式子,但是由于有多组数据， \\ 因此我们需要整除分块 \end{gathered}$$¹

$上面的结论，我们可以换一种方式来推导，这里使用了莫比乌斯反演$²
$我们设f(d)=gcd(i,j)为d的个数$
$$f(d)=\sum _{i=1}^a\sum _{j=1}^b[gcd(i,j)=d]$$
$再设g(d)为gcd(i,j)为d的倍数的个数$
$$g(d)=\sum _{i=1}^{min(a/d,b/d)}f(di)$$
$同时，根据定义$
$$g(d)=\lfloor \frac{a}{d}\rfloor \lfloor \frac{b}{d}\rfloor$$
$根据莫比乌斯反演的第二种形式，可以从g推出f$
$$f(d)=\sum _{i=1}^{min(a/d,b/d)}\mu (i)g(di)$$
$设a<b,整理式子得到$
$$f(d)=\boxed{\sum _{i=1}^{a/d}\ast \mu (i)\cdot \lfloor \frac{a}{di}\rfloor \lfloor \frac{b}{di}\rfloor }$$

$用i替换为k后，和上面的式子是相同的$

$代码$

#include <bits/stdc++.h>
using namespace std;

typedef long long LL;

const int N = 50020;
int n, mu[N+50], p[N+50];
LL sum[N+50];
int tot;
bool used[N+50];

//预处理出素数 莫比乌斯函数
void init() {
	mu[1] = 1;
	for (int i = 2; i <= N; i++) {
		if (!used[i]) {
			mu[i] = -1;
			p[++tot] = i;
		}
		for (int j = 1; j <= tot && p[j] * i <= N; j++) {
			used[p[j] * i] = true;
			if (i % p[j] == 0) {
				mu[i * p[j]] = 0;
				break;
			}
			mu[i * p[j]] = mu[i] * (-1);
		}
	}
	for (int i = 1; i <= N; i++) sum[i] = mu[i] + sum[i - 1];
	
}

LL cal(int a, int b, int d) {
	if (a > b) swap(a, b);
	a /= d;
	b /= d;
	LL ans = 0;
	//暴力计算 会TLE
	// for (int i = 1; i <= a; i++) {
		// ans += 1LL * mu[i] * (a / i) * (b / i);
	// }
	
	//整除分块
	for (int i = 1, j; i <= a; i = j + 1) {
		//这里是两个整除分块，取最小值，可以拿(n = 10, m = 20 从i = 6开始想一下)
		j = min(a / (a / i), b / (b / i));
		ans += 1LL * (sum[j] - sum[i - 1]) * (a / i) * (b / i);
	}
	return ans;
}

int main() {
	int T, a, b, d;
	scanf("%d", &T);
	init();
	while (T--) {
		scanf("%d%d%d", &a, &b, &d);
		printf("%lld\n", cal(a, b, d));
	}
	return 0;
}

---

1. 整除分块 ↩︎

2. 莫比乌斯反演 ↩︎