BZOJ
Luogu
题意:
给定n,m,求(sum_{i=1}^{n}sum_{j=1}^{m}d(ij)),其中(d(x))表示x的约数个数。多组数据,n,m<=50000,T<=50000
sol
首先我们大胆猜想,
[d(ij)=sum_{u|i}sum_{v|j}[gcd(u,v)==1]
]
(证明晚上再补。。。)
好然后我们看,我们已知
[ans=sum_{i=1}^{n}sum_{j=1}^{m}sum_{u|i}sum_{v|j}[gcd(u,v)==1]
]
显然这四个(sum)很不好搞,所以我们考虑枚举(u,v),计算每一对((u,v))的贡献。
[ans=sum_{u=1}^{n}sum_{v=1}^{m}[gcd(u,v)==1]lfloor frac nu
floorlfloor frac mv
floor
]
写出这个形式那就好办了,
[f(d)=sum_{u=1}^{n}sum_{v=1}^{m}[gcd(u,v)==d]lfloor frac nu
floorlfloor frac mv
floor
]
[F(d)=sum_{u=1}^{n}sum_{v=1}^{m}[d|gcd(u,v)]lfloor frac nu
floorlfloor frac mv
floor
]
在(F(d))的表达式中显然(u)和(v)都是(d)的倍数,所以我们可以令(u=id,v=jd)然后
[F(d)=sum_{i=1}^{n/d}sum_{j=1}^{m/d}lfloor frac {n}{id}
floorlfloor frac {m}{jd}
floor=sum_{i=1}^{n/d}lfloor frac {n/d}{i}
floor * sum_{j=1}^{m/d}lfloor frac {m/d}{j}
floor
]
注意上面的那个结构,形如(sum_{i=1}^{n}frac ni),我们把它记作(sum(n))。如果你做过这道题[AHOI2005]约数研究就应该不难知道这是啥。
(sum(n)=sum_{i=1}^{n} frac ni)表示1~n中每个数的约数个数和
所以就是对每个数求一下约数个数再取前缀和即可。
对一个数求约数个数使用唯一分解定理,复杂度(O(nsqrt n))(预处理)
别忘了答案式
[ans=f(1)=sum_{d=1}^{n}mu(d)F(d)=sum_{d=1}^{n}mu(d)sum(lfloor frac {n}{d}
floor)sum(lfloor frac {m}{d}
floor)
]
(O(n))处理出(mu(d))的前缀和然后直接数论分块一波
复杂度(O(Tsqrt n))
code
#include<cstdio>
#include<algorithm>
using namespace std;
#define ll long long
const int N = 50000;
int gi()
{
int x=0,w=1;char ch=getchar();
while ((ch<'0'||ch>'9')&&ch!='-') ch=getchar();
if (ch=='-') w=0,ch=getchar();
while (ch>='0'&&ch<='9') x=(x<<3)+(x<<1)+ch-'0',ch=getchar();
return w?x:-x;
}
int mu[N+5],pri[N+5],tot,zhi[N+5];
ll s[N+5],f[N+5];
void Mobius()
{
zhi[1]=mu[1]=1;
for (int i=2;i<=N;i++)
{
if (!zhi[i]) pri[++tot]=i,mu[i]=-1;
for (int j=1;j<=tot&&i*pri[j]<=N;j++)
{
zhi[i*pri[j]]=1;
if (i%pri[j]) mu[i*pri[j]]=-mu[i];
else {mu[i*pri[j]]=0;break;}
}
}
for (int i=1;i<=N;i++)
s[i]=s[i-1]+mu[i];
}
int Divide(int x)
{
int p[10]={0},k[10]={0},t=0;
for (int i=2;i*i<=x;i++)
if (x%i==0)
{
p[++t]=i;
while (x%i==0) k[t]++,x/=i;
}
if (x>1) p[++t]=x,k[t]=1;
int res=1;
for (int i=1;i<=t;i++)
res*=k[i]+1;
return res;
}
int main()
{
Mobius();
for (int i=1;i<=N;i++)
f[i]=f[i-1]+Divide(i);
int T=gi();
while (T--)
{
int n=gi(),m=gi();
if (n>m) swap(n,m);
int i=1;ll ans=0;
while (i<=n)
{
int j=min(n/(n/i),m/(m/i));
ans+=(s[j]-s[i-1])*f[n/i]*f[m/i];
i=j+1;
}
printf("%lld
",ans);
}
return 0;
}