【BZOJ4305】数列的GCD

Description
给出一个长度为N的数列{a[n]}，1<=a[i]<=M(1<=i<=N)。
现在问题是，对于1到M的每个整数d，有多少个不同的数列b[1], b[2], …, b[N]，满足：
(1)1<=b[i]<=M(1<=i<=N)；
(2)gcd(b[1], b[2], …, b[N])=d；
(3)恰好有K个位置i使得a[i]<>bi
注：gcd(x1,x2,…,xn)为x1, x2, …, xn的最大公约数。
输出答案对1,000,000,007取模的值。

Input
第一行包含3个整数，N，M，K。
第二行包含N个整数：a[1], a[2], …, a[N]。

Output
输出M个整数到一行，第i个整数为当d=i时满足条件的不同数列{b[n]}的数目mod 1,000,000,007的值。

Sample Input
3 3 3
3 3 3
Sample Output
7 1 0
HINT

当d=1，{b[n]}可以为：(1, 1, 1), (1, 1, 2), (1, 2, 1), (1, 2, 2), (2, 1, 1), (2, 1, 2), (2, 2, 1)。

当d=2，{b[n]}可以为：(2, 2, 2)。

当d=3，因为{b[n]}必须要有k个数与{a[n]}不同，所以{b[n]}不能为(3, 3, 3)，满足条件的一个都没有。

对于100%的数据，1<=N,M<=300000, 1<=K<=N, 1<=a[i]<=M。

考虑正难则反,既然是有k个位置不同,那么我们可以试着求另外n-k个位置相同
根据题意枚举d
对于每个d,枚举倍数(<=m)计算a数列中有多少个数满足是d的倍数.
我们假设共有sum个数是d的倍数,那么对于当前的d,可以通过排列组合和容斥原理计算得到答案为:

$$ans[i]=C_{sum}^{n-k}*(\lfloor \frac m d\rfloor-1)^{sum-(n-k)}*\lfloor \frac m d \rfloor^{n-sum}-\sum_{i=1}^{\lfloor \frac m d \rfloor } ans[i*d]$$
膜Claris讲解…膜Po姐题解…

#include<iostream>
#include<cstdio>
#include<cstring>
#include<cmath>
#include<algorithm>
#define LL long long
#define MAXN 300100
#define P 1000000007
#define GET (ch>='0'&&ch<='9')
using namespace std;
int n,m,k;
int a[MAXN];
LL fac[MAXN],inv[MAXN];
int cnt[MAXN],ans[MAXN];
LL Pow(LL x,int _)
{
    LL ret=1;
    for (int i=_;i;i>>=1,x=x*x%P)   if (i&1)    ret*=x,ret%=P;
    return ret;
}
void in(int &x)
{
    char ch=getchar();x=0;
    while (!GET)    ch=getchar();
    while (GET) x=x*10+ch-'0',ch=getchar();
}
LL C(int n,int m)
{
    return fac[n]*inv[m]%P*inv[n-m]%P;
}
int main()
{
    in(n);in(m);in(k);k=n-k;fac[0]=1;int sum=0;
    for (int i=1;i<=n;i++)  in(a[i]),cnt[a[i]]++;
    for (int i=1;i<=n;i++)  fac[i]=fac[i-1]*i%P;
    inv[1]=1;for (int i=2;i<=n;i++) inv[i]=(P-P/i)*inv[P%i]%P;
    inv[0]=1;for (int i=1;i<=n;i++) inv[i]=inv[i]*inv[i-1]%P;
    for (int i=1;i<=m;i++)
    {
        sum=0;
        for (int j=i;j<=m;j+=i) sum+=cnt[j];
        if (sum<k)  {ans[i]=0;continue;}
        ans[i]=C(sum,k)*Pow(m/i-1,sum-k)%P*Pow(m/i,n-sum)%P;
    }
    for (int i=m;i;i--)
        for (int j=i+i;j<=m;j+=i)   ans[i]=(ans[i]+(P-ans[j]))%P;
    for (int i=1;i<m;i++)   printf("%d ",ans[i]);printf("%d\n",ans[m]);
}