[题解] LuoguP2257 YY的GCD

传送门

给(n,m)，让你求

[sumlimits_{i=1}^n sumlimits_{j=1}^m [gcd(i,j) in prime] ]

有(T)组询问((T le 10^4,n,mle 10^7))。

枚举质数(p)，然后柿子变成

[sumlimits_{p} sumlimits_{i=1}^nsumlimits_{j=1}^m [gcd(i,j)=p] ]

等价于

[sumlimits_{p} sumlimits_{i=1}^{leftlfloor n/p ight floor}sumlimits_{j=1}^{leftlfloor m/p ight floor} [gcd(i,j)=1] ]

因为把(leftlfloor n/p ight floor,leftlfloor m/p ight floor)以内互质的一对数乘上(p)就是(gcd=p)的数了。

令(S(n,m)=sumlimits_{i=1}^nsumlimits_{j=1}^m [gcd(i,j)=1])，由于(mu)的性质([n=1]=sumlimits_{dmid n} mu(d))，所以

[egin{aligned}S(n,m)&=sumlimits_{i=1}^n sumlimits_{j=1}^m sumlimits_{dmid gcd} mu(d) \&= sumlimits_{d}mu(d)sumlimits_{i=1}^n [dmid i]sumlimits_{j=1}^m [dmid j] \&= sumlimits_{d} mu(d)leftlfloor frac{n}{d} ight floorleftlfloor frac{m}{d} ight floorend{aligned} ]

又因为(leftlfloorfrac{leftlfloorfrac{n}{p} ight floor}{d} ight floor=leftlfloorfrac{n}{dp} ight floor)，所以柿子是

[sumlimits_{p} sumlimits_{d} mu(d) leftlfloorfrac{n}{dp} ight floor leftlfloorfrac{m}{dp} ight floor ]

可以枚举(k=dp)，有

[sumlimits_{k} leftlfloorfrac{n}{k} ight floor leftlfloorfrac{m}{k} ight floorsumlimits_{pmid k} muleft(frac{k}{p} ight) ]

令(f(n)=sumlimits_{pmid n} mu(frac{n}{p}))，预处理出(f)的前缀和，数论分块就好了。

并不会算预处理的复杂度qwq...

#include <bits/stdc++.h>
using namespace std;
#define rep(i,a,b) for (int i=(a);i<(b);++i)
#define per(i,a,b) for (int i=(a)-1;i>=(b);--i)
#define mp make_pair
#define pb push_back
#define fi first
#define se second
#define all(x) (x).begin(),(x).end()
#define SZ(x) ((int)(x).size())
typedef double db;
typedef long long ll;
typedef pair<int,int> PII;
typedef vector<int> VI;

const int maxn=1e7,N=maxn+10;
int vis[N],p[N],pn,mu[N],sum[N];

#define ss(l,r) (sum[r]-sum[l-1])

void init(int n) {
    mu[1]=1;
    rep(i,2,n+1) {
        if(!vis[i]) {p[pn++]=i;mu[i]=-1;}
        for(int j=0;j<pn&&i*p[j]<=n;j++) {
            vis[i*p[j]]=1;
            if(i%p[j]==0) {mu[i*p[j]]=0;break;}
            else mu[i*p[j]]=-mu[i];
        }
    }
    rep(i,0,pn) for(int j=p[i];j<=n;j+=p[i])
        sum[j]+=mu[j/p[i]];
    rep(i,1,n+1) sum[i]+=sum[i-1];
}

ll solve(int n,int m) {
    int tn=min(n,m); ll ans=0;
    for(int l=1,r=0;l<=tn;l=r+1) {
        r=min(n/(n/l),m/(m/l));
        ans+=1ll*(n/l)*(m/l)*ss(l,r);
    }
    return ans;
}

int main() {
#ifdef LOCAL
    freopen("a.in","r",stdin);
#endif
    init(maxn);
    int _,n,m;for(scanf("%d",&_);_;_--) {
        scanf("%d%d",&n,&m);
        printf("%lld
",solve(n,m));
    }
    return 0;
}