GuGuFishtion HDU

GuGuFishtion (HDU - 6390)

题意：

定义(G_u (a,b)=frac{phi(ab)}{phi(a)phi(b)})。

求((sumlimits_{a=1}^msumlimits_{b=1}^nG_u (a,b))pmod p)。

题解：

考虑(phi(x) = x*(1-frac{1}{p_1})*(1-frac{1}{p_2})...*(1-frac{1}{p_n}))。

将(G_u (a,b))的分子与分母按上述分解、约分后，得(frac{1}{(1-frac{1}{p_1})(1-frac{1}{p_2})(1-frac{1}{p_3})...(1-frac{1}{p_k})})，其中(p_1*p_2*p_3...p_k=gcd(a,b))。

上下同时乘以(gcd(a,b))，则(G_u (a,b)=frac{gcd(a,b)}{phi(gcd(a,b))})。

问题变成了求((sumlimits_{a=1}^msumlimits_{b=1}^nfrac{gcd(a,b)}{phi(gcd(a,b))})pmod p)。

上式 (=(sum_d^{min(n, m)}frac{d}{phi(d)}sumlimits_{a=1}^msumlimits_{b=1}^n[gcd(a,b)==d])pmod p)。

可以枚举(d)，对于(sumlimits_{a=1}^msumlimits_{b=1}^n[gcd(a,b)==d])，即(sumlimits_{a=1}^{m/d}sumlimits_{b=1}^{n/d}[gcd(a,b)==1])，用莫比乌斯反演求。

代码：

#include <bits/stdc++.h>
#define fopi freopen("in.txt", "r", stdin)
#define fopo freopen("out.txt", "w", stdout)
using namespace std;
typedef long long LL;
typedef long double ld;
const int maxn = 1e6 + 10;

int check[maxn], phi[maxn], mu[maxn], prime[maxn];
LL inv[maxn], sum[maxn];

void init(int N) {
    memset(check, false, sizeof(check));
    mu[1] = 1;
    phi[1] = 1;
    int tot = 0;
    for (int i = 2; i <= N; i++) {
        if (!check[i]) {
            prime[tot++] = i;
            mu[i] = -1;
            phi[i] = i-1;
        }
        for (int j = 0; j < tot; j++) {
            if (i * prime[j] > N) break;
            check[i * prime[j]] = true;
            if (i % prime[j] == 0) {
                mu[i * prime[j]] = 0;
                phi[i * prime[j]] = phi[i] * prime[j];
                break;
            }
            else {
                mu[i * prime[j]] = -mu[i];
                phi[i * prime[j]] = phi[i] * (prime[j]-1);
            }
        }
    }
}

LL solve(int n, int m, int p) {
    LL res = 0;
    for (int i = 1, la = 0; i <= m; i = la+1) {
        la = min(n/(n/i), m/(m/i));
        res = (res + 1ll * (sum[la] - sum[i-1]) * (n/i) % p * (m/i) % p) % p;
    }
    return res;
}

int T, n, m, p;
int main() {
    init(maxn - 10);
    
    scanf("%d", &T);
    while(T--) {
        scanf("%d%d%d", &n, &m, &p);
        if (n < m) swap(n, m);

        inv[1] = 1;
        for (int i = 2; i <= m; i++)
            inv[i] = 1ll * inv[p % i] * (p - p/i) % p;
        for (int i = 1; i <= m; i++)
            sum[i] = (sum[i-1] + mu[i] + p) % p;

        LL ans = 0;
        for (int i = 1; i <= m; i++) {
            ans = (ans + 1ll * i * inv[phi[i]] % p * solve(n/i, m/i, p) % p) % p;
        }

        printf("%lld
", ans);
    }
}