[bzoj3994] [SDOI2015]约数个数和

Description

设d(x)为x的约数个数，给定N、M，求 \(\sum_{i=1}^N\sum_{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\)个的情况：

• 若\(q\leqslant 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}^{\lfloor\frac{n}{d^\prime}\rfloor}\lfloor\frac{n}{d_1d^\prime}\rfloor\sum_{d_2=1}^{\lfloor\frac{m}{d^\prime}\rfloor}\lfloor\frac{m}{d_2d^\prime}\rfloor \]
令\(f\)为式子后面那一块，即：
\[ f(n)=\sum_{i=1}^{n}\lfloor\frac{n}{i}\rfloor \]
然后考虑下这个函数的性质：

对于枚举的\(i\)，\(\lfloor\frac{n}{i}\rfloor\)实质上就是\(n\)以内能整除\(i\)的数的个数。

反过来想，对于每个数，它的每一个约数都对答案有1的贡献，所以\(f(n)\)等价于\(n\)以内的所有数的约数个数和。

然后线筛下约数个数和莫比乌斯函数，求下前缀和，整除分块下，就做完了。

时间复杂度\(O(n+q\sqrt{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('\n');}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;
}