【bzoj4305】数列的GCD

拿来练练基础的组合数学思维还是可以的

实际上题目可以转化成：给一个序列A，修改其中的K个成为序列B，对$\forall i\in [1,m]$求出使得序列B满足$gcd\{x\mid x\in B\}==i$的方案数。

既然每个i都要check一次，那就看看对于某个i怎么计算答案$f_{i}$好了。
首先注意注意到，对于$\{x\mid x\in A \verb'&'\verb'&'i\nmid x\}$里的数必须要修改，否则gcd必定不为i。
令cnt为A内i的倍数的个数，那么称另外$n-cnt$个数为必修改数，必修改数修改方案则有$\lfloor\frac{m}{i}\rfloor^{n-cnt}$种。而对于这cnt个非必须修改数，我们可以选$n-k$个让其保持不变，然后剩余的cnt-n+k个强制改变，这样就有$\binom{n-k}{cnt}(\lfloor\frac{m}{i}\rfloor-1)^{cnt-n+k}$种选择。
还要注意到，这样子会出现gcd是i的倍数的情况，最后筛一下去掉即可。
也就是说对于$i$，答案为

$f_{i}=\binom{n-k}{cnt}(\lfloor\frac{m}{i}\rfloor-1)^{cnt-n+k}\lfloor\frac{m}{i}\rfloor^{n-cnt}-\sum_{i|d}f_{d}$

筛筛筛！
时间$O(n\log m)$，空间$O(n+m)$。

#include <bits/stdc++.h>
#define rep(i,a,b) for(int i=a;i<=b;i++)
#define per(i,a,b) for(int i=a;i>=b;i--)
#define maxn 300007

inline int rd() {
    char c = getchar();
    while (!isdigit(c)) c = getchar() ; int x = c - '0';
    while (isdigit(c = getchar())) x = x * 10 + c - '0';
    return x;
}

typedef long long ll;
typedef int arr[maxn];

const int mod = 1000000007;

inline int add(int a , int b) { ll x = (ll) a + b ; if (x >= mod) x %= mod ; return x ; }
inline int mul(int a , int b) { ll x = (ll) a * b ; if (x >= mod) x %= mod ; return x ; }
inline int dec(int a , int b) { ll x = (ll) a - b ; if (x < 0) x += mod ; if (x >= mod) x %= mod ; return x ; }

arr frac , _inv , f , mark;
int n , m , K;

inline int Pow(int a , int b) {
    int ret = 1;
    while (b) {
        if (b & 1) ret = (ll) ret * a % mod;
        a = (ll) a * a % mod , b >>= 1;
    }
    return ret;
}

inline int C(int n , int m) {
    return mul(frac[n] , mul(_inv[n - m] , _inv[m]));
}

void input() {
    n = rd() , m = rd() , K = rd();
    rep(i , 1 , n) mark[rd()] ++;
    frac[0] = 1 , _inv[0] = 1;
    rep(i , 1 , n) frac[i] = mul(frac[i - 1] , i) , _inv[i] = mul(_inv[i - 1] , Pow(i , mod - 2));
}

void solve() {
    int s = n - K;
    rep(i , 1 , m) {
        int cnt = 0;
        rep(j , 1 , m) if ((ll) j * i > m) break;
        else cnt += mark[j * i];
        if (cnt < s) continue;
        f[i] = C(cnt , s);
        f[i] = mul(f[i] , Pow(m / i - 1 , cnt - s));
        f[i] = mul(f[i] , Pow(m / i , n - cnt));
    }
    per(i , m , 1) rep(j , 2 , m) if ((ll)j * i > m) break;
    else f[i] = dec(f[i] , f[j * i]);
    rep(i , 1 , m) printf("%d%c" , f[i] , i == m ? '\n' : ' ');
}

int main() {
    #ifndef ONLINE_JUDGE
        freopen("data.txt" , "r" , stdin);
    #endif
    input();
    solve();
    return 0;
}