自动搬运

查看原文

来自洛谷，原作者为

	FangZeLi 宁愿重伤也不愿悲伤，让伤痕变成了我的徽章，刺在我心脏，永远不忘。

搬运于2025-08-24 22:16:56，当前版本为作者最后更新于2020-02-02 17:41:01，作者可能在搬运后再次修改，您可在原文处查看最新版

自动搬运只会搬运当前题目点赞数最高的题解，您可前往洛谷题解查看更多

以下是正文

---

作为一道立志于打破套路的套路数论题，这题其实不是很难。

这道题有两种思路，一种很暴力，需要很扎实的基本功。一种很巧妙，需要很强的思维能力。

现在，我们以推出正解式子的过程，来讲解一下这道题目。

法一：暴力

$$\sum_{i=1}^N\sum_{j=1}^N\sum_{p=1}^{\lfloor\frac{N}{j}\rfloor}\sum_{q=1}^{\lfloor\frac{N}{j}\rfloor}[\gcd(i,j)==1][\gcd(p,q)==1] $$

拿到这个式子，我们先进行最最基本的变形：

$$=\sum_{i=1}^N\sum_{j=1}^N[\gcd(i,j)==1]\sum_{p=1}^{\lfloor\frac{N}{j}\rfloor}\sum_{q=1}^{\lfloor\frac{N}{j}\rfloor}[\gcd(p,q)==1] $$

假如在这时候停笔枚举，就可以过第一个测试点。10 pts get！（然而复杂度我不会算233）

将后面拎出来看：

$$\sum_{p=1}^{\lfloor\frac{N}{j}\rfloor}\sum_{q=1}^{\lfloor\frac{N}{j}\rfloor}[\gcd(p,q)==1] $$

这不正是我们最最熟悉的吗？

接下来有两种做法：

1. 我会莫比乌斯反演！

   套路的化简后半部分：

$$=\sum_{p=1}^{\lfloor\frac{N}{j}\rfloor}\sum_{q=1}^{\lfloor\frac{N}{j}\rfloor}\sum_{d\mid\gcd(j,k)}\mu(d) $$$$=\sum_{d=1}^{\lfloor\frac{N}{j}\rfloor}\lfloor\frac{\lfloor\frac{N}{j}\rfloor}{d}\rfloor\lfloor\frac{\lfloor\frac{N}{j}\rfloor}{d}\rfloor\mu(d) $$

代回原式，得：

$$=\sum_{i=1}^N\sum_{j=1}^N[\gcd(i,j)==1]\sum_{d=1}^{\lfloor\frac{N}{j}\rfloor}\lfloor\frac{\lfloor\frac{N}{j}\rfloor}{d}\rfloor\lfloor\frac{\lfloor\frac{N}{j}\rfloor}{d}\rfloor\mu(d) $$

直接数论分块后按公式计算，复杂度$O(n^2\sqrt{n})$，20 pts get！

3. 我精通数论函数！

   假如你做过SDOI2008仪仗队，就应该会这种做法，化简后半部分：

$$=2\sum_{p=1}^{\lfloor\frac{N}{j}\rfloor}\sum_{q=1}^p[\gcd(p,q)==1]-1 $$

我们看这一部分

$$\sum_{q=1}^p[\gcd(p,q)==1]$$

这恰好是$\varphi$函数的定义式。

因此：

$$=2\sum_{p=1}^{\lfloor\frac{N}{j}\rfloor}\varphi(p)-1 $$

记$\varphi^s(n)=\sum_{i=1}^n\varphi(i)$，得：

$$=2\varphi^s(\lfloor\frac{N}{j}\rfloor)-1$$

回代入原式：

$$=\sum_{i=1}^N\sum_{j=1}^N[\gcd(i,j)==1](2\varphi^s(\lfloor\frac{N}{j}\rfloor)-1) $$

用线性筛预处理出$\varphi$及其前缀和，再按公式计算即可。

复杂度$O(n^2)$。20 pts get！

显然，我们绝不能满足于这样的复杂度，所以继续前进。

注意到后边的部分已经做的很好了，很难有继续优化的可能，在这样的情况下，前面的$[\gcd(i,d)==1]$就成了我们的眼中钉，肉中刺。

但这一部分却不能像之前一样优化，因为后面的部分依赖于此处枚举的$j$。

那没办法，莫比乌斯反演总行吧，大不了把$d$留着呗。

$$=\sum_{g=1}^N\sum_{i=1}^{\lfloor\frac{N}{g}\rfloor}\sum_{j=1}^{\lfloor\frac{N}{g}\rfloor}\mu(g)(2\varphi^s(\lfloor\frac{N}{jg}\rfloor)-1) $$$$=\sum_{g=1}^N\lfloor\frac{N}{g}\rfloor\mu(g)\sum_{j=1}^{\lfloor\frac{N}{g}\rfloor}(2\varphi^s(\lfloor\frac{N}{jg}\rfloor)-1) $$

接下来，又有两种方法了：

1. 我精通套路！

   发现原式关于$g$和$j$的枚举其实就是在枚举乘积不大于$n$的二元组，于是令$T=gj$，得：

$$=\sum_{T=1}^N(2\varphi^s(\lfloor\frac{N}{T}\rfloor)-1)\sum_{g\mid T}\lfloor\frac{N}{g}\rfloor\mu(g) $$

假如我们令$f(T)=\sum_{d\mid T}\lfloor\frac{N}{d}\rfloor\mu(d)$，$h(T)=2\varphi^s(T)-1$，原式就化为：

$$=\sum_{T=1}^Nh(\lfloor\frac{N}{T}\rfloor)f(T)$$

这是一个除法分块的形式，可以在预处理后$O(\sqrt{n})$的时间内求出。

问题就在于如何预处理。

$g$是很好求的，而$f$很棘手，$f$甚至不是一个积性函数。

那，既然无法线性筛，我们就枚举每一个数$p$及其倍数，暴筛吧。

复杂度是调和级数$O(n\log n)$，60 pts get！

3. 我精通数论分块！

   我们观察原式，将其直接分为两个部分：

$$=\sum_{g=1}^N\mu(g)\times\lfloor\frac{N}{g}\rfloor\sum_{j=1}^{\lfloor\frac{N}{g}\rfloor}(2\varphi^s(\lfloor\frac{N}{jg}\rfloor)-1) $$

我们重新定义函数$f$：

$$ f(T)=T\sum_{i=1}^T(2\varphi^s(\lfloor\frac{T}{i}\rfloor)-1) $$

显然，这个函数可以数论分块$O(\sqrt{n})$求。

代回原式：

$$\sum_{g=1}^N\mu(g)f(\lfloor\frac{N}{g}\rfloor)$$

这个式子也是可以数论分块$O(\sqrt{n})$求的。

但我们再一次遇到了预处理的难题，这次如果再暴力，复杂度就是$O(n\sqrt{n})$的，还不如上一种方法。

答案是：我们根本不需要预处理！

再对第二个式子数论分块时，我们求出所需的$f$值即可。

这样做，复杂度是$O(n^\frac{3}{4})$的，考虑到预处理的$O(n)$复杂度，总复杂度依然是$O(n)$的。80 pts get！

如何改进复杂度呢？

预处理使用杜教筛，即可使复杂度瓶颈变为计算的过程。新复杂度是$O(n^\frac{3}{4})$的。90 pts get！

文末会有该档部分分的代码。

法二：思维

不进行任何变形：

$$\sum_{i=1}^N\sum_{j=1}^N\sum_{p=1}^{\lfloor\frac{N}{j}\rfloor}\sum_{q=1}^{\lfloor\frac{N}{j}\rfloor}[\gcd(i,j)==1][\gcd(p,q)==1] $$

仔细观察这个式子，你会发现，关于$j$的枚举，很像是在枚举$p$，$q$的$\gcd$。

我们来试一试。

$$=\sum_{i=1}^N\sum_{p=1}^N\sum_{q=1}^N[\gcd(i,\gcd(p,q))==1] $$$$=\sum_{i=1}^N\sum_{p=1}^N\sum_{q=1}^N[\gcd(i,p,q)==1] $$

写到这里，你可以停笔$O(n^3)$的去枚举（10 pts get），但其实你离另一种正解只有两步之遥。

$$=\sum_{d=1}^N\sum_{i=1}^{\lfloor\frac{N}{d}\rfloor}\sum_{p=1}^{\lfloor\frac{N}{d}\rfloor}\sum_{q=1}^{\lfloor\frac{N}{d}\rfloor}\mu(d) $$$$=\sum_{d=1}^N\mu(d)\lfloor\frac{N}{d}\rfloor\lfloor\frac{N}{d}\rfloor\lfloor\frac{N}{d}\rfloor $$

这个式子可以在$O(n)$的预处理后，$O(\sqrt{n})$的求得，总复杂度$O(n)$。80 pts get！

和上一种方法一样，我们使用杜教筛优化。瓶颈依旧是预处理，但复杂度降至$O(n^\frac{2}{3})$，100 pts get！

文末会有该标程的代码。

接下来是代码时间：

90 pts的代码：

#define _CRT_SECURE_NO_WARNINGS

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<map>

using namespace std;

#define _N 1000001
#define _P 998244353
#define _INV 499122177

int n;

int pricnt = 0;
bool visited[_N];
int prime[_N];
int mu[_N];
int phi[_N];

map<int, int> ansmu, ansphi;

long long ans = 0;

void sieve()
{
	int range = n + 1 < _N ? n + 1 : _N;
	visited[1] = true;
	mu[1] = phi[1] = 1;
	for (int i = 1; i < range; i++)
	{
		if (!visited[i])
		{
			prime[++pricnt] = i;
			mu[i] = _P - 1;
			phi[i] = i - 1;
		}
		for (int j = 1; j <= pricnt; j++)
		{
			if (i * prime[j] >= range)
			{
				break;
			}
			visited[i * prime[j]] = true;
			if (i % prime[j])
			{
				mu[i * prime[j]] = (_P - mu[i]) % _P;
				phi[i * prime[j]] = phi[i] * phi[prime[j]] % _P;
			}
			else
			{
				phi[i * prime[j]] = phi[i] * prime[j] % _P;
				break;
			}
		}
	}
	for (int i = 1; i < range; i++)
	{
		mu[i] = (mu[i - 1] + mu[i]) % _P;
		phi[i] = (phi[i - 1] + phi[i]) % _P;
	}
}

int getmu(int x)
{
	if (x < _N)
	{
		return mu[x];
	}
	if (ansmu[x])
	{
		return ansmu[x];
	}
	long long res = 1;
	for (int l = 2, r; l <= x; l = r + 1)
	{
		r = x / (x / l);
		res = (res - ((long long)r - l + 1) * getmu(x / l) % _P + _P) % _P;
	}
	return ansmu[x] = res;
}

int getphi(int x)
{
	if (x < _N)
	{
		return phi[x];
	}
	if (ansphi[x])
	{
		return ansphi[x];
	}
	long long res = (long long)x * ((long long)x + 1) % _P * _INV % _P;
	for (int l = 2, r; l <= x; l = r + 1)
	{
		r = x / (x / l);
		res = (res - ((long long)r - l + 1) * getphi(x / l) % _P + _P) % _P;
	}
	return ansphi[x] = res;
}

long long f(int n)
{
	long long res = 0;
	for (int l = 1, r = 0; l <= n; l = r + 1)
	{
		r = n / (n / l);
		res = (res + (r - l + 1) * (2 * getphi(n / l) % _P - 1) % _P) % _P;
	}
	res = res * n % _P;
	return res;
}

int main()
{
	scanf("%d", &n);
	sieve();
	for (int l = 1, r = 0; l <= n; l = r + 1)
	{
		r = n / (n / l);
		ans = (ans + ((long long)getmu(r) - getmu(l - 1) + _P) % _P * f(n / l) % _P) % _P;
	}
	printf("%lld\n", ans);
	return 0;
}

100 pts的代码：

#define _CRT_SECURE_NO_WARNINGS

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<map>

using namespace std;

#define _N 1000001
#define _P 998244353

int n;

int pricnt = 0;
bool visited[_N];
int prime[_N];
int mu[_N];

map<int, int> ansmu;

long long ans = 0;

void sieve()
{
	int range = n + 1 < _N ? n + 1 : _N;
	visited[1] = true;
	mu[1] = 1;
	for (int i = 1; i < range; i++)
	{
		if (!visited[i])
		{
			prime[++pricnt] = i;
			mu[i] = _P - 1;
		}
		for (int j = 1; j <= pricnt; j++)
		{
			if (i * prime[j] >= range)
			{
				break;
			}
			visited[i * prime[j]] = true;
			if (i % prime[j])
			{
				mu[i * prime[j]] = (_P - mu[i]) % _P;
			}
			else
			{
				break;
			}
		}
	}
	for (int i = 1; i < range; i++)
	{
		mu[i] = (mu[i - 1] + mu[i]) % _P;
	}
}

int getmu(int x)
{
	if (x < _N)
	{
		return mu[x];
	}
	if (ansmu[x])
	{
		return ansmu[x];
	}
	long long res = 1;
	for (int l = 2, r; l <= x; l = r + 1)
	{
		r = x / (x / l);
		res = (res - ((long long)r - l + 1) * getmu(x / l) % _P + _P) % _P;
	}
	return ansmu[x] = res;
}

int main()
{
	scanf("%d", &n);
	sieve();
	for (int l = 1, r = 0; l <= n; l = r + 1)
	{
		r = n / (n / l);
		ans = (ans + ((long long)getmu(r) - getmu(l - 1) + _P) % _P * ((long long)n / l) % _P * ((long long)n / l) % _P * ((long long)n / l) % _P) % _P;
	}
	printf("%lld\n", ans);
	return 0;
}

Update At 2024/11/19: 修复了文中的 $\LaTeX$ 公式显示。