题面
求(sum_{i=1}^nsum_{j=1}^md(ij))
(leq 50000)组数据，(1leq n,mleq 50000)。
题目分析
首先，你需要知道一个结论：
你可以认为(x,y)表示你选择的因数为(frac i x cdot y)，即：(x)表示(i)中不要的部分，(y)表示(j)中要的部分。
如果(gcd(x,y)==p_i)，那么(frac i x)表示在约数中拿掉(p_i)，(y)表示在约数中加入(p_i)，这样一拿一加，自然会在答案中重复。
那么，现在我们的问题转化为求
这样还是无法计算，所以我们把枚举因数提前
现在看起来就可以反演了，设(f(x))表示(gcd(i,j)==x)时的答案，(g(x))表示(gcd(i,j)==kx,kin Z)时的答案，则：
比较巧的一点是：(sumlimits_{i=1}^nlfloorfrac n i floor)可以表示(1 sim n)的约数个数的前缀和。
约数个数可以在线性筛中预处理，原理如下：
对于(x=p_1^{a_1}p_2^{a_2}p_3^{a_3}...p_n^{a_n})，(x)的约数个数为((a_1+1)cdot(a_2+1)cdot(a_3+1)cdot...cdot(a_n+1))
由于在线性筛中，每个数只会被它最小的质因子更新，所以：
如果(i\%prime[j]==0)，说明(i)中含有(prime[j])，此时(x)中(prime[j])的个数为(i)中(prime[j])的个数(+1)，(x)的约数个数=(i)的约数个数/((i)中(prime[j])的个数)*((i)中(prime[j])的个数(+1))；
如果(i\%prime[j]!=0)，说明(prime[j])在(x)中只有(1)个，(x)的约数个数=(i)的约数个数*(2)。
这样一来(g(x))可以进行预处理，然后(O(1))计算。
反演得(f(x)=sumlimits_{x|d}^nmu(frac dx)g(d))，为了针对多组数据，整除分块即可。
代码实现
#include<iostream>
#include<cstring>
#include<cmath>
#include<algorithm>
#include<cstdio>
#include<iomanip>
#include<cstdlib>
#define MAXN 0x7fffffff
typedef long long LL;
const int N=50005;
using namespace std;
inline int Getint(){register int x=0,f=1;register char ch=getchar();while(!isdigit(ch)){if(ch=='-')f=-1;ch=getchar();}while(isdigit(ch)){x=x*10+ch-'0';ch=getchar();}return x*f;}
int mu[N],prime[N];
bool vis[N];
int ys[N],lw[N],g[N];
int main(){
mu[1]=g[1]=1;
for(int i=2;i<=50000;i++){
if(!vis[i])prime[++prime[0]]=i,mu[i]=-1,ys[i]=2,lw[i]=1;
for(int j=1;j<=prime[0]&&i*prime[j]<=50000;j++){
vis[i*prime[j]]=1;
if(i%prime[j]==0){
ys[i*prime[j]]=ys[i]/(lw[i]+1)*(lw[i]+2);
lw[i*prime[j]]=lw[i]+1;
break;
}
mu[i*prime[j]]=-mu[i];
ys[i*prime[j]]=ys[i]*2,lw[i*prime[j]]=1;
}
mu[i]+=mu[i-1],g[i]=g[i-1]+ys[i];
}
int T=Getint();
while(T--){
int n=Getint(),m=Getint();
if(n>m)swap(n,m);
LL ans=0;
for(int l=1,r;l<=n;l=r+1){
r=min(n/(n/l),m/(m/l));
ans+=1ll*(mu[r]-mu[l-1])*g[n/l]*g[m/l];
}
cout<<ans<<'
';
}
return 0;
}