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]<>b[i](1<=i<=N) 

 

注：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。 

 

 

---

 

 

 

考虑当 d 为公倍数，但不一定为最小公倍数时的方案：

因为有 k 个数要不同

设 s = n-k ，表示有 s 个数要与原数列相同

当 d 为公倍数时，设 a 数列中有 sum 个数为 d 的倍数

因为对于与原数列相同的那 s 个数，它们必须都为 d 的倍数

所有一定要从 sum 个 d 的倍数的数中取出 s 个数

从 sum 个 d 的倍数的数中取出 s 个数的方案数为 C(s,cum) （C(n,m)表示从m个数中取出n个数的方案数）

然后考虑与原数列不相同的数，显然每一个数只要取 d 的倍数就好了

那么每个数取的方案有 [m/d] 种 （[a/b] 表示 a/b 向下取整），（可以取每一个 ≤ m 的 d 的倍数）

但是要注意，对于那些原来就是 d 的倍数，但没被选的数，它们不能取原来的值

所以对于这些数，它们的取数方案只有 [m/d]-1 种 （不能取自己原来的值）

所以总方案数 ans [ d ] = C(s,sum) * ([m/d]-1)^(sum-s)  *  ([m/d])^(n-sum) 

即

原数列相同的那 s 个数的方案数 * 对于那些原来就是 d 的倍数但没被选的数的方案数 * 原本不是d的倍数的数的方案数

 

然后再考虑当 d 为最小公倍数时的情况：

设 Ans[ i ] 表示 i 为最小公倍数时的方案

那么 Ans [ d ] = ans[ d ] - Ans[ d*i ]  （2<= i && d*i <= m）

这个用反向的递推可以全部求出来

关于 C 的计算可以先预处理一波阶乘

那么 C(n,m) = m!/( n! * (m-n)! )

但是因为要取模质数 = 1e9+7

可以用费马小定理求出 x 模P意义下的逆元 = x^(P-2) （P是质数）

然后就可以搞了，注意开long long

复杂度大概为 O(nlogP)

具体看代码：

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cstring>
#include<cmath>
using namespace std;
typedef long long ll;
const int N=5e5+7,mo=1e9+7;
int n,m,k;
ll f[N],ans[N];//f是阶乘,ans是答案
int cnt[N];//cnt[i]表示值为i的数的个数
inline ll ksm(ll x,ll y)//快速幂
{
    ll res=1;
    while(y)
    {
        if(y&1) res=(res*x)%mo;
        x=x*x%mo; y>>=1;
    }
    return res;
}
inline ll C(ll n,ll m){ return f[m]*ksm(f[n],mo-2)%mo*ksm(f[m-n],mo-2)%mo; }//求C(n,m)
int main()
{
    int a;
    cin>>n>>m>>k;
    for(int i=1;i<=n;i++) scanf("%d",&a),cnt[a]++;
    f[0]=1;//注意f[0]=1
    for(int i=1;i<=n;i++) f[i]=f[i-1]*i%mo;
    int s=n-k;
    for(int i=m;i;i--)//从后往前算
    {
        int sum=0;
        for(int j=i;j<=m;j+=i) sum+=cnt[j];//计算sum
        if(sum<s){ ans[i]=0; continue; }//如果原数列还没有s个数为d的倍数，那么方案不存在
        ans[i]=C(s,sum)*ksm(m/i-1,sum-s)%mo*ksm(m/i,n-sum)%mo;//先算出d为公倍数但不一定是最小公倍数时的情况
        //然后再减去不是最小公倍数的情况
        for(int j=i*2;j<=m;j+=i) ans[i]=(ans[i]+mo-ans[j])%mo;//这时后面的ans都是处理好的
    }
    for(int i=1;i<=m;i++) printf("%lld ",ans[i]);
    return 0;
}