【BZOJ 4305】数列的GCD

传送门

---

problem

给出一个长度为 $n$ 的数列 { a n } \{a_n\} {an​}，$1\le a_i\le m(1\le i\le n)$。

现在问题是，对于 $1$ 到 $m$ 的每个整数 $d$，有多少个不同的数列 $b_1,b_2,...,b_n$，满足：

• $1\le b_i\le m(1\le i\le n)$；
• $\gcd (b_1,b_2,...,b_n)=d$；
• 恰好有 $k$ 个位置 $i$ 使得 $a_i̸=b_i(1\le i\le n)$。

注：$\gcd (x_1,x_2,...,x_n)$ 为 $x_1,x_2,...,x_n$ 的最大公约数。

输出答案对 $10^9+7$ 取模的值。

数据范围：$1\le n,m\le 3\times 10^5$，$1\le k\le n$，$1\le a_i\le m$。

---

solution

设 $f(d)$ 表示 $\gcd (b_1,...,b_n)=d$ 的答案，$F(d)={\sum }_{d|D}f(D)$。

那么我们可以先计算出 $F(d)$，然后用莫比乌斯反演得到 $f(d)$。

我们用 $cn{t}_d$ 表示在序列 { a n } \{a_n\} {an​} 中，是 $d$ 倍数的数的个数。先把 $F(d)$ 的表达式放上来：

$$F(d)=\left(\genfrac{}{}{0ex}{}{cn{t}_d}{n-k}\right)(\lfloor \frac{m}{d}\rfloor -1{)}^{cn{t}_d-(n-k)}(\lfloor \frac{m}{d}\rfloor {)}^{n-cn{t}_d}$$

接下来就是对这个式子的理解。

首先，由于序列 { b n } \{b_n\} {bn​} 的 $\gcd =d$，所以每个数都必须是 $d$ 的倍数。即如果不加任何限制，每个位置的选法有 $\lfloor \frac{m}{d}\rfloor$ 种。

我们要求要有恰好 $k$ 个位置使 $a_i̸=b_i$，就要选恰好 $(n-k)$ 个位置使 $a_i=b_i$。由于只有 $cn{t}_d$ 个位置上可能有 $a_i=b_i$，所以我们可以从 $cn{t}_d$ 中选中 $n-k$ 个，钦定这些位置相等，那就是 $\left(\genfrac{}{}{0ex}{}{cn{t}_d}{n-k}\right)$。

剩下的 $cn{t}_d-(n-k)$ 个位置上，不能有位置相等，所以只有 $\lfloor \frac{m}{d}\rfloor -1$ 种选法，就是 $(\lfloor \frac{m}{d}\rfloor -1{)}^{cn{t}_d-(n-k)}$。

最后的 $n-cn{t}_d$ 个位置，$a_i$ 和 $b_i$ 已经没有机会相等了，这部分的贡献就是 $(\lfloor \frac{m}{d}\rfloor {)}^{n-cn{t}_d}$。

那么最后反演一下，得到：

$$f(x)=\sum _{x|d}\mu (\frac{d}{x})F(d)$$

然后枚举倍数计算答案即可。

时间复杂度 $O(n\log n)$。

---

code

#include<cmath>
#include<cstdio>
#include<cstring>
#include<algorithm>
#define N 300005
#define P 1000000007
using namespace std;
int n,m,k,Max;
int num[N],cnt[N],fac[N],ifac[N],F[N];
int sum,prime[N],mark[N],mu[N];
int add(int x,int y)  {return x+y>=P?x+y-P:x+y;}
int dec(int x,int y)  {return x-y< 0?x-y+P:x-y;}
int mul(int x,int y)  {return 1ll*x*y%P;}
int power(int a,int b,int ans=1){
	for(;b;b>>=1,a=mul(a,a))
		if(b&1)  ans=mul(ans,a);
	return ans;
}
void linear_sieves(){
	mu[1]=1;
	for(int i=2;i<N;++i){
		if(!mark[i])  prime[++sum]=i,mu[i]=-1;
		for(int j=1;j<=sum&&i*prime[j]<N;++j){
			mark[i*prime[j]]=1;
			if(i%prime[j]==0)  {mu[i*prime[j]]=0;break;}
			mu[i*prime[j]]=-mu[i];
		}
	}
}
void prework(){
	fac[0]=fac[1]=1;
	for(int i=2;i<N;++i)  fac[i]=mul(fac[i-1],i);
	ifac[N-1]=power(fac[N-1],P-2);
	for(int i=N-2;~i;--i)  ifac[i]=mul(ifac[i+1],i+1);
}
int C(int n,int m){
	if(n<m)  return 0;
	return mul(fac[n],mul(ifac[m],ifac[n-m]));
}
int calc(int d){
	if(n-k>cnt[d])  return 0;
	return mul(C(cnt[d],n-k),mul(power(m/d-1,cnt[d]-(n-k)),power(m/d,n-cnt[d])));
}
int main(){
	scanf("%d%d%d",&n,&m,&k);
	linear_sieves(),prework();
	for(int i=1,x;i<=n;++i){
		scanf("%d",&x);
		num[x]++,Max=max(Max,x);
	}
	for(int i=1;i<=Max;++i)
		for(int j=i;j<=Max;j+=i)
			cnt[i]+=num[j];
	for(int i=1;i<=Max;++i)  F[i]=calc(i);
	for(int i=1;i<=m;++i){
		int ans=0;
		for(int j=1;i*j<=m;++j){
			if(!mu[j])  continue;
			ans=(mu[j]==1)?add(ans,F[i*j]):dec(ans,F[i*j]);
		}
		printf("%d ",ans);
	}
	return 0;
}