【数论】YY的GCD（P2257）

正题

P2257

题目大意

给你T组询问，每组询问给出n,m，让你求 1≤x≤n,1≤y≤m1\leq x\leq n,1\leq y\leq m1≤x≤n,1≤y≤m 且 gcd(x,y)=primegcd(x,y)=primegcd(x,y)=prime 的方案数

解题思路

根据题意，有

∑i=1n∑j=1m[gcd(i,j)=prime]\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=prime] i=1∑nj=1∑m[gcd(i,j)=prime]

考虑把prime拿出来枚举

∑p∈primen∑i=1n/p∑j=1m/p[gcd(i,j)=1]∑p∈primen∑i=1n/p∑j=1m/p∑d∣i,d∣jμ(d)∑p∈primemin(n,m)∑d=1min(n,m)/pμ(d)×⌊ndp⌋×⌊mdp⌋\sum_{p\in prime}^n\sum_{i=1}^{n/p}\sum_{j=1}^{m/p}[gcd(i,j)=1]\\ \sum_{p\in prime}^n\sum_{i=1}^{n/p}\sum_{j=1}^{m/p}\sum_{d|i,d|j}\mu(d)\\ \sum_{p\in prime}^{min(n,m)}\sum_{d=1}^{min(n,m)/p}\mu(d)\times\left\lfloor\frac{n}{dp}\right\rfloor\times\left\lfloor\frac{m}{dp}\right\rfloor\\ p∈prime∑ni=1∑n/pj=1∑m/p[gcd(i,j)=1]p∈prime∑ni=1∑n/pj=1∑m/pd∣i,d∣j∑μ(d)p∈prime∑min(n,m)d=1∑min(n,m)/pμ(d)×⌊dpn⌋×⌊dpm⌋

然后可以考虑改变枚举变量

可以发现后面的两个整除都与dp有关，那么考虑枚举dp，那么有

∑k=1min(n,m)∑p∈prime,p∣kμ(kp)×⌊nk⌋×⌊mk⌋∑k=1min(n,m)⌊nk⌋⌊mk⌋∑p∈prime,p∣kμ(kp)\sum_{k=1}^{min(n,m)}\sum_{p\in prime,p|k}\mu(\frac{k}{p})\times\left\lfloor\frac{n}{k}\right\rfloor\times\left\lfloor\frac{m}{k}\right\rfloor\\ \sum_{k=1}^{min(n,m)}\left\lfloor\frac{n}{k}\right\rfloor\left\lfloor\frac{m}{k}\right\rfloor\sum_{p\in prime,p|k}\mu(\frac{k}{p})\\ k=1∑min(n,m)p∈prime,p∣k∑μ(pk)×⌊kn⌋×⌊km⌋k=1∑min(n,m)⌊kn⌋⌊km⌋p∈prime,p∣k∑μ(pk)

对于后面的∑p∈prime,p∣kμ(kp)\sum_{p\in prime,p|k}\mu(\frac{k}{p})∑p∈prime,p∣kμ(pk)，可以预处理，先枚举质数，然后枚举倍数，然后计算贡献

对于前面的一段再整除分块就好了

时间复杂度 O(Tn)O(T\sqrt{n})O(Tn)

code

#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#define ll long long
#define N 10000010
using namespace std;
ll T,n,m,ans,w,p[N],g[N],mu[N],prime[N];
const ll MX=1e7;
void work()
{p[1]=mu[1]=1;for(ll i=2;i<=MX;++i){if(!p[i]){mu[i]=-1;prime[++w]=i;}for(ll j=1;j<=w&&i*prime[j]<=MX;++j){p[i*prime[j]]=1;if(i%prime[j]==0)break;else mu[i*prime[j]]=-mu[i];}}for(ll i=1;i<=w;++i)for(ll j=1;j*prime[i]<=MX;++j)g[j*prime[i]]+=mu[j];//对每个数提出一个质因数计算mu的贡献for(ll i=2;i<=MX;++i)g[i]+=g[i-1];return;
}
int main()
{scanf("%lld",&T);work();while(T--){scanf("%lld%lld",&n,&m);ans=0;for(ll l=1,r=0;l<=min(n,m);l=r+1){r=min(n/(n/l),m/(m/l));ans+=(n/l)*(m/l)*(g[r]-g[l-1]);}printf("%lld\n",ans);}return 0;
}