原题
求x<=a,y<=b,并且gcd(x,y)=d的x,y对数
(sum^{n}_{x=1}sum^{m}_{y=1}[gcd(x,y)==d])
(=sum^{lfloor frac{n}{d}
floor}_{x=1}sum^{lfloor frac{m}{d}
floor}_{y=1}[gcd(x,y)==1])
根据莫比乌斯反演,可以把式子化为
(=sum^{lfloor frac{n}{d}
floor}_{x=1}sum^{lfloor frac{m}{d}
floor}_{y=1}sum_{d|gcd(x,y)}mu(d))
因为$d|gcd(x,y) $ => (d|x,d|y)
所以式子变为 (sum^{n/d}_{i=1}mu(d)*lfloor frac{lfloor frac{n}{d}
floor}{d}
floor*lfloor frac{lfloor frac{m}{d}
floor}{d}
floor)
用n/d的取值分一个段,然后对mu求一个前缀和即可解决……
分段:
因为我们知道(lfloor frac{n}{d}
floor)的取值有(sqrt n)种取值,所以对于(lfloor frac{n}{d}
floor)和(lfloor frac{m}{d}
floor)一共有(sqrt n + sqrt m)种情况
经典枚举方法:
for (int i=1,last=0;i<=a;i=last+1)
{
last=min(a/(a/i),b/(b/i));
ans+=(ll)(a/i)*(b/i)*(sum[last]-sum[i-1]);
}
对mu求前缀和:
利用欧拉筛法:
void init()
{
miu[1]=sum[1]=1;
for (int i=2;i<=N;i++)
{
if (!notprime[i]) miu[i]=-1,prime[++tot]=i;
for (int j=1;j<=tot && i*prime[j]<=N;j++)
{
notprime[i*prime[j]]=1;
if (i%prime[j]) miu[i*prime[j]]=-miu[i];
else
{
miu[i*prime[j]]=0;
break;
}
}
sum[i]=sum[i-1]+miu[i];
}
}
AC代码:
#include<cstdio>
#include<algorithm>
#define N 50000
typedef long long ll;
using namespace std;
int t,a,b,d,miu[N+10],prime[N+10],tot,sum[N];
bool notprime[N+10];
ll ans;
int read()
{
int ans=0,fu=1;
char j=getchar();
for (;j<'0' || j>'9';j=getchar()) if (j=='-') fu=-1;
for (;j>='0' && j<='9';j=getchar()) ans*=10,ans+=j-'0';
return ans*fu;
}
void init()
{
miu[1]=sum[1]=1;
for (int i=2;i<=N;i++)
{
if (!notprime[i]) miu[i]=-1,prime[++tot]=i;
for (int j=1;j<=tot && i*prime[j]<=N;j++)
{
notprime[i*prime[j]]=1;
if (i%prime[j]) miu[i*prime[j]]=-miu[i];
else
{
miu[i*prime[j]]=0;
break;
}
}
sum[i]=sum[i-1]+miu[i];
}
}
int main()
{
init();
t=read();
while (t--)
{
a=read();b=read();d=read();
a/=d;
b/=d;
ans=0;
if (a>b) swap(a,b);
for (int i=1,last=0;i<=a;i=last+1)
{
last=min(a/(a/i),b/(b/i));
ans+=(ll)(a/i)*(b/i)*(sum[last]-sum[i-1]);
}
printf("%lld
",ans);
}
return 0;
}