HDU 4675 - GCD of Sequence

解题思路参考：http://www.cnblogs.com/kuangbin/archive/2013/08/13/3255943.html

如果gcd(b[1...n])=d，那么b[1...n]中每个b[i]都应是d的倍数，

而b[1...n]是从a[1...n]修改k个数进而变过来的，

那么，假设a[1...n]中有cnt个a[i]本来就是d的倍数，

那么首先 ( n-cnt )个不是d的倍数的必须修改，

然后，剩下的 cnt 个本来就是d的倍数的a[i]中，我们再挑选 k - ( n - cnt ) 个修改（k个要修改的数，已经修改了 n-cnt 个），就可以满足条件，

而在小于m的情况下，d的倍数有 l = floor(m/d) 个，

所以

又

所以

接下来是关于逆元的知识：http://www.cnblogs.com/linyujun/p/5194184.html

#include<cstdio>
#include<cstring>
#define MOD 1000000007
#define MAX 300005
using namespace std;
typedef long long ll;
ll C[MAX],ans[MAX];
int a[MAX],b[MAX],num[MAX];
ll pow(ll a,ll b){
	ll r=1,base=a%MOD;
	while(b){
		if(b&1) r*=base , r%=MOD;
		base*=base;
		base%=MOD;
		b>>=1;
	}
	return r;
}
ll inv(ll a,ll p){return pow(a, p-2);}//费马求a关于b的逆元
int main()
{
	int n,m,k;
    while(scanf("%d%d%d",&n,&m,&k)!=EOF)
    {
        memset(num,0,sizeof(num));
        C[n-k]=1;//当cnt=n-k时，C[cnt][n-k]=C[n-k][n-k]=1 
        for(int i = n-k+1;i <= n;i++)//当cnt = n-k+1 to n 时，计算C[cnt][n-k]
        {
            C[i] = C[i-1]*i%MOD*inv(i-(n-k),MOD)%MOD;
        }
        for(int i=1;i<=n;i++)
        {
            scanf("%d",&a[i]);
            num[a[i]]++;
        }
        for(int d=m;d>=1;d--)
        {
            int cnt=0;
            ll accu=0;
            for(int times=1;times * d <= m;times++)
            {
                cnt+=num[times * d];
                if(times > 1) accu = (accu + ans[times *d])%MOD;
            }
            int l=m/d;
            if(l==1)
            {
                if(k - (n - cnt) == 0) ans[d]=1;//如果正好有k个不是d的倍数 
                else ans[d]=0;//如果k>(n-cnt)，那么剩下的就要在cnt里面找数进行修改，但是这时l==1，因此修改不了；或者k<(n-cnt)，这时修改了n-cnt个数，已经超过了k个，不满足 
            }
            else
            {
                if(k < n - cnt) ans[d]=0;
                else
                {
                    ll ans_tmp=1;
                    ans_tmp=ans_tmp*C[cnt]%MOD;
                    ans_tmp=ans_tmp*pow(l-1,k-n+cnt)%MOD;
                    ans_tmp=ans_tmp*pow(l,n-cnt)%MOD;
                    ans[d]=(ans_tmp-accu+MOD)%MOD;
                }
            }
        }
        for(int d=1;d<=m;d++)
        {
            if(d>1) printf(" ");
            printf("%I64d",ans[d]);
        }
        printf("\n");
    }
}