HDU 6053 - TrickGCD | 2017 Multi-University Training Contest 2

/*

HDU 6053 - TrickGCD [ 莫比乌斯函数,筛法分块 ]  | 2017 Multi-University Training Contest 2

题意：

	给出数列 A[N]，问满足：

		1 <= B[i] <= A[i] ;

		对任意(l, r) (1<=l<=r<=n) ，使得 gcd(bl,...br) >= 2 ;

	的 B[N] 数列的个数

分析：

	设 gcd(b1,...bn) = k (k >= 2)，此时 k 对答案的贡献为 (a1/k)*(a2/k)*(a3/k)*...*(an/k)

	根据容斥原理，ans = +[k=一个素数之积 时对答案的贡献]

						-[k=两个素数之积 时对答案的贡献]

						+[k=三个素数之积 时对答案的贡献]

						...

	故任意k对答案的贡献系数 μ(k) = 0 , k是完全平方数的倍数

								 = (-1)^(n-1) , k = p1*p2*p3*...*pn ,p是素数

	贡献系数可以O(nsqrt(n)) 或者 O(nlogn) 预处理，再或者可以看出μ(k) 是莫比乌斯函数的相反数

	现在枚举k需要O(n)的时间，计算k对答案的贡献必须在O(sqrt(n))的时间之内

	将a[]处理成权值数组，并求前缀和，设为 sum[]

	对于每个k，对sum[]进行埃式筛法的分块，即根据k的倍数分块

	此时每个k的贡献 = 1^(sum[2k-1]-sum[k-1]) * 2^(sum[3k-1]-sum[2k-1]) * 3^(sum[4k-1]-sum[3k-1]) ...

	就做到 O(n(logn)^2)

编码时长： 40分钟(0)

*/

#include  <bits/stdc++.h>

using namespace std;

#define LL long long

const LL MOD = 1e9+7;

const int N = 1e5+4;

bool notp[N];

int prime[N], pnum, mu[N];

void Mobius() {

	memset(notp, 0, sizeof(notp));

	mu[1] = 1;

	for (int i = 2; i < N; i++) {

		if (!notp[i]) prime[++pnum] = i, mu[i] = -1;

		for (int j = 1; prime[j]*i < N; j++) {

			notp[prime[j]*i] = 1;

			if (i%prime[j] == 0) {

				mu[prime[j]*i] = 0;

				break;

			}

			mu[prime[j]*i] = -mu[i];

		}

	}

	for (int i = 0; i < N; i++) mu[i] = -mu[i];

}

LL PowMod(LL a, int m)

{

    if (a == 1 || m == 0) return 1;

    if (a == 0) return 0;

    LL res = 1;

    while (m)

    {

        if (m&1) res = res*a % MOD;

        a = a*a % MOD;

        m >>= 1;

    }

    return res;

}

int a[N];

int t, n;

int sum[N];

int Max, Min;

LL ans;

void solve()

{

    int i, j, k, p;

    ans = 0;

    for (i = 2; i <= Min; ++i)

    {

        if (!mu[i]) continue;

        LL res = 1;

        j = min(i, Max), k = min((i<<1)-1, Max);

        for (p = 1; ; ++p)

        {

            if (sum[k] - sum[j-1])

                res = res*PowMod(p, sum[k] - sum[j-1]) % MOD;

            if (k == Max) break;

            j += i;

            k += i;

            if (k > Max) k = Max;

        }

        ans += mu[i]*res;

        if (ans > MOD) ans -= MOD;

        if (ans < 0) ans += MOD;

    }

}

int main()

{

    int i;

    Mobius();

    scanf("%d", &t);

    for (int tt = 1; tt <= t; ++tt)

    {

        scanf("%d", &n);

        for (i = 0; i < N; ++i) sum[i] = 0;

        for (i = 1; i <= n; ++i)

        {

            scanf("%d", &a[i]);

            ++sum[a[i]];

        }

        Max = Min = a[1];

        for (i = 2; i <= n; ++i)

        {

            Max = max(Max, a[i]);

            Min = min(Min, a[i]);

        }

        for (i = 1; i <= Max; ++i) sum[i] += sum[i-1];

        solve();

        printf("Case #%d: %lld\n", tt, ans);

    }

}