題解 P3327 【[SDOI2015]約數個數和】
首先有個東西是這題解題的關鍵：
\[{\rm{d}}(ij)=\sum_{x|i}\sum_{y|i}[{\rm{gcd}}(x,y)=1] \]
有位dalao的題解證明了這個式子，可以去看看。
然後就可以開始推式子了：
\[\sum_{i=1}^n\sum_{j=1}^md(ij) =\sum_{i=1}^n\sum_{j=1}^m\sum_{x|i}\sum_{y|i}[{\rm{gcd}}(x,y)=1]=\sum_{x=1}^n\sum_{y=1}^m[{\rm{gcd}}(x,y)=1]\lfloor\frac{n}{x}\rfloor\lfloor\frac{m}{y}\rfloor \]
接下來用莫比烏斯函式替換掉 \([{\rm{gcd}}(x,y)=1]\) ：
\[\sum_{d|x,d|y}\mu(d)=[{\rm{gcd}}(x,y)=1] \]
代入：
\[\begin{align} \sum_{x=1}^n\sum_{y=1}^m[{\rm{gcd}}(x,y)=1]\lfloor\frac{n}{x}\rfloor\lfloor\frac{m}{y}\rfloor & =\sum_{x=1}^n\sum_{y=1}^m \big(\sum_{d|x,d|y}\mu(d)\big)\lfloor\frac{n}{x}\rfloor\lfloor\frac{m}{y}\rfloor \\ \end{align} \]
將這個式子由列舉 \(x,y\) 的形式轉變為列舉 \(d\) ，即用 \(xd,yd\) 代替 \(x,y\) ：
\[= \sum_{d=1}^{{\rm{min}}(n,m)}\mu(d)\sum_{x=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{y=1}^{\lfloor\frac{m}{d}\rfloor}\lfloor\frac{n}{xd}\rfloor\lfloor\frac{m}{yd}\rfloor = \sum_{d=1}^{{\rm{min}}(n,m)}\mu(d)\sum_{x=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{y=1}^{\lfloor\frac{m}{d}\rfloor}\lfloor\frac{\lfloor\frac{n}{d}\rfloor}{x}\rfloor\lfloor\frac{\lfloor\frac{m}{d}\rfloor}{y}\rfloor \]
這裡列舉的 \(x,y\) 是上式中的 \(\frac{x}{d},\frac{y}{d}\) 。
設 \(f(n)=\sum_{i=1}^{n}\lfloor\frac{n}{i}\rfloor\) ，將其代入上式：
\[\sum_{i=1}^n\sum_{j=1}^md(ij)=\sum_{d=1}^{{\rm{min}}(n,m)}\mu(d)f(\lfloor\frac{n}{d}\rfloor)f(\lfloor\frac{m}{d}\rfloor) \]
到這裡，很多題解就直接採用整除分塊來計算 \(f\) 了，複雜度 \(O(n\sqrt{n})\) 。但是其實還有一種更快的做法能做到 \(O(n)\) 預處理。
不難發現 \(f(n)\) 其實是計算了 \(1 \sim n\) 每個數小於等於 \(n\) 的倍數個數之和，也就是說一個數對 \(f(n)\) 做出的貢獻就是它的約數個數。那麼有：
\[f(n)=\sum_{i=1}^{n}\lfloor\frac{n}{i}\rfloor=\sum_{i=1}^n{\rm{d}}(i) \]
由於 \({\rm{d}}\) 是積性函式，可以使用線性篩做到 \(O(n)\) 預處理。
關於線性篩 \({\rm{d}}\) ，可以看這個。
\({\rm{Code}}:\)
#include <iostream>
#include <cstring>
#include <cstdio>
#include <algorithm>
#include <cmath>
#define maxn 50005
#define Rint register int
#define INF 0x3f3f3f3f
using namespace std;
typedef long long lxl;
template <typename T>
inline T read()
{
T x=0,f=1;char ch=getchar();
while(ch<'0'||ch>'9') {if(ch=='-') f=-1;ch=getchar();}
while(ch>='0'&&ch<='9') {x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
return x*f;
}
int prime[maxn],cnt;
bool flag[maxn];
lxl d[maxn],f[maxn],mu[maxn],sum[maxn],num[maxn];
inline void sieve()
{
d[1]=1;
mu[1]=1;
for(int i=2;i<maxn;++i)
{
if(!flag[i]) prime[++cnt]=i,d[i]=2,num[i]=1,mu[i]=-1;
for(int j=1;j<=cnt&&i*prime[j]<maxn;++j)
{
flag[i*prime[j]]=true;
if(i%prime[j])
{
num[i*prime[j]]=1;
d[i*prime[j]]=d[i]*2;
}
else
{
num[i*prime[j]]=num[i]+1;
d[i*prime[j]]=d[i]/num[i*prime[j]]*(num[i*prime[j]]+1);
break;
}
mu[i*prime[j]]=-mu[i];
}
}
for(int i=1;i<maxn;++i)
f[i]=f[i-1]+d[i],sum[i]=sum[i-1]+mu[i];
}
inline lxl calcu(int n,int m)
{
lxl res=0;
if(n>m) swap(n,m);
for(int l=1,r=0;l<=n;l=r+1)
{
r=min(n/(n/l),m/(m/l));
res+=f[n/l]*f[m/l]*(sum[r]-sum[l-1]);
}
return res;
}
int main()
{
// freopen("P3327.in","r",stdin);
sieve();
int T=read<int >();
while(T--)
{
int n=read<int >(),m=read<int >();
printf("%lld\n",calcu(n,m));
}
return 0;
}