Description
设d(x)为x的约数个数，给定N、M，求 (sum_{i=1}^Nsum_{j=1}^Md(ij))
Input
输入文件包含多组测试数据。
第一行，一个整数T，表示测试数据的组数。
接下来的T行，每行两个整数N、M。
Output
T行，每行一个整数，表示你所求的答案。
Sample Input
2
7 4
5 6
Sample Output
110
121
HINT
1<=N, M<=50000
1<=T<=50000
solution
前置知识：莫比乌斯反演
首先对于(d)函数有一个结论：
[d(ij)=sum_{d|i}sum_{d|j}[gcd(i,j)=1]
]
我就是因为不知道这个结论推了半个小时无果QAQ
证明大致如下：对于(ij)的每一个质因数(x)，在(i)里出现了(a)次，在(j)里出现了(b)次，则总共有(a+b+1)中情况，
然后进行转换，对于任意一个要选(q)个的情况：
- 若(qleqslant a)，就在(i)里选出(q)个，(j)里不选。
- 否则就在(j)里选出(q-a)个，这里的在(j)里选实质上是在(i)里选满了，然后再在(j)里选的，为了不重复计数，在(j)里选(z)个实质上是选了(z+a)个。
对于每一个质因数都这么考虑，然后只要保证(i,j)互质即为一种方案。
然后把这个玩意带到题目给的式子里去，大力推一下，可得：
[ans=sum_{d^prime=1}^{min(n,m)}mu(d^prime)sum_{d_1=1}^{lfloorfrac{n}{d^prime}
floor}lfloorfrac{n}{d_1d^prime}
floorsum_{d_2=1}^{lfloorfrac{m}{d^prime}
floor}lfloorfrac{m}{d_2d^prime}
floor
]
令(f)为式子后面那一块，即：
[f(n)=sum_{i=1}^{n}lfloorfrac{n}{i}
floor
]
然后考虑下这个函数的性质：
对于枚举的(i)，(lfloorfrac{n}{i} floor)实质上就是(n)以内能整除(i)的数的个数。
反过来想，对于每个数，它的每一个约数都对答案有1的贡献，所以(f(n))等价于(n)以内的所有数的约数个数和。
然后线筛下约数个数和莫比乌斯函数，求下前缀和，整除分块下，就做完了。
时间复杂度(O(n+qsqrt{n}))。
#include<bits/stdc++.h>
using namespace std;
#define int long long
void read(int &x) {
x=0;int f=1;char ch=getchar();
for(;!isdigit(ch);ch=getchar()) if(ch=='-') f=-f;
for(;isdigit(ch);ch=getchar()) x=x*10+ch-'0';x*=f;
}
void print(int x) {
if(x<0) putchar('-'),x=-x;
if(!x) return ;print(x/10),putchar(x%10+48);
}
void write(int x) {if(!x) putchar('0');else print(x);putchar('
');}
const int maxn = 5e4+1;
int d[maxn],pri[maxn],tot,p[maxn],vis[maxn],mu[maxn];
void sieve() {
d[1]=mu[1]=1;
for(int i=2;i<maxn;i++) {
if(!vis[i]) pri[++tot]=i,d[i]=2,p[i]=2,mu[i]=-1;
for(int t,j=1;j<=tot&&i*pri[j]<maxn;j++) {
vis[t=i*pri[j]]=1;
if(!(i%pri[j])) {
p[t]=p[i]+1;d[t]=d[i]/p[i]*p[t];
mu[t]=0;break;
}
p[t]=2,d[t]=d[i]*p[t];mu[t]=-mu[i];
}
}
for(int i=1;i<maxn;i++) mu[i]=mu[i]+mu[i-1];
for(int i=1;i<maxn;i++) d[i]=d[i-1]+d[i];
}
signed main() {
sieve();int t;read(t);
while(t--) {
int n,m;read(n),read(m);
int T=1,ans=0;if(n>m) swap(n,m);
while(T<=n) {
int pre=T;T=min(n/(n/T),m/(m/T));
ans+=(mu[T]-mu[pre-1])*d[n/T]*d[m/T];T++;
}write(ans);
}
return 0;
}