HDU 6053 TrickGCD (莫比乌斯函数+分块筛)

来源:互联网 发布:淘宝权为什么被禁赛 编辑:程序博客网 时间:2024/06/02 03:15

思路:

找到 ai 中最小的那个数,枚举他的因数作为gcd,然后根据这个gcd我们可知 ai/gcd 即为这个这个以这个数为gcd的序列对答案的贡献,然而这么算会有重复。比如我们在算2和3的时候已经把6的序列算了,所以我们需要容斥一下。

容斥的方法:莫比乌斯函数。
考虑一个数n,假设他能分解成三个质因数相乘:a*b*c,那么在算gcd为a、b、c时分别算了一遍,所以我们要减去gcd为a*b、a*c、b*c的,再加上gcd为a*b*c的。推广:我们需要加上质因数为奇数个的,减去质因数为偶数个的。当某个质因数的指数大于等于2,那么这个数我们之前肯定算过,所以这个不计入最后答案。
而莫比乌斯函数正好符合我们的需要(取个相反数就ok了)。

我们枚举gcd的时候已经花了n的时间了,所以接下来在算 ai/gcd 的时间就一定要控制在logn内。所以我们把这个式子化简:maxa/gcdj=1im

我来解释一下这个式子:就是把枚举n个 ai 转化成了枚举 ai/gcd 的值,这个点值一定在1和maxa/gcd之间,这个化简大概优化到了logn。然后就是m的含义,m代表n个 ai 中除以gcd为i的数量,这个m我们可以通过求一遍 ai 的前缀和得到(就是开个数组,把 ai 当作下标,记录 ai 出现了多少次,然后求这个数组的前缀和,这样就可以O(1)的求出一个范围内 ai 出现了多少次了,也就得到了m的值)

#include <iostream>#include <cstdio>#include <queue>#include <string.h>#include <algorithm>#define inf 0x3f3f3f3ftypedef long long int lli;using namespace std;const int mod = 1e9+7;const int maxn = 100200;bool isprime[maxn];//0代表是素数int prime[maxn],miu[maxn];void moblus(){    int cnt = 0;    miu[1] = 1;    for(lli i = 2; i < maxn; i++){        if(!isprime[i]){            prime[cnt++] = i;            miu[i] = -1;        }        for(lli j = 0; j < cnt && i*prime[j] < maxn; j++){            isprime[i*prime[j]] = 1;            if(i % prime[j])                miu[i*prime[j]] = -miu[i];            else{                miu[i*prime[j]] = 0;                break;            }        }    }}inline lli qp(lli a,lli x){    lli ans = 1;    for(;x;x>>=1){        if(x&1) ans = ans*a%mod;        a = a*a % mod;    }    return ans;}int a[102000],suma[202000];int main(){    int cas,ncase = 0,n,v,maxa,mina;    scanf("%d",&cas);    moblus();    while(cas--){        ncase++;        lli ans = 0;maxa = -100200,mina = 102000;        scanf("%d",&n);        memset(a,0,sizeof(a));        for(int i = 1;i <= n;i++){            scanf("%d",&v);            a[v]++;            mina = min(mina,v);maxa = max(maxa,v);        }        for(int i = 1;i <= maxa*2;i++){            suma[i] = suma[i-1] + a[i];        }        lli tempans = 1;        for(int i = 2;i <= mina;i++){            tempans = 1;            for(int j = 1;j*i <= maxa;j++){                tempans *= qp(j,(suma[i*j+i-1]-suma[i*j-1]));                tempans %= mod;            }            ans -= miu[i]*tempans;ans %= mod;        }        printf("Case #%d: %lld\n",ncase,(ans+mod)%mod);    }}