【题解】Luogu-P4240 毒瘤之神的考验

P4240 毒瘤之神的考验

Description

• 多测，\(t\) 组数据。

• 每次给定两个整数 \(n, m\)，请求出

  \[\left[ \sum_{i = 1}^n \sum_{j = 1}^m \varphi(ij) \right] \bmod 998244353 \]

• \(1\le t\le 10^4\)，\(1\le n, m\le 10^5\)。

Solution

引理：

\[\varphi(ij) = \dfrac{\varphi(i) \varphi(j) \gcd(i, j)}{\varphi(\gcd(i, j))} \]

证明：

\[ \begin{aligned} \varphi(i) \varphi(j) \gcd(i, j) & = i \prod_{p_1\in \mathbb{P}, p_1\mid i} \dfrac{p_1 - 1}{p_1} \cdot j \prod_{p_2\in \mathbb{P}, p_2\mid j} \dfrac{p_2 - 1}{p_2} \cdot \gcd(i, j) \\ & = ij \gcd(i, j) \left(\prod_{p_1\in \mathbb{P}, p_1\mid i} \dfrac{p_1 - 1}{p_1} \prod_{p_2\in \mathbb{P}, p_2\mid j} \dfrac{p_2 - 1}{p_@} \right) \end{aligned} \]

你会发现所有 \(\gcd(i, j)\) 的质因数都被算了 \(2\) 次。

\[ \begin{aligned} & = \left(ij\prod_{p_1\in \mathbb{P}, p_1\mid ij} \dfrac{p_1 - 1}{p_1} \right) \cdot \left(\gcd(i, j) \prod_{p_2\in \mathbb{P}, p_2\mid \gcd(i, j)} \dfrac{p_2 - 1}{p_2} \right) \\ & = \varphi(ij) \varphi(\gcd(i, j)) \end{aligned} \]

把 \(\varphi(\gcd(i, j))\) 除过去即可。

不妨设 \(n\le m\)，根据引理，有

\[\begin{aligned} ans & = \sum_{i = 1}^n \sum_{j = 1}^m \dfrac{\varphi(i) \varphi(j) \gcd(i, j)}{\varphi(\gcd(i, j))} \\ & = \sum_{d = 1}^n \sum_{i = 1}^n \sum_{j = 1}^m \dfrac{\varphi(i) \varphi(j) d}{\varphi(d)} [\gcd(i, j) = d] \\ & = \sum_{d = 1}^n \dfrac{d}{\varphi(d)} \sum_{i = 1}^{\left\lfloor\frac{n}{d}\right\rfloor} \sum_{j = 1}^{\left\lfloor\frac{m}{d}\right\rfloor} \varphi(id) \varphi(jd) [\gcd(i, j) = 1] \\ & = \sum_{d = 1}^n \dfrac{d}{\varphi(d)} \sum_{i = 1}^{\left\lfloor\frac{n}{d}\right\rfloor} \varphi(id) \sum_{j = 1}^{\left\lfloor\frac{m}{d}\right\rfloor} \varphi(jd) [\gcd(i, j) = 1] \\ & = \sum_{d = 1}^n \dfrac{d}{\varphi(d)} \sum_{i = 1}^{\left\lfloor\frac{n}{d}\right\rfloor} \varphi(id) \sum_{j = 1}^{\left\lfloor\frac{m}{d}\right\rfloor} \varphi(jd) \sum_{k\mid \gcd(i, j)} \mu(k) \\ & = \sum_{d = 1}^n \dfrac{d}{\varphi(d)} \sum_{k = 1}^{\left\lfloor\frac{n}{d}\right\rfloor} \mu(k) \sum_{i = 1}^{\left\lfloor\frac{n}{dk}\right\rfloor} \varphi(idk) \sum_{j = 1}^{\left\lfloor\frac{m}{dk}\right\rfloor} \varphi(jdk) \\ & = \sum_{T = 1}^n \sum_{d\mid T} \dfrac{d}{\varphi(d)} \mu\left(\dfrac{T}{d}\right) \sum_{i = 1}^{\left\lfloor\frac{n}{T}\right\rfloor} \varphi(iT) \sum_{j = 1}^{\left\lfloor\frac{m}{T}\right\rfloor} \varphi(jT) \end{aligned} \]

其中

\[f(n) = \sum_{d\mid n} \dfrac{d\, \mu\left(\dfrac{n}{d}\right)}{\varphi(d)} \]

不是积性函数（它甚至不一定是整数），直接预处理逆元 \(\Theta(n\ln n)\) 计算。

对于

\[\sum_{i = 1}^{\left\lfloor\frac{n}{T}\right\rfloor} \varphi(iT) \sum_{j = 1}^{\left\lfloor\frac{m}{T}\right\rfloor} \varphi(jT) \]

发现这东西带两个参数，令

\[g(k, n) = \sum_{i = 1}^n \varphi(ik) \]

你有递推式

\[g(k, n) = g(k, n - 1) + \varphi(nk) \]

而此处的 \(nk\le n\textsf{(题目中的 } n \textsf{)}\)，也可以 \(\Theta(n\ln n)\) 预处理，\(g\) 数组用 vector 存。

代回去

\[ans = \sum_{T = 1}^n f(T) g\left(T, \left\lfloor\dfrac{n}{T}\right\rfloor \right) g\left(T, \left\lfloor\dfrac{m}{T}\right\rfloor \right) \]

这个 \(g\) 带一个 \(T\)，无法整除分块，只能暴力算。

此时是 \(\Omicron(n\ln n)\) 的预处理，\(\Omicron(n)\) 的单次回答，过不去。

---

\[h(n, m, l) = \sum_{T = 1}^l f(T) g\left(T, n \right) g\left(T, m \right) \]

这样在整除分块中

\[h\left(\left\lfloor\dfrac{n}{l}\right\rfloor, \left\lfloor\dfrac{m}{l}\right\rfloor, r \right) - h\left(\left\lfloor\dfrac{n}{l}\right\rfloor, \left\lfloor\dfrac{m}{l}\right\rfloor, l - 1 \right) \]

又有递推式

\[h(n, m, l) = h(n, m, l - 1) + f(l) g(l, n) g(l, m) \]

有 \(\max(n, m) l\le n \textsf{(题目中的 } n \textsf{)}\)。

\(h\) 数组同样用 vector 存。

此时是

\[\begin{aligned} \Omicron\left(n \sum_{i = 1}^n \dfrac{n}{i} \right) & = \Omicron\left(n \cdot \int_1^n \dfrac{n}{x} \, dx \right) \\ & = \Omicron(n^2 \ln n) \end{aligned} \]

的预处理，\(\Omicron(\sqrt{n})\) 的单次回答，预处理会 \(\text{MLE} + \text{TLE}\)。

---

既然两种方法都不行，那么自然想到 用根号分治来平衡时间。

设阈值为 \(k\)，表示对于 \(\le k\) 的部分预处理。

于是当 \(\left\lfloor\dfrac{m}{l}\right\rfloor \le k\)，即 \(l > \left\lfloor\dfrac{m}{k}\right\rfloor\) 的部分用已经预处理过的 \(h\) 整除分块算，当 \(l \le \left\lfloor\dfrac{m}{k}\right\rfloor\) 时直接暴力算。

• 暴力：\(\Omicron(n\ln n) + \Omicron\left(\dfrac{tn}{k} \right)\)
• 分块：\(\Omicron(nk \ln k) + \Omicron(t\sqrt{n})\)

总时间复杂度为 \(\Omicron\left(nk \ln k + t\left(\dfrac{n}{k} + \sqrt{n} \right) \right)\)，提个公因数就是 \(k\ln k + \dfrac{t}{k}\) 最小。

写个程序去算

	int t = 1e4;
	int ans = 0x3f3f3f3f, pos = 0;
	for (int k = 1; k <= t; k++)
	{
		int res = k * log(k) + t / k;
		if (res < ans)
		{
			ans = res, pos = k;
		}
	}
	printf("%d %d\n", ans, pos);

输出

392 47

所以取 \(k = 47\)。

空间复杂度为 \(\Omicron(nk\ln k)\)。

然而实测 \(k\) 取 \(60\sim 70\) 更快（

Code

// 18 = 9 + 9 = 18.
#include 
#include 
#include 
#define Debug(x) cout << #x << "=" << x << endl
typedef long long ll;
using namespace std;

const int MAXN = 1e5 + 5;
const int N = 1e5;
const int MAXK = 47 + 5;
const int K = 47;
const int MOD = 998244353;

typedef int arr[MAXN];

int add(int a, int b)
{
	return (a + b) % MOD;
}

int sub(int a, int b)
{
	return (a - b + MOD) % MOD;
}

int mul(int a, int b)
{
	return (ll)a * b % MOD;
}

int qpow(int a, int b)
{
	int base = a, ans = 1;
	while (b)
	{
		if (b & 1)
		{
			ans = (ll)ans * base % MOD;
		}
		base = (ll)base * base % MOD;
		b >>= 1;
	}
	return ans;
}

int inv(int a)
{
	return qpow(a, MOD - 2);
}

arr p, mu, phi, phi_pro, phi_pro_inv, phi_inv, f;
bool vis[MAXN];
vector g[MAXN], h[MAXK][MAXK];

void pre()
{
	mu[1] = phi[1] = 1;
	for (int i = 2; i <= N; i++)
	{
		if (!vis[i])
		{
			p[++p[0]] = i;
			mu[i] = MOD - 1;
			phi[i] = i - 1;
		}
		for (int j = 1; j <= p[0] && i * p[j] <= N; j++)
		{
			vis[i * p[j]] = true;
			if (i % p[j] == 0)
			{
				mu[i * p[j]] = 0;
				phi[i * p[j]] = phi[i] * p[j];
				break;
			}
			mu[i * p[j]] = MOD - mu[i];
			phi[i * p[j]] = phi[i] * phi[p[j]];
		}
	}
	
	phi_pro[0] = phi_inv[1] = 1;
	for (int i = 1; i <= N; i++)
	{
		phi_pro[i] = mul(phi_pro[i - 1], phi[i]);
	}
	phi_pro_inv[N] = inv(phi_pro[N]);
	for (int i = N - 1; i >= 1; i--)
	{
		phi_pro_inv[i] = mul(phi_pro_inv[i + 1], phi[i + 1]);
		phi_inv[i + 1] = mul(phi_pro_inv[i + 1], phi_pro[i]);
	}
	
	for (int i = 1; i <= N; i++)
	{
		for (int j = 1; i * j <= N; j++)
		{
			f[i * j] = add(f[i * j], mul(mul(i, mu[j]), phi_inv[i]));
		}
	}
	
	for (int k = 1; k <= N; k++)
	{
		g[k].resize(N / k + 5);
		for (int n = 1; n * k <= N; n++)
		{
			g[k][n] = add(g[k][n - 1], phi[n * k]);
		}
	}
	
	for (int n = 1; n <= K; n++)
	{
		for (int m = n; m <= K; m++)
		{
			h[n][m].resize(N / m + 5);
			for (int l = 1; m * l <= N; l++)
			{
				h[n][m][l] = add(h[n][m][l - 1], mul(f[l], mul(g[l][n], g[l][m])));
			}
		}
	}
}

int block(int n, int m)
{
	int res = 0;
	for (int i = 1; i <= m / K; i++)
	{
		res = add(res, mul(f[i], mul(g[i][n / i], g[i][m / i])));
	}
	for (int l = m / K + 1, r; l <= n; l = r + 1)
	{
		int k1 = n / l, k2 = m / l;
		r = min(n / k1, m / k2);
		res = add(res, sub(h[k1][k2][r], h[k1][k2][l - 1]));
	}
	return res;
}

int main()
{
	pre();
	int t;
	scanf("%d", &t);
	while (t--)
	{
		int n, m;
		scanf("%d%d", &n, &m);
		if (n > m)
		{
			swap(n, m);
		}
		printf("%d\n", block(n, m));
	}
	return 0;
}