hdu 6053 TrickGCD

Problem

acm.hdu.edu.cn/showproblem.php?pid=6053

莫比乌斯反演介绍

莫比乌斯反演（by PoPoQQQ）
莫比乌斯反演入门
莫比乌斯反演

Reference

2017 Multi-University Training Contest 2 solutions BY 电子科技大学（1009题）
HDU 6053 TrickGCD （莫比乌斯函数）
解题报告：HDU_6053 TrickGCD 莫比乌斯反演

Meaning

给定一个序列 $\{A_{n}\}$，要构造一个序列 $\{B_{n}\}$，满足：

• $1\leq B_{i}\leq A_{i}$
• $gcd（B_{1}，…，B_{n}）\geq 2$

问构造的方案数。

Analysis

暴力的做法，就是枚举 gcd，每一个 gcd 产生的答案是 $\prod^{n}_{i=1} \frac{A_{i}}{gcd}$，总的答案就是它们的求和，但要去重，可以用莫比乌斯反演。
…好吧其实现在并不太懂莫比乌斯反演…按官方题解，定义：
F(x)：gcd 是 x 倍数的方案数；
我猜相应的就有：
f(x)：gcd 是 x 的方案数；
于是：
$F(x) = \sum_{x | d}f(d) \Rightarrow f(x) = \sum_{x | d} \mu (\frac{d}{x})F(d)$
而 $F(x) = \prod^{n}_{i=1} \frac{A_{i}}{x}$，就可以做了。
按 dalao 的说法，因为枚举 gcd 时是从 2 开始的，f(1) 没有用到，所以在乘 $\mu (i)$的时候，要多乘一个 -1，于是就有了第二个参考博客的反演式。
在处理那个连乘式的时候，要用一点优化处理。
好吧还是很迷…

Code

#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;
const int N = 100000, A = N, MOD = 1000000007;

int a[N], sum[N<<1|1];
int prime[A+1];
int mu[A+1];

void Mobius()
{
    mu[1] = 1;
    memset(prime, -1, sizeof prime);
    for(int i = 2, top = 0; i <= A; ++i)
    {
        if(prime[i])
        {
            prime[top++] = i;
            mu[i] = -1;
        }
        for(int j = 0; j < top && i * prime[j] <= A; ++j)
        {
            prime[i * prime[j]] = 0;
            if(i % prime[j])
                mu[i * prime[j]] = -mu[i];
            else
            {
                mu[i * prime[j]] = 0;
                break;
            }
        }
    }
}

int fast_pow(long long a, int b)
{
    long long x = 1;
    for( ; b; b >>= 1, a = a * a % MOD)
        if(b & 1)
            x = x * a % MOD;
    return x % MOD;
}

int main()
{
    Mobius();
    int T;
    scanf("%d", &T);
    for(int kase = 1; kase <= T; ++kase)
    {
        int n, small = A, big = 0;
        scanf("%d", &n);
        memset(sum, 0, sizeof sum);
        for(int i = 0; i < n; ++i)
        {
            scanf("%d", a+i);
            ++sum[a[i]];
            small = min(small, a[i]);
            big = max(big, a[i]);
        }
        for(int i = 1; i <= big << 1; ++i)
            sum[i] += sum[i-1];
        long long ans = 0, tmp;
        for(int d = 2; d <= small; ++d)
            if(tmp = -mu[d]) // mu[d] != 0
            {
                for(int i = d; i <= big; i += d)
                    tmp = tmp * fast_pow(i/d, sum[i+d-1]-sum[i-1]) % MOD;
                ans = (ans + tmp + MOD) % MOD;
            }
        printf("Case #%d: %I64d\n", kase, ans);
    }
    return 0;
}