莫比乌斯反演1006 HDU 4675 好题 推荐

题意:
给一个长度为n的数组,每个数<=m,从中选出k个数
改变k个数变为数组bi
求gcd(bi)==(1~m)的方案数
思路:
求gcd的对数
Ans=∑f(d)
f(d)是gcd(bi)==d的方案数
F(d)是d|gcd(bi)的方案数
p=ai数组里拥有d因子的数的个数
F(d)=C(p,k-(n-p))((m/d)-1)^k-(n-p)(m/d)^(n-p)
因为k个数一定是选定的,我们要从拥有d因子的数里边选择数来改变
所有才有的((m/d)-1)^k-(n-p) * (m/d)^(n-p)

#include<stdio.h>
#include<string.h>
#include<iostream>
#include<algorithm>
#include<math.h>
#include<queue>
#include<stack>
#include<string>
#include<vector>
#include<map>
#include<set>
using namespace std;
#define lowbit(x) (x&(-x))
typedef long long LL;
const int maxn = 300005;
const LL MOD = 1000000007;
const int inf=(1<<28)-1;
#define maxp 1000005
bool notprimes[maxp];
int primes[maxp];
int mu[maxp];
void get_mu()
{
    memset(notprimes,false,sizeof(notprimes));
    primes[0]=0;
    mu[1]=1;
    for(int i=2;i<maxp;++i)
    {
        if(!notprimes[i])
        {
            primes[++primes[0]]=i;
            mu[i]=-1;
        }
        for(int j=1;j<=primes[0];++j)
        {
            if((LL)primes[j]*i>=maxp) break;
            notprimes[i*primes[j]]=true;
            if(i%primes[j])
                mu[i*primes[j]]=-mu[i];
            else
            {
                mu[i*primes[j]]=0;
                break;//代表i不是素数,mu[i*primes[j]]必然是0 
            }
        }
    }
}

LL A[maxn],Cnt[maxn];
LL F[maxn];
LL fac[1000010];
LL quick_pow(LL a,LL b)
{
    LL res=1,tmp=a;
    while(b)
    {
        if(b&1) res=(res*tmp)%MOD;
        tmp=(tmp*tmp)%MOD;
        b/=2;
    }
    return res;
}
LL C(LL n,LL m)
{
    if(m>n|m<0) return 0;
    LL s1=fac[n],s2=fac[n-m]*fac[m]%MOD;
    return s1*quick_pow(s2,MOD-2)%MOD;
}
void init()
{
    fac[0]=1;
    for(int i=1;i<1e6+10;++i)
    fac[i]=(fac[i-1]*i)%MOD;
}

int main()
{
    init();
    get_mu();
    int n,m,k;
    while(~scanf("%d%d%d",&n,&m,&k))
    {
        memset(Cnt,0,sizeof(Cnt));
        for(int i=1;i<=n;++i)
        {
            LL x;
            scanf("%lld",&x);
            Cnt[x]++;
        }
        for(int i=1;i<=m;++i)
        {
            for(int j=i+i;j<=m;j+=i)
                Cnt[i]+=Cnt[j];
        }
        for(int i=1;i<=m;++i)
        {
            if(k<n-Cnt[i])
                F[i]=0;
            else
            {
                LL p=Cnt[i];
                F[i]=((C(p,k-(n-p))*quick_pow(m/i-1,k-(n-p)))%MOD*quick_pow(m/i,n-p))%MOD;
            }
        }
        for(int i=1;i<=m;++i)
        {
            if(i!=1) printf(" ");
            if(F[i]==0)
            {
                printf("0");
                continue;
            }
            LL Ans=0;
            for(int j=i;j<=m;j+=i)
            {
                Ans+=(LL)mu[j/i]*F[j];
                Ans%=MOD;
            }
            Ans=(Ans+MOD)%MOD;
            printf("%lld",Ans);
        }
        printf("\n");
    }
    return 0;
}