HDU 4675 GCD of Sequence（莫比乌斯反演 + 打表注意事项）题解

Posted 2022-09-14 kirinsb

题意：

给出\\(M\\)和\\(a数组\\)，询问每一个\\(d\\in[1,M]\\)，有多少组数组满足：正好修改\\(k\\)个\\(a\\)数组里的数使得和原来不同，并且要\\(\\leq M\\)，并且\\(gcd(a_1,a_2,\\dots,a_n)=d\\)。

思路：

对于每一个\\(d\\)，即求\\(f(d)\\)：修改\\(k\\)个后\\(gcd(a_1,a_2,\\dots,a_n)=d\\)的对数。
那么假设\\(F(d)\\)：修改\\(k\\)个后\\(gcd(a_1,a_2,\\dots,a_n)\\)是\\(d\\)倍数的对数。故：
\\[ f(k) = \\sum_k|d\\mu(\\fracdk)F(d) \\]
打表求\\(F(d)\\)即可。假设\\(num[d]\\)为\\(a\\)中是\\(d\\)倍数的数量，则
\\[ F(d)=(\\fracMd)^n-num[d]*C_num[d]^k-(n-num[d])*(\\fracMd-1)^k-(n-num[d]) \\]
然后\\(nlogn\\)打出\\(num\\)数组即可。

思考：

这样的打表法是\\(nlogn\\)的：
证明 O(n/1+n/2+…+n/n)=O(nlogn)

for(int i = 1; i <= n; i++)
     scanf("%d", &a[i]);
     num[a[i]]++;
 
 for(int i = 1; i <= m; i++)
     for(int j = i + i; j <= m; j += i)
         num[i] += num[j];
     
 

这样是\\(n\\sqrt n\\)的

for(int i = 1; i <= n; i++)
     scanf("%d", &a[i]);
     for(int j = 1; j <= sqrt(a[i]); j++)
         if(a[i] % j == 0)
             num[j]++;
             if(j * j != a[i]) num[a[i] / j]++;
         
     
 

代码：

#include<map>
#include<set>
#include<queue>
#include<stack>
#include<ctime>
#include<cmath>
#include<cstdio>
#include<string>
#include<vector>
#include<cstring>
#include<sstream>
#include<iostream>
#include<algorithm>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
const int maxn = 3e5 + 5;
const int INF = 0x3f3f3f3f;
const ull seed = 131;
const ll MOD = 1000000007;
using namespace std;

int mu[maxn], vis[maxn];
int prime[maxn], cnt;
ll fac[maxn], inv[maxn];
ll ppow(ll a, ll b)
    ll ret = 1;
    while(b)
        if(b & 1) ret = ret * a % MOD;
        a = a * a % MOD;
        b >>= 1;
    
    return ret;

void init(int n)
    memset(vis, 0, sizeof(vis));
    memset(mu, 0, sizeof(mu));
    cnt = 0;
    mu[1] = 1;
    for(int i = 2; i <= n; i++) 
        if(!vis[i])
            prime[cnt++] = i;
            mu[i] = -1;
        
        for(int j = 0; j < cnt && prime[j] * i <= n; j++)
            vis[prime[j] * i] = 1;
            if(i % prime[j] == 0) break;
            mu[i * prime[j]] = -mu[i];
        
    

    fac[0] = inv[0] = 1;
    for(int i = 1; i <= n; i++) fac[i] = fac[i - 1] * i % MOD;
    inv[n] = ppow(fac[n], MOD - 2);
    for(int i = n - 1; i >= 1; i--) inv[i] = (i + 1LL) * inv[i + 1] % MOD;

ll C(int n, int m)
    return fac[n] * inv[m] % MOD * inv[n - m] % MOD;

int num[maxn], a[maxn];
//num[i]:是i的倍数的个数
ll F[maxn], f[maxn];
int main()
    init(3e5);
    int n, m, k;
    while(~scanf("%d%d%d", &n, &m, &k))
        memset(num, 0, sizeof(num));
        for(int i = 1; i <= n; i++)
            scanf("%d", &a[i]);
            num[a[i]]++;
        
        for(int i = 1; i <= m; i++)
            for(int j = i + i; j <= m; j += i)
                num[i] += num[j];
            
        
        for(int i = 1; i <= m; i++)
            int no = n - num[i];
            if(no > k) F[i] = 0;
            else
                F[i] = ppow(m / i, no) * C(num[i], k - no) % MOD * ppow(m / i - 1, k - no) % MOD;
            
        

        for(int i = 1; i <= m; i++)
            f[i] = 0;
            for(int j = i; j <= m; j += i)
                f[i] += mu[j / i] * F[j];
                f[i] %= MOD;
            
            printf("%lld%c", (f[i] % MOD + MOD) % MOD, i == m? '\\n' : ' ');
        

    
    return 0;