约数个数和

题目链接：luogu P3327

题目大意

给你 n,m 要你求 ∑i=1~n∑j=1~m d(ij)，d 是约数个数。
多组数据。

思路

首先由一个比较神奇的东西：
\(d(ij)=\sum\limits_{x|i}\sum\limits_{y|j}[gcd(x,y)=1]\)

这里证明一下：
先看我们可以只看一个质数的部分。
假设它的因子的出现次数是 \(x\)，\(i\) 中有 \(a\) 个，\(j\) 中有 \(b\) 个，那 \(ij\) 中就是 \(a+b\) 个。

如果 \(x\leqslant a\)，那它在 \(i\) 中会出现被数到。
那如果 \(x>a\)，我们可以在 \(j\) 中数一个 \(x-a\) 就表示满了之后再我那个这边的个数。

所以是对的。

然后就来化式子：
\(\sum\limits_{i=1}^n\sum\limits_{j=1}^md(ij)\)
\(=\sum\limits_{i=1}^n\sum\limits_{j=1}^m\sum\limits_{x|i}\sum\limits_{y|j}[gcd(x,y)=1]\)
\(=\sum\limits_{x=1}^n\sum\limits_{y=1}^m\sum\limits_{i=1}^{\frac{n}{x}}\sum\limits_{j=1}^{\frac{m} {y}}[gcd(x,y)=1]\)
\(=\sum\limits_{x=1}^n\sum\limits_{y=1}^m\left\lfloor\dfrac{n}{x}\right\rfloor\left\lfloor\dfrac{m} {y}\right\rfloor[gcd(x,y)=1]\)
\(=\sum\limits_{i=1}^n\sum\limits_{j=1}^m\left\lfloor\dfrac{n}{i}\right\rfloor\left\lfloor\dfrac{m} {j}\right\rfloor[gcd(i,j)=1]\)

然后莫反。
\(f_x=\sum\limits_{i=1}^n\sum\limits_{j=1}^m\left\lfloor\dfrac{n}{i}\right\rfloor\left\lfloor\dfrac{m} {j}\right\rfloor[gcd(i,j)=x]\)

\(g_i=\sum\limits_{x|i}f_i=\sum\limits_{i=1}^n\sum\limits_{j=1}^m\left\lfloor\dfrac{n}{i}\right\rfloor\left\lfloor\dfrac{m} {j}\right\rfloor[x|gcd(i,j)]\)
\(=\sum\limits_{i=1}^n\sum\limits_{j=1}^m\left\lfloor\dfrac{n}{i}\right\rfloor\left\lfloor\dfrac{m} {j}\right\rfloor\sum\limits_{x|i,x|j}1\)
\(=\sum\limits_{x}\sum\limits_{i=1}^{\left\lfloor\frac{n}{x}\right\rfloor}\sum\limits_{j=1}^{\left\lfloor\frac{m}{x}\right\rfloor}\left\lfloor\dfrac{n}{ix}\right\rfloor\left\lfloor\dfrac{m} {jx}\right\rfloor\)

\(f_i=\sum\limits_{i|x}\mu(\dfrac{x}{i})g_x\)
\(f_1=\sum\limits_{1|x}\mu(\dfrac{x}{1})g_x=\sum\limits_{i=1}^{\min(n,m)}\mu(i)g_i\)

那询问我们可以预处理 \(\mu\) 的前缀和然后用数论分块求 \(g_i\)。
（用数论分块预处理 \(\sum\limits_{i=1}^x\left\lfloor\dfrac{x}{i}\right\rfloor\)）
（然后用数论分块搞 \(\sum\limits_{x}\) 的部分）

然后每次询问你就跑一次数论分块 \(O(\sqrt{n})\)，所以询问的复杂度是 \(O(T\sqrt{n})\) 可以的。

代码

#include<cstdio>
#include<iostream>
#define ll long long

using namespace std;

int T, n, m, prime[50001];
ll f[50001], ans, miu[50001];
bool np[50001];

int main() {
	miu[1] = 1;
	for (int i = 2; i <= 50000; i++) {
		if (!np[i]) prime[++prime[0]] = i, miu[i] = -1;
		for (int j = 1; j <= prime[0] && i * prime[j] <= 50000; j++) {
			np[i * prime[j]] = 1;
			if (i % prime[j] == 0) {miu[i * prime[j]] = 0; break;}
				else miu[i * prime[j]] = -miu[i];
		}
	}
	for (int i = 1; i <= 50000; i++) miu[i] += miu[i - 1];
	
	for (int i = 1; i <= 50000; i++) {
		int r;
		for (int l = 1; l <= i; l = r + 1) {
			r = i / (i / l);
			f[i] += 1ll * (r - l + 1) * (i / l);
		}
	}
	
	scanf("%d", &T);
	while (T--) {
		scanf("%d %d", &n, &m);
		ans = 0;
		
		if (n > m) swap(n, m);
		
		int r;
		for (int l = 1; l <= n; l = r + 1) {
			r = min(n / (n / l), m / (m / l));
			ans += (miu[r] - miu[l - 1]) * f[n / l] * f[m / l];
		}
		
		printf("%lld\n", ans);
	}
	
	return 0;
}