设
[f(n)=sum_{i=1}^Nsum_{j=1}^M[(i,j)=n]
]
我们的答案显然是
[ans=sum_{pin prime}f(p)
]
设
[F(n)=sum_{i=1}^Nsum_{j=1}^M[n|(i,j)]
]
即有多少个数对的最大公约数是(n)的倍数
显然(F(n)=left lfloor frac{N}{n} ight floor imesleft lfloor frac{M}{n} ight floor)
同时还存在
[F(n)=sum_{n|d}f(d)
]
看起来并不能反演，但是我们大胆猜测会存在这样的性质
[f(n)=sum_{n|d}mu(frac{d}{n})F(d)
]
看起来很靠谱啊，那就证明一下吧
[sum_{n|d}mu(frac{d}{n})F(d)=sum_{n|d}mu(frac{d}{n})sum_{d|i}f(i)
]
考虑把(f(i))提前，因为(n|d,d|i)，所以(n|i)
[sum_{n|d}mu(frac{d}{n})sum_{d|i}f(i)=sum_{n|i}f(i)sum_{d|i,n|d}mu(frac{d}{n})
]
设(d=kn,i=cd)，则有(i=ckn)
则(frac{d}{n}=k,frac{i}{n}=ck)，所以其实是在(frac{i}{n})的约数
所以可以写成
[sum_{n|i}f(i)sum_{d|i,n|d}mu(frac{d}{n})=sum_{n|i}f(i)sum_{d|frac{i}{n}}mu(d)
]
所以只有在(i=n)的时候(sum_{d|frac{i}{n}}mu(d)=1)，所以这个柿子的值是成立的
所以有一种新的反演形式
[F(n)=sum_{n|d}f(d)
]
就有
[f(n)=sum_{n|d}mu(frac{d}{n})F(d)
]
之后我们的柿子变成了
[ans=sum_{pin prime}sum_{n|p}mu(frac{p}{n})F(n)=sum_{pin prime}sum_{n|p}mu(frac{p}{n}) imes left lfloor frac{N}{n}
ight
floor imesleft lfloor frac{M}{n}
ight
floor
]
于是现在得到了一个复杂度非常玄学的做法，就是枚举(p)之后枚举(p)的倍数
暴力就写好了
#include<iostream>
#include<cstring>
#include<cstdio>
#define re register
#define maxn 10000005
#define LL long long
#define min(a,b) ((a)<(b)?(a):(b))
#define max(a,b) ((a)>(b)?(a):(b))
inline int read()
{
char c=getchar();
int x=0;
while(c<'0'||c>'9') c=getchar();
while(c>='0'&&c<='9')
x=(x<<3)+(x<<1)+c-48,c=getchar();
return x;
}
int mu[maxn],p[maxn],f[maxn];
int n,m,T;
int main()
{
scanf("%d",&T);
f[1]=mu[1]=1;
for(re int i=2;i<=10000000;i++)
{
if(!f[i]) p[++p[0]]=i,mu[i]=-1;
for(re int j=1;j<=p[0]&&p[j]*i<=10000000;j++)
{
f[p[j]*i]=1;
if(i%p[j]==0) break;
mu[p[j]*i]=-1*mu[i];
}
}
while(T--)
{
scanf("%d%d",&n,&m);
LL ans=0;
for(re int i=1;i<=p[0]&&p[i]<=min(n,m);i++)
{
for(re int j=1;j*p[i]<=min(n,m);j++)
ans+=mu[j]*(n/(j*p[i]))*(m/(j*p[i]));
}
printf("%lld
",ans);
}
return 0;
}
把柿子化到这里显然还不够啊
我们需要继续搞一搞
设(k imes p=n)
那么
[sum_{pin prime}sum_{n|p}mu(frac{p}{n}) imes left lfloor frac{N}{n}
ight
floor imesleft lfloor frac{M}{n}
ight
floor=sum_{pin prime}sum_{k=1}^{left lfloorfrac{N}{p}
ight
floor}mu(k) imes left lfloor frac{N}{kp}
ight
floor imesleft lfloor frac{M}{kp}
ight
floor
]
设(T=kp)
于是就有
[sum_{pin prime}sum_{k=1}^{left lfloorfrac{N}{p}
ight
floor}mu(k) imes left lfloor frac{N}{kp}
ight
floor imesleft lfloor frac{M}{kp}
ight
floor=sum_{T=1}^{N}sum_{t|T,tin prime}mu(frac{T}{t}) imes left lfloor frac{N}{T}
ight
floor imesleft lfloor frac{M}{T}
ight
floor
]
[=sum_{T=1}^{N} left lfloor frac{N}{T}
ight
floor imesleft lfloor frac{M}{T}
ight
floorsum_{t|T,tin prime}mu(frac{T}{t})
]
发现好像前面那两个东西可以两个整除分块一起上，后面这个(sum_{t|T,tin prime}mu(frac{T}{t}))看起来好像需要一个前缀和
于是就可以啦
代码
#include<iostream>
#include<cstring>
#include<cstdio>
#define re register
#define maxn 10000005
#define LL long long
#define min(a,b) ((a)<(b)?(a):(b))
#define max(a,b) ((a)>(b)?(a):(b))
inline int read()
{
char c=getchar();
int x=0;
while(c<'0'||c>'9') c=getchar();
while(c>='0'&&c<='9')
x=(x<<3)+(x<<1)+c-48,c=getchar();
return x;
}
int mu[maxn],f[maxn],p[maxn>>2];
LL pre[maxn];
int T,n,m;
int main()
{
f[1]=mu[1]=1;
for(re int i=2;i<=10000000;i++)
{
if(!f[i]) p[++p[0]]=i,mu[i]=-1;
for(re int j=1;j<=p[0]&&p[j]*i<=10000000;j++)
{
f[p[j]*i]=1;
if(i%p[j]==0) break;
mu[p[j]*i]=-1*mu[i];
}
}
for(re int i=1;i<=p[0];i++)
for(re int j=1;j*p[i]<=10000000;j++) pre[j*p[i]]+=mu[j];
for(re int i=1;i<=10000000;i++) pre[i]+=pre[i-1];
scanf("%d",&T);
while(T--)
{
LL ans=0;
n=read(),m=read();
if(n>m) std::swap(n,m);
for(re int l=1,r;l<=n;l=r+1)
{
r=min(n/(n/l),m/(m/l));
ans+=(LL)(n/l)*(m/l)*(pre[r]-pre[l-1]);
}
printf("%lld
",ans);
}
return 0;
}