YY的GCD
原题连接c++
这应该是我作的第一道莫比乌斯反演的题目。
题目描述
- 神犇YY虐完数论后给傻×kAc出了一题
- 给定N, M,求1<=x<=N, 1<=y<=M且gcd(x, y)为质数的(x, y)有多少对
kAc这种傻×必然不会了,因而向你来请教……
- 多组输入
输入输出格式
输入格式:
- 第一行一个整数T 表述数据组数
- 接下来T行,每行两个正整数,表示N, M
输出格式:
说明
解题思路
- 显然,题目要求的\(Ans\)实际上就是\(\sum_{i=1}^{n}\sum_{j=1}^{m}[gcd(x,y)=prim]\)的值
- 接下来,咱们就开始进行欢乐的推式子了
- 对于这种与\(gcd\)有关的莫比乌斯反演,通常咱们都是套路的去设\(f(d)\)为\(gcd(i,j)=d\)的个数,\(F(n)\)为\(gcd(i,j)=n\)和\(n\)的倍数的个数,即:
\[f(d)=\sum_{i=1}^{n}\sum_{j=1}^{m}[gcd(i,j)=d]\]
\[F(n)=\sum_{n|d}f(d)=\lfloor\frac{N}{n}\rfloor\lfloor\frac{M}{n}\rfloor\]
\[f(n)=\sum_{n|d}\mu(\lfloor\frac{d}{n}\rfloor)F(d)\]
- 这样,咱们即可以开心的化简这个式子了!
\[Ans=\sum_{p\in prim}\sum_{i=1}^{n}\sum_{j=1}^{m}[gcd(i,j)=p]\]
将\(f(p)\)带入得:
\[Ans=\sum_{p\in prim}f(p)\]
而后就莫比乌斯反演一下
\[Ans=\sum_{p\in prim}\sum_{p|d}\mu(\lfloor\frac{d}{p}\rfloor)F(d)\]
咱们换一个枚举项,咱们枚举\(\lfloor\frac{d}{p}\rfloor\)
\[Ans=\sum_{p\in prim}\sum_{d=1}^{min(\lfloor\frac{n}{p}\rfloor,\lfloor\frac{m}{p}\rfloor)}\mu(d)F(dp)=\sum_{p\in prim}\sum_{d=1}^{min(\lfloor\frac{n}{p}\rfloor,\lfloor\frac{m}{p}\rfloor)}\mu(d)\lfloor\frac{n}{dp}\rfloor\lfloor\frac{m}{dp}\rfloor\]
这个\(dp\)一看就很不爽,因而咱们把它换成\(T\)
\[Ans=\sum_{T=1}^{min(n,m)}\sum_{t|T,t\in prime}\mu(\lfloor\frac{T}{t}\rfloor)\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor\]
\[Ans=\sum_{T=1}^{min(n,m)}\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor(\sum_{t|T,t\in prime}\mu(\lfloor\frac{T}{t}\rfloor))\]
- 推到这里,咱们就能够开始作了。若是是单组询问,咱们就直接\(O(n)\)作。
(不过好像通常这种题,都不会让你直接处理。)若是是多组数据的话,咱们就只要在打一个简单的整除分块就能够了。后面的\(\mu\)函数能够线筛出来。因为整除分块的缘故,咱们就只须要记一个前缀和就能够了。
- 下面贴一个完整的代码吧。
#include<bits/stdc++.h>
#define N 10000100
using namespace std;
inline void read(int &x)
{
x=0;
static int p;p=1;
static char c;c=getchar();
while(!isdigit(c)){if(c=='-')p=-1;c=getchar();}
while(isdigit(c)) {x=(x<<1)+(x<<3)+(c-48);c=getchar();}
x*=p;
}
inline void print(long long x)
{
static int cnt;
static int a[15];
cnt=0;
do
{
a[++cnt]=x%10;
x/=10;
}while(x);
for(int i=cnt;i>=1;i--)putchar(a[i]+'0');
puts("");
}
bool vis[N];
long long sum[N];
int prim[N];
int mu[N],g[N];
int cnt;
void get_mu(int n)
{
mu[1]=1;
for(int i=2;i<=n;i++)
{
if(!vis[i]){mu[i]=-1;prim[++cnt]=i;}
for(int j=1;j<=cnt&&prim[j]*i<=n;j++)
{
vis[i*prim[j]]=1;
if(i%prim[j]==0)break;
else mu[prim[j]*i]=-mu[i];
}
}
for(int j=1;j<=cnt;j++)
for(int i=1;i*prim[j]<=n;i++)g[i*prim[j]]+=mu[i];
for(int i=1;i<=n;i++)sum[i]=sum[i-1]+(long long)g[i];
}
int n,m;
int main()
{
#ifndef ONLINE_JUDGE
freopen("P2257.in","r",stdin);
freopen("P2257.out","w",stdout);
#endif
int t;
read(t);
get_mu(10000000);
while(t--)
{
read(n);read(m);
if(n>m)swap(n,m);
static long long ans;ans=0;
for(int l=1,r;l<=n;l=r+1)
{
r=min(n/(n/l),m/(m/l));
ans+=1ll*(n/l)*(m/l)*(sum[r]-sum[l-1]);
}
print(ans);
}
return 0;
}