Description:
共\(T \le 5 \times 10^4\)组询问, 每组询问给定\(n\)和\(m\), 请你求出
\[ \sum_{i = 1}^n \sum_{j = 1}^m \sigma_0 (ij) \]
Solution:
先给出一个结论:
\[ \sigma_0(ij) = \sum_{a | i} \sum_{b | j} [\gcd(a, b) = 1] \]
证明: 我们令\(i = p_1^{a_1} p_2^{a_2} \cdots\), \(j = p_1^{b_1} p_2^{b_2} \cdots\), \(d | ij\)且\(d = p_1^{c_1} p_2^{c_2} \cdots\), 则\(c_n \le a_n + b_n\).
考虑如何不重复地统计每一个\(d\): 令\(c_n = A_n + B_n\), 其中\(A_n\)和\(B_n\)分别为\(i\)和\(j\)对\(c_n\)的贡献, 则我们要求
\[ \begin{cases} B_n = 0 & A_n < a_n \\ B_n \ge 0 & A_n = a_n \end{cases} \]
这样一来, \(c_n\)的表示形式就变成唯一的了, 因而不会被重复统计. 我们再考虑如何统计这样的\(A_n\)和\(B_n\): 我们令\(A_n' = a_n - A_n\), 则约束条件变为
\[ \begin{cases} B_n = 0 & A_n' \ne 0 \\ B_n \ge 0 & A_n' = 0 \end{cases} \]
等价于\(\gcd(A_n', B_n) = 1\).
因此得证.
好吧, 假如看不懂上面的这一些证明, 就这么想吧: \(i\)表示\(a\)中不取多少, \(j\)表示\(b\)中取多少, 只要保证\(\gce(a, b) = 1\), 即不会重复统计.
因此我们考虑原题的式子
\[ \begin{aligned} \sum_{i = 1}^n \sum_{j = 1}^m \sigma_0(ij) &= \sum_{i = 1}^n \sum_{j = 1}^m \sum_{a | i} \sum_{b | j} [\gcd(a, b) = 1] \\ &= \sum_{i = 1}^n \sum_{j = 1}^m \sum_{a | i} \sum_{b | j} \sum_{d | \gcd(a, b)} \mu(d) \\ &= \sum_{i = 1}^n \sum_{j = 1}^m \sum_{d | \gcd(i, j)} \mu(d) \sigma_0(\frac i d) \sigma_0(\frac j d) \\ &= \sum_{d = 1}^n \sum_{i = 1}^{\lfloor \frac n d \rfloor} \sum_{j = 1}^{\lfloor \frac m d \rfloor} \mu(d) \sigma_0(i) \sigma_1(j) \\ &= \sum_{d = 1}^n \mu(d) \sum_{i = 1}^{\lfloor \frac n d \rfloor} \sigma_0(i) \sum_{j = 1}^{\lfloor \frac m d \rfloor} \sigma_0(j) \end{aligned} \]
分块处理后半部分即可.
时间复杂度: 预处理\(O(n)\), 单次询问\(O(n^\frac 1 2)\)
#include <cstdio>
#include <cctype>
#include <cstring>
#include <algorithm>
using namespace std;
namespace Zeonfai
{
inline int getInt()
{
int a = 0, sgn = 1; char c;
while (! isdigit(c = getchar())) if (c == '-') sgn *= -1;
while (isdigit(c)) a = a * 10 + c - '0', c = getchar();
return a * sgn;
}
}
const int N = (int)5e4, MOD = (int)1e9;
typedef int arr[N + 7];
typedef long long Larr[N + 7];
int tot;
arr isNotPrime, prm, mu, minDivisor, minDivisorDegree, sgm;
Larr a, b, c;
inline void initialize()
{
memset(isNotPrime, 0, sizeof(isNotPrime));
tot = 0;
sgm[1] = mu[1] = 1;
for (int i = 2; i <= N; ++ i)
{
if (! isNotPrime[i])
{
prm[tot ++] = i;
mu[i] = -1;
minDivisor[i] = i;
minDivisorDegree[i] = 1;
sgm[i] = 2;
}
for (int j = 0; j < tot && i * prm[j] <= N; ++ j)
{
int x = i * prm[j]; isNotPrime[x] = 1;
if (i % prm[j])
{
mu[x] = - mu[i];
minDivisor[x] = prm[j];
minDivisorDegree[x] = 1;
sgm[x] = sgm[i] * 2;
}
else
{
mu[x] = 0;
minDivisor[x] = minDivisor[i] * prm[j];
minDivisorDegree[x] = minDivisorDegree[i] + 1;
sgm[x] = sgm[i / minDivisor[i]] * (minDivisorDegree[x] + 1);
}
}
}
a[0] = b[0] = c[0] = 0;
for (int i = 1; i <= N; ++ i) a[i] = a[i - 1] + sgm[i], b[i] = a[i] * a[i], c[i] = c[i - 1] + mu[i];
}
int main()
{
using namespace Zeonfai;
initialize();
int T = getInt();
for (int cs = 0; cs < T; ++ cs)
{
int n = getInt(), m = getInt();
long long ans = 0;
int L = 1;
while (L <= min(n, m))
{
int R = min(n / (n / L), m / (m / L));
ans = ans + a[n / L] * a[m / L] * (c[R] - c[L - 1]);
L = R + 1;
}
printf("%lld\n", ans);
}
}