[题解] POI2007 ZAP-Queries

题目描述

Byteasar the Cryptographer works on breaking the code of BSA (Byteotian Security Agency). He has alreadyfound out that whilst deciphering a message he will have to answer multiple queries of the form”for givenintegers aa, bb and dd, find the number of integer pairs $(x, y)​$ satisfying the following conditions:

$1\le x\le a​$, $1\le y\le b​$, $\gcd(x,y)=d​$, where $\gcd(x,y)​$ is the greatest common divisor of $x​$ and $y​$.

Byteasar would like to automate his work, so he has asked for your help.

TaskWrite a programme which:

reads from the standard input a list of queries, which the Byteasar has to give answer to, calculates answers to the queries, writes the outcome to the standard output.

输入格式

The first line of the standard input contains one integer $n​$ ($1\le n\le 50\ 000​$),denoting the number of queries.

The following $n$ lines contain three integers each: $a$, $b$ and $d$ ($1\le d\le a,b\le 50\ 000$), separated by single spaces.

Each triplet denotes a single query.

输出格式

Your programme should write $n$ lines to the standard output. The $i$’th line should contain a single integer: theanswer to the $i$’th query from the standard input.

样例输入

2
4 5 2
6 4 3

样例输出

3
2

解题思路

$$
\begin{aligned}
& \sum_{x=1}^a \sum_{y=1}^{b} [\gcd(x, y) = d] \\
= & \sum_{x=1}^{a/d} \sum_{y=1}^{b/d} [\gcd(x, y) = 1] \\
= & \sum_{x=1}^{a/d} \sum_{y=1}^{b/d} \sum_{k \mid \gcd(x, y)} \mu(k) \\
= & \sum_{x=1}^{a/d} \sum_{y=1}^{b/d} \sum_{k \mid x} \sum_{k \mid y} \mu(k) \\
= & \sum_{x=1}^{a/d} \sum_{y=1}^{b/d} \sum_{k=1}^{\min(a, b)/d} [k \mid x] [k \mid y] \mu(k) \\
= & \sum_{k=1}^{\min(a, b)/d} \mu(k) \sum_{x=1}^{a/d} \sum_{y=1}^{b/d} [k \mid x] [k \mid y] \\
= & \sum_{k=1}^{\min(a, b)/d} \mu(k) \lfloor \frac{a/d}{k} \rfloor \lfloor \frac{b/d}{k} \rfloor
\end{aligned}
$$
其中第二步到第三步是依据
$$
[x = 1] = \boldsymbol{1} * \boldsymbol {\mu} = \sum_{d \mid x} \mu(x)
$$

代码

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm>

using namespace std;
const int MAXN = 50000 + 10;
const int SIZ = 50000;

int pri[MAXN], mu[MAXN], cnt;
long long sum[MAXN];
bool isNotPrime[MAXN];

inline void init () {
    mu[1] = 1;

    for (register int i=2; i<=SIZ; i++) {
        if (!isNotPrime[i]) {
            pri[++cnt] = i;
            mu[i] = -1;
        }

        for (register int j=1; j<=cnt; j++) {
            int m = i * pri[j];     

            if (m > SIZ) {
                break;
            }

            isNotPrime[m] = true;

            if (i % pri[j] == 0) {
                mu[m] = 0;
                break;
            } else {
                mu[m] = -mu[i];
            }
        }
    }

    for (register int i=1; i<=SIZ; i++) {
        sum[i] = sum[i-1] + mu[i];
    }
}

inline long long solve (int a, int b, int d) {
    int r = 0;

    a /= d;
    b /= d;

    int minnum = min(a, b);
    long long ans = 0;

    for (register int l = 1; l <= minnum; l = r + 1) {
        r = min(a / (a / l), b / (b / l));
        ans += 1LL * (a / l) * (b / l) * (sum[r] - sum[l-1]);
    }

    return ans;
}

int n;
int a, b, d;

int main(){
    init();
    scanf("%d", &n);

    for (register int i=1; i<=n; i++){
        scanf("%d%d%d", &a, &b, &d);
        printf("%lld\n", solve(a, b, d));
    }    

    return 0;
}