GCD 题解
FangZeLi
·
·
题解
作为一道立志于打破套路的套路数论题，这题其实不是很难。
这道题有两种思路，一种很暴力，需要很扎实的基本功。一种很巧妙，需要很强的思维能力。
现在，我们以推出正解式子的过程，来讲解一下这道题目。
法一：暴力
\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]
这不正是我们最最熟悉的吗？
接下来有两种做法：
-
我会莫比乌斯反演！
套路的化简后半部分：
=\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！
-
我精通数论函数！
假如你做过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)
接下来，又有两种方法了：
-
我精通套路！
发现原式关于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})的时间内求出。
问题就在于如何预处理。
那，既然无法线性筛，我们就枚举每一个数$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的代码：
```cpp
#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的代码：
```cpp
#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$ 公式显示。