hdu 6053TrickGCD（线性筛+莫比乌斯函数+前缀和）

TrickGCD Time Limit: 5000/2500 MS (Java/Others)    Memory Limit: 262144/262144 K (Java/Others) Total Submission(s): 746    Accepted Submission(s): 293 Problem Description You are given an array  A  , and Zhu wants to know there are how many different array  B  satisfy the following conditions? *  1≤Bi≤Ai * For each pair( l , r ) ( 1≤l≤r≤n ) ,  gcd(bl,bl+1...br)≥2   Input The first line is an integer T( 1≤T≤10 ) describe the number of test cases. Each test case begins with an integer number n describe the size of array  A . Then a line contains  n  numbers describe each element of  A You can assume that  1≤n,Ai≤105   Output For the  k th test case , first output "Case #k: " , then output an integer as answer in a single line . because the answer may be large , so you are only need to output answer  mod   109+7   Sample Input 1 4 4 4 4 4   Sample Output Case #1: 17

首先要了解莫比乌斯函数参考网上一个简单易懂的PPT这里就不展开讲了

https://wenku.baidu.com/view/fbec9c63ba1aa8114431d9ac.html

莫比乌斯函数完整定义的通俗表达：

1）莫比乌斯函数μ(n)的定义域是N

2）μ(1)=1

3）当n存在平方因子时，μ(n)=0

4）当n是素数或奇数个不同素数之积时，μ(n)=-1

5）当n是偶数个不同素数之积时，μ(n)=1  （参考百度百科）

其次要莫比乌斯函数的求解其实是依靠线性筛来完成的，两者代码的区别仅仅莫比乌斯函数代码多了一个mu[]数组

不同线性筛的可以参考这两个博客

http://blog.csdn.net/leolin_/article/details/6642126

http://www.cnblogs.com/grubbyskyer/p/3852421.html

#include<cstdio>
#include<cstring>
using namespace std;

const int N = 1e5+4;
bool vis[N];
int prime[N],cnt,mu[N]; 

void Init()  //线性筛求莫比乌斯
{
    memset(vis,0,sizeof(vis));
    mu[1] = 1;
    cnt = 0;
    for(int i=2; i<N; i++)
    {
        if(!vis[i])   //vis用于标记是否是非素数，若是素数则为false
        {
            prime[cnt++] = i;   //i=d
            mu[i] = -1;
        }
        for(int j=0; j<cnt&&i*prime[j]<N; j++)
        {
            vis[i*prime[j]] = 1;
            if(i%prime[j]) mu[i*prime[j]] = -mu[i];   //如果这个合数(暂且说成合数，质数的情况同理)不能被prime[j]整除,prime[j]肯定不是i的因子(换句话说i中没有与prime[j]相同的因子)
            else    //如果i%prime[j]==0那么i肯定可以分解成prime[j]*(某一个数),这样到之后k=i*prime[j+1]时,k肯定可以从prime[j+1]*(某一个数)*prime[j]得到(一个更大的合数*一个更小的质数得到你)
            {
                mu[i*prime[j]] = 0;   //如果这个合数(质数)能被prime[j]整除，那么prime[j]就是i的因子，i*prime[j]中肯定有相同的质因子
                break;         //这样相比与普通筛法O(n*logn*logn)就避免了重复，使得复杂度可以降到O(n)
            }
        }
    }
}

int main()
{
	Init();
	printf("1\n");
}

了解完莫比乌斯函数，这道题的题解参考

http://www.cnblogs.com/nicetomeetu/p/7248040.html

整体思想就是因为b中每一个bi都要小于ai,有因为b要任意两个都不互质，所以通过枚举gcd来求解(因为ai并不是特别大),并且gcd一定要小于a中的最小值（易证）

(a1/k) 表示在1..a1中k的倍数的个数。[9/3=3：3，6，9]

这样对每一个gcd都算贡献= (a1/k)*(a2/k)*(a3/k)*...*(an/k),之后就是要去重的情况，因为你在算2是已经把2的倍数都考虑到了，之后再算就可能会有重复

莫比乌斯函数就是来处理这个的，根据了容斥原理

这里主要优化是前缀和sum[]数组，sum[i]表示a中不大于i的数的个数，这样将一个k的贡献值=(a1/k)*(a2/k)*(a3/k)*......转换成求1^(sum[2k-1]-sum[k-1])*2^(sum[3k-1]-sum[2k-1])*.....

就是将a中的数放进一个个k的倍数区间里，这里sum[3k-1]-sum[2k-1]表示a中整除k答案为2的数的个数，表示在原式中有sum[3k-1]-sum[2k-1]个2相乘，之后的同理。

#include<stdio.h>
#include<string.h>
#include<iostream>
#include<algorithm>
using namespace std;

typedef long long int ll;

const int MAXN = 1e5+10;
const int INF = 1999999999;
const int MOD = 1e9+7;

bool vis[MAXN];
ll prime[MAXN],cnt,mu[MAXN];  
ll a[MAXN],n,sum[MAXN];  //sum[i]表示a中不大于i的数的个数
ll Max,Min;

void Mobius()  //根据线性筛法来求莫比乌斯函数
{
    memset(vis,0,sizeof(vis));
    cnt=0;
    mu[1]=1;
    for(int i=2;i<MAXN;i++)
    {
        if(!vis[i])
        {
            prime[cnt++]=i;
            mu[i]=-1;
        }
        for(int j=0;j<cnt&&i*prime[j]<MAXN;j++)
        {
            vis[i*prime[j]]=1;
            if(i%prime[j]) mu[i*prime[j]]=-mu[i];
            else
            {
                mu[i*prime[j]]=0;
                break;
            }
        }
    }
    for(int i=2;i<MAXN;i++) mu[i]=-mu[i];
}

ll quickpow(ll b,ll k)  //快速幂
{
    ll  ans=1;
    while(k)
    {
        if(k&1)ans=(ans*b)%MOD;
        b=(b*b)%MOD;
        k>>=1;
    }
    return ans;
}


void solve()
{
    ll ans=0;
    for(ll i=2;i<=Min;i++)  //枚举k
    {
        if(!mu[i]) continue;  //符号为0是跳过
        ll j=i-1;
        ll k=2*i-1;

        ll res=1;
        for(ll p=1;;p++)    //计算每一个i对答案的贡献
        {                   //通过计算每个i的p倍区间内,a[]有几个数在这个区间内这就表示p的指数
            if(sum[k]-sum[j])    //表示在i的p倍区间内a中有sum[k]-sum[j]个数
                res=(res*quickpow(p,sum[k]-sum[j]))%MOD;  // = 1^(sum[2k-1]-sum[k-1]) * 2^(sum[3k-1]-sum[2k-1]) * 3^(sum[4k-1]-sum[3k-1]) ...
            if(k>=Max)break;
            j+=i;
            k+=i;
            if(k>Max)k=Max;
        }
        ans=(ans+mu[i]*res)%MOD;   //mu[i]表示符号
    }
    if(ans<0)ans=ans+MOD;
    printf("%lld\n",ans%MOD);

}

int main()
{
    int t,q;
    scanf("%d",&t);
    Mobius();    //求每一个k的符号
    q=0;
    while(t--)
    {
        q++;
        scanf("%lld",&n);
        Max=-1;
        Min=INF;
        memset(sum,0,sizeof(sum));
        for(ll i=1;i<=n;i++)
        {
            scanf("%lld",&a[i]);
            Max=max(a[i],Max);    //这个a中的最大值
            Min=min(a[i],Min);
            sum[a[i]]++;
        }
        sum[0]=0;
        for(ll i=1;i<=Max;i++)    //算出每一个前缀和来求贡献
            sum[i]+=sum[i-1];
        printf("Case #%d: ",q);
        solve();
    }
}