LuoguP3455 [POI2007]ZAP - Queries

题目描述

给定$n,m,k$,求以下式子的值:

$$\sum_{i=1}^n \sum_{j=1}^m [\gcd(i,j) = k]$$

Input & Output's examples

Input's examples

2
4 5 2
6 4 3

Output's examples

3
2

分析

比较显然的Mobius反演题目，没学过莫反的别尝试暴力了

我们来一起推一下式子(注意：以下的分数默认向下取整)。

$$\sum_{i=1}^n \sum_{j=1}^m [\gcd(i,j) = k]$$

同时除一下$k$:

$$\sum_{i=1}^{\frac nk} \sum_{j=1}^{\frac mk} [\gcd(i,j) = 1]$$

由Mobius函数的定义式：

$$\sum_{d | n} \mu(d) = [n = 1]$$

我们把$gcd(i,j)$代入可得。

$$\sum_{d|gcd(i,j)} \mu(d) = [\gcd(i,j) = 1]$$

所以，原式可以简化为：

$$\sum_{i=1}^{\frac nk} \sum_{j=1}^{\frac mk} \sum_{d | \gcd(i,j)} \mu(d)$$

我们变换枚举顺序，首先枚举$d$，将它的求和符号提前，后面内容则变为布尔表达式

$$\sum_{d=1}^{\min \{ \frac nk,\frac mk \}} \sum_{i=1}^{\frac nk} \sum_{j=1}^{\frac mk} [d | \gcd(i,j)] \mu(d)$$

因为$d | \gcd(i,j)$，那么我们只需要枚举$d$的倍数即可，同时除一下$d$:

$$\sum_{d=1}^{\min \{ \frac nk,\frac mk \}} \sum_{i=1}^{\frac {n}{kd}} \sum_{j=1}^{\frac {m}{kd}} [\gcd(i,j) \in Z] \mu(d)$$

显然$\gcd(i,j) \in Z$恒成立，那么直接去掉。

$$\sum_{d=1}^{\min \{ \frac nk,\frac mk \}} \sum_{i=1}^{\frac {n}{kd}} \sum_{j=1}^{\frac {m}{kd}} \mu(d)$$

$i,j$已经消失，那么我们直接根据乘法分配率：

$$\sum_{d=1}^{\min \{ \frac nk,\frac mk \}} \lfloor \frac {n}{kd} \rfloor \lfloor \frac {m}{kd} \rfloor \mu(d)$$

因为$\mu(d)$是不完全积性函数，我们可以$O(n)$求和。

但是本题要求多组数据,于是我们可以数论分块一下，然后预处理$\mu$前缀和。

代码

/* Headers */
#include<cstdio>
#include<cstring>
#include<cmath>
#include<cctype>
#include<algorithm>
#include<vector>
#include<queue>
#include<stack>
#include<climits>
#include<iostream>
#include<map>
#define FOR(i,a,b,c) for(int i=(a);i<=(b);i+=(c))
#define ROF(i,a,b,c) for(int i=(a);i>=(b);i-=(c))
#define FORL(i,a,b,c) for(long long i=(a);i<=(b);i+=(c))
#define ROFL(i,a,b,c) for(long long i=(a);i>=(b);i-=(c))
#define FORR(i,a,b,c) for(register int i=(a);i<=(b);i+=(c))
#define ROFR(i,a,b,c) for(register int i=(a);i>=(b);i-=(c))
#define lowbit(x) x&(-x)
#define LeftChild(x) x<<1
#define RightChild(x) (x<<1)+1
#define RevEdge(x) x^1
#define FILE_IN(x) freopen(x,"r",stdin);
#define FILE_OUT(x) freopen(x,"w",stdout);
#define CLOSE_IN() fclose(stdin);
#define CLOSE_OUT() fclose(stdout);
#define IOS(x) std::ios::sync_with_stdio(x)
#define Dividing() printf("-----------------------------------\n");
namespace FastIO{
    const int BUFSIZE = 1 << 20;
    char ibuf[BUFSIZE],*is = ibuf,*its = ibuf;
    char obuf[BUFSIZE],*os = obuf,*ot = obuf + BUFSIZE;
    inline char getch(){
        if(is == its)
            its = (is = ibuf)+fread(ibuf,1,BUFSIZE,stdin);
        return (is == its)?EOF:*is++;
    }
    inline int getint(){
        int res = 0,neg = 0,ch = getch();
        while(!(isdigit(ch) || ch == '-') && ch != EOF)
            ch = getch();
        if(ch == '-'){
            neg = 1;ch = getch();
        }
        while(isdigit(ch)){
            res = (res << 3) + (res << 1)+ (ch - '0');
            ch = getch();
        }
        return neg?-res:res;
    }
    inline void flush(){
        fwrite(obuf,1,os-obuf,stdout);
        os = obuf;
    }
    inline void putch(char ch){
        *os++ = ch;
        if(os == ot)    flush();
    }
    inline void putint(int res){
        static char q[10];
        if(res==0)  putch('0');
        else if(res < 0){putch('-');res = -res;}
        int top = 0;
        while(res){
            q[top++] = res % 10 + '0';
            res /= 10;
        }
        while(top--)    putch(q[top]);
    }
    inline void space(bool x){
        if(!x) putch('\n');
        else putch(' ');
    }
}
inline void read(int &x){
    int rt = FastIO::getint();
    x = rt;
}
inline void print(int x,bool enter){
    FastIO::putint(x);
    FastIO::flush();
    FastIO::space(enter);
}
/* definitions */
const int MAXN = 5e4 + 10;
int n,m,d,u[MAXN];
int T,prime[MAXN],cnt;
bool isprime[MAXN];
/* functions */
inline long long solve(){
    long long ans = 0;
    scanf("%d%d%d",&n,&m,&d);
    if(n > m) std::swap(n,m);
    n /= d, m /= d;
    for(int i=1,j;i <= n;i = j + 1){
        j = std::min(n / (n / i), m / (m / i));
        ans += (long long) (u[j] - u[i - 1]) * (n / i) * (m / i);
    }
    return ans;
}
inline void init(){
    u[1] = 1;
    FOR(i,2,MAXN,1){
        if(!isprime[i]) u[i] = -1, prime[++cnt] = i;
        for(int j=1;j <= cnt && (long long) i * prime[j] <= MAXN;++j){
            isprime[i * prime[j]] = true;
            if(i % prime[j] == 0) break;
            u[i * prime[j]] = -u[i];
        }
        u[i] +=  u[i-1];
    }
}
int main(int argc,char *argv[]){
    init();
    scanf("%d",&T);
    while(T--) printf("%lld\n",solve());
    return 0;
}

THE END