题目
对于正整数n，定义f(n)为n所含质因子的最大幂指数。例如f(12250)=f(2^1 * 5^3 * 7^2)=3, f(10007)=1, f(1)=0。
给定正整数a,b，求(sumlimits_{i=1}^{a}{sumlimits_{j=1}^{b}{f(gcd(i,j))}})。
题解
[sumlimits_{i=1}^{a}{sumlimits_{j=1}^{b}{f(gcd(i,j))}}
]
[sumlimits_dsumlimits_{i=1}^{frac{a}{d}}{sumlimits_{j=1}^{frac{b}{d}}{f(d)[gcd(i,j)=1]}}
]
[sumlimits_dsumlimits_{i=1}^{frac{a}{d}}{sumlimits_{j=1}^{frac{b}{d}}{f(d)sumlimits_{d_1|gcd(i,j)}{mu(d_1)}}}
]
[sumlimits_d{sumlimits_{d_1}{mu(d_1)f(d)lfloorfrac{a}{dd_1}
floorlfloorfrac{b}{dd_1}
floor}}
]
令(T=dd_1)
[sumlimits_T^{min(a, b)}{lfloorfrac{a}{T}
floorlfloorfrac{b}{T}
floorsumlimits_{d|T}{f(d)mu(frac{T}{d})}}
]
现在关键是求(sumlimits_{d|T}{f(d)mu(frac{T}{d})})里的(f)。观察可发现，(f)有类似积性函数的性质，即当(gcd(i,j)=1)时，有(f(ij)=max(f(i),f(j)))，可以使用线性筛。多维护一个最小质因子的次数的函数num，在i%p==0情况时容易想到有(f(i*p)=max(f(i), num(i *p)))
然后就是常规操作数论分块了。
#include <bits/stdc++.h>
#define endl '
'
#define IOS std::ios::sync_with_stdio(0); cin.tie(0); cout.tie(0)
#define mp make_pair
#define seteps(N) fixed << setprecision(N)
typedef long long ll;
using namespace std;
/*-----------------------------------------------------------------*/
ll gcd(ll a, ll b) {return b ? gcd(b, a % b) : a;}
#define INF 0x3f3f3f3f
const int N = 1e7 + 10;
const double eps = 1e-5;
int mu[N];
ll f[N];
int prime[N];
int cnt;
int vis[N];
int num[N];
int mx[N];
void init() {
mu[1] = num[1] = 1;
for(int i = 2; i < N; i++) {
if(!vis[i]) {
mu[i] = -1;
num[i] = 1;
mx[i] = 1;
prime[cnt++] = i;
}
for(int j = 0; j < cnt; j++) {
int p = prime[j];
if(i * p > N) break;
vis[i * p] = 1;
if(i % p == 0) {
mu[i * p] = 0;
num[i * p] = num[i] + 1;
mx[i * p] = max(num[i * p], mx[i]);
break;
}
mu[i * p] = -mu[i];
num[i * p] = 1;
mx[i * p] = mx[i];
}
}
num[1] = 0;
for(int i = 1; i < N; i++) {
for(int j = i; j < N; j += i) {
f[j] += mx[i] * mu[j / i];
}
}
for(int i = 2; i < N; i++) f[i] += f[i - 1];
}
int main() {
IOS;
init();
int t;
cin >> t;
while(t--) {
ll ans = 0;
int n, m;
cin >> n >> m;
int cur = 1;
while(cur <= min(n, m)) {
int p1 = n / (n / cur);
int p2 = m / (m / cur);
int nt = min(min(n, m), min(p1, p2));
ans += 1ll * (n / cur) * (m / cur) * (f[nt] - f[cur - 1]);
cur = nt + 1;
}
cout << ans << endl;
}
}