题目大意
对于正整数 \(n\)，定义 \(f(n)\) 为 \(n\) 所包含质因子的最大幂指数。例如 \(f(1960)=f(2^3 \times 5^1 \times 7^2)=3\)，\(f(10007)=1\)，\(f(1)=0\)。
给定正整数 \(a,b\)，求下式的值：
\[\sum^{a}_{i=1}\sum^{b}_{j=1}f(\gcd(i,j)) \]推导
首先记 \(n=\min(a,b),m=\max(a,b)\)。
\[\begin{aligned} & \ \ \ \ \ \ \sum^{a}_{i=1}\sum^{b}_{j=1}f(\gcd(i,j)) \\ &= \sum^{n}_{d=1}f(d)\sum^{n}_{i=1}\sum^{m}_{j=1}\left[ \gcd(i,j)=d \right] \\ &= \sum^{n}_{d=1}f(d)\sum^{\left\lfloor \frac{n}{d} \right\rfloor}_{i=1}\sum^{\left\lfloor \frac{m}{d} \right\rfloor}_{j=1}\left[ \gcd(i,j)=1 \right] \\ &= \sum^{n}_{d=1}f(d)\sum^{\left\lfloor \frac{n}{d} \right\rfloor}_{i=1}\sum^{\left\lfloor \frac{m}{d} \right\rfloor}_{j=1}\sum_{t \mid \gcd(i,j)}\mu(t) \\ &= \sum^{n}_{d=1}f(d)\sum^{\left\lfloor \frac{n}{d} \right\rfloor}_{i=1}\sum^{\left\lfloor \frac{m}{d} \right\rfloor}_{j=1}\sum_{t \mid i \land t \mid j}\mu(t) \\ &= \sum^{n}_{d=1}f(d) \sum^{\left\lfloor \frac{n}{d} \right\rfloor}_{t =1}\mu(t) \sum^{\left\lfloor \frac{n}{d} \right\rfloor}_{i=1} \left[ t \mid i \right] \sum^{\left\lfloor \frac{m}{d} \right\rfloor}_{j=1} \left[ t \mid j \right] \\ &= \sum^{n}_{d=1}f(d)\sum^{\left\lfloor \frac{n}{d} \right\rfloor}_{t=1}\mu(t) \left\lfloor \frac{n}{dt} \right\rfloor \left\lfloor \frac{m}{dt} \right\rfloor \\ \end{aligned} \]设 \(T=dt\)（十分常用的技巧），那么有
\[= \sum^{n}_{T=1} \left\lfloor \frac{n}{T} \right\rfloor \left\lfloor \frac{m}{T} \right\rfloor \sum_{d \mid T}\mu(\frac{T}{d})f(d) \]记 \(h=f*\mu\)，那么有
\[= \sum^{n}_{T=1} \left\lfloor \frac{n}{T} \right\rfloor \left\lfloor \frac{m}{T} \right\rfloor h(T) \]那么，现在的问题是如何求 \(h\)？考虑求 \(h(n)\)。
我们可以把 \(n\) 进行质因数分解：
\[n=\prod^{m}_{i=1}p_i^{c_i} \]考虑 \(h(T)\) 的原本写法：
\[h(T)=\sum_{d \mid T}\mu(\frac{T}{d})f(d) \]记 \(T=\prod^{m}_{i=1}p_i^{a_i}\)，\(d=\prod^{m}_{i=1}p_i^{b_i}\)。
考虑 \(\mu(n)\) 的定义，当 \(n\) 中含有平方质因子时 \(\mu(n)=0\)，即不会在 \(h\) 中产生贡献，可以不考虑。
所以 \(\dfrac{T}{d}\) 的质因数要满足最大幂指数小于 \(2\)，即 \(a_i-b_i \in \left\{ 0,1 \right\}\)。
所以将 \(h(T)\) 改写为如下形式时
\[\displaystyle h(T) = \sum_{ab = T} f(a) \mu(b) \]有
\[b=\prod^{m}_{i=1}p_i^{d_i}(d_i \in \left\{ 0,1 \right\}) \]记 \(l=\max^m\limits_{i=1}c_i,k=\sum^{m}\limits_{i=1} \left[c_i=l \right]\)。即 \(l\) 为 \(n\) 所包含质因子的最大幂指数，\(k\) 为 \(n\) 所包含质因子幂指数中最大幂指数的个数。
我们发现 \(f(a)\) 的取值只有 \(l\) 和 \(l-1\) 两种可能（\(b\) 可能把最大幂指数的都取走，导致 \(f(a)\) 的取值少了 \(1\)）。
我们先按 \(k=m\) 和 \(k \neq m\) 两种情况分类讨论，在每一项讨论中，我们再分 \(f(a)=l\) 和 \(f(a)=l-1\) 两种子情况。
当 \(k=m\) 时，贡献为（加号前为 \(l\) 的情况，加号后为 \(l-1\) 的情况）
\[\sum^{m-1}_{s=0}(l) \times (-1)^s \times {m \choose s} + (l-1) \times (-1)^m \times {m \choose m} \]\(f(a)=l\) 的情况不会把最大幂指数的质因数都取走。
所以我们枚举到 \(m-1\)，中间的 \((-1)^s\) 和 \((-1)^m\) 实质上就是莫比乌斯函数，最右边的二项式系数是枚举选取的方案。
将上面的式子加号右边的部分拆开，把带 \(l\) 的部分合并到左面，得到
\[\sum^{m}_{s=0}(l) \times (-1)^s \times {m \choose s} + 1 \times (-1)^m \]二项式定理，得
\[-1 \times (-1)^m=(-1)^{m+1} \]当 \(k \neq m\) 时，
当 \(f(a)=l\) 时，设在 \(k\) 个 \(c_i=l\) 的质数中选了 \(t\) 个，另外 \(m-k\) 个质数中选了 \(s\) 个，贡献为：
\[\sum^{m-k}_{s=0}\sum^{k-1}_{t=0} {k \choose t} \times l \times (-1)^{s+t} \times {m-k \choose s} \]中间的 \((-1)^{s+t}\) 仍是莫比乌斯函数。因为 \(f(a)=l\)，所以 \(k\) 个质数不能被选完，因此枚举到 \(k-1\)。
左右拆开，分别二项式定理，得
\[\begin{aligned} \sum^{k-1}_{t=0} {k \choose t} \times l \times (-1)^t \times \sum^{m-k}_{s=0}(-1)^s \times 1 ^{m-k-s} \times {m - k \choose s} =0 \end{aligned} \]当 \(f(a)=l-1\) 时，\(k\) 个 \(c_i=l\) 的质数一定全部被选择，设在另外 \(m-k\) 个质数中选择 \(s\) 个，贡献为：
\[\begin{aligned} & \ \ \ \ \ \sum^{m-k}_{s=0}(l-1) \times (-1)^k \times (-1)^s \times {m-k \choose s} \\ &= \sum^{m-k}_{s=0}(l-1) \times (-1)^k \times (-1)^s \times 1^{m-k-s} \times {m-k \choose s} \\ &= (l-1) \times (-1)^k \times \sum^{m-k}_{s=0} (-1)^s \times 1^{m-k-s} \times {m-k \choose s} \\ &= 0 \end{aligned} \]所以最终 \(h(n)\) 的表达式为
\[h(n)= \begin{cases} (-1)^{m + 1} & k = m \\ 0 & k \neq m \end{cases} \]Code
#include <bits/stdc++.h>
using namespace std;
#define int long long
const int N = 1e7;
const int M = 10500;
int T,n,m;
bool p[N + M];
int pri[N >> 2],cnt,low[N + M];
int h[N + M],a[N + M];
void Seive() {
p[1] = 1;
for(int i = 2;i <= N; i++) {
if(!p[i]) {
cnt ++;
low[i] = pri[cnt] = i;
a[i] = h[i] = 1;
}
for(int j = 1;j <= cnt && i * pri[j] <= N; j++) {
p[i * pri[j]] = 1;
if(i % pri[j] == 0) {
a[i * pri[j]] = a[i] + 1;
low[i * pri[j]] = low[i] * pri[j];
if(i == low[i])
h[i * pri[j]] = 1;
else if(a[i / low[i]] == a[i * pri[j]])
h[i * pri[j]] = -h[i / low[i]];
else
h[i * pri[j]] = 0;
break;
}
a[i * pri[j]] = 1;
low[i * pri[j]] = pri[j];
if(a[i] == 1)
h[i * pri[j]] = -h[i];
else
h[i * pri[j]] = 0;
}
}
for(int i = 1;i <= N; i++)
h[i] += h[i - 1];
return ;
}
signed main() {
ios::sync_with_stdio(false);
Seive();
cin >> T;
while(T --> 0) {
cin >> n >> m;
if(n > m)
swap(n,m);
int l = 1,r = 0,ans = 0;
while(l <= n) {
r = min(n / (n / l),m / (m / l));
ans += (n / l) * (m / l) * (h[r] - h[l - 1]);
l = r + 1;
}
cout << ans << "\n";
}
return 0;
}