BZOJ 2820:YY的GCD 莫比乌斯反演

BZOJ2301的强化版

上面这题是求gcd(i,j)==1的个数

这题是求gcd(i,j) ==p的个数 p是素数

我就在上面那个式子中额外枚举个素数p

#include <bits/stdc++.h>
#define FOR(i,s,t) for(int i=(s);i<=(t);i++)
#define ROF(i,s,t) for(int i=(s);i>=(t);i--)
#define pb push_back
#define mp make_pair
#define eb emplace_back
#define fi first
#define se second
#define endl '\n'
using namespace std;
typedef unsigned long long ull;
typedef long long ll;
const int maxn = 1e7+1;
const ll mod = 1e9 + 7;
const int inf = 0x3f3f3f3f;
const ll INF = 0x3f3f3f3f3f3f3f3f;
int readInt(){
    int x=0;
    bool sign=false;
    char c=getchar();
    while(!isdigit(c)){
        sign=c=='-';
        c=getchar();
    }
    while(isdigit(c)){
        x=x*10+c-'0';
        c=getchar();
    }
    return sign?-x:x;
}
ll readLong(){
    ll x=0;
    bool sign=false;
    char c=getchar();
    while(!isdigit(c)){
        sign=c=='-';
        c=getchar();
    }
    while(isdigit(c)){
        x=x*10+c-'0';
        c=getchar();
    }
    return sign?-x:x;
}
string readString(){
    string s;
    char c=getchar();
    while(isspace(c)){
        c=getchar();
    }
    while(!isspace(c)){
        s+=c;
        c=getchar();
    }
    return s;
}

int pri[maxn];
int vis[maxn +1];
int mu[maxn + 1];
int sum[maxn + 1];

void getPrime(){
    mu[1] = 1;
    for (int i = 2; i < maxn; i++){
        if (!vis[i]){
            pri[++pri[0]] = i;
            mu[i] = -1;
        }
        for (int j = 1; j <= pri[0] &&i * pri[j] < maxn; j++){
            vis[i * pri[j]] = 1;
            if (i % pri[j] == 0){
                mu[i * pri[j]] = 0;
                break;
            }
            mu[i * pri[j]] = -mu[i];
        }
    }
    for(int i = 1; i <= pri[0]; i++)
        for (int j = pri[i]; j <maxn; j += pri[i]){
            sum[j] += mu[j / pri[i]];
        }
    sum[0] = 0;
    for (int i = 1; i < maxn; i++){
        sum[i] = sum[i-1] + sum[i];
    }
}

ll solve(int n, int m){
    ll res = 0;
    if (n > m) swap(n, m);
    int pos = 0;
    for (int i = 1; i <= n; i = pos + 1){
        pos = min(n/(n/i), m/(m/i));
        res += 1ll * (sum[pos] - sum[i - 1]) * (n / i) * (m / i);
        //cout << i << "~" << pos << "---" <<res << endl;
    }
    return res;
}

int main(){
    getPrime();
    int T = readInt();
    while (T--){
        int n = readInt();
        int m = readInt();
        printf("%lld\n",solve(n, m));
    }
    return 0;
}