hdu4675 GCD of Sequence
题目链接:http://acm.hdu.edu.cn/showproblem.php?pid=4675
题意:
给定一个长度为n的序列a,且 1<=a[i]<=m,求分别有多少个序列b,使得GCD(b[1],b[2],...b[n])=x (1<=x<=m),且正好有k个b[i]!=a[i]。
分析:
莫比乌斯反演,主要是确定F(x)。
用F(x)表示gcd为x的倍数的方案数,f(x)表示gcd为x的方案数。
先考虑F(d)怎么计算。可以把a数组中的数分成两类,第一类是必须对应下标不等的,即a[i]不是d的倍数),其他的就是第二类。
假设第二类的数量是p,第一类的数量就是n−p,因为要选择k个不同的,第一类必须不同,所以需要在第二类中选择k−n+p个,而b数组中每一个数都有[m/p]种选择。
所以最终的结果就是F(d)=C(p,k-n+p) ∗ ( ([m/d]−1)^(k−n+p) )∗( [m/d]^(n−p) ),然后暴力计算每一个f(d)就好了。总的复杂度是n(logn)。
代码:
1 #include<iostream> 2 #include<cstdio> 3 #include<cmath> 4 #include<cstring> 5 6 using namespace std; 7 const int maxn=300010; 8 const long long mod=1e9+7; 9 10 int n,m,k; 11 int a[maxn]; 12 int mu[maxn]; 13 int vis[maxn]; 14 int prime[maxn]; 15 int cnt; 16 int num[maxn]; 17 long long jie[maxn],ni[maxn]; 18 long long F[maxn]; 19 long long res[maxn]; 20 21 void init() 22 { 23 memset(vis,0,sizeof(vis)); 24 mu[1]=1; 25 cnt=0; 26 for(int i=2;i<maxn;i++) 27 { 28 if(!vis[i]) 29 { 30 prime[cnt++]=i; 31 mu[i]=-1; 32 } 33 for(int j=0;j<cnt&&i*prime[j]<maxn;j++) 34 { 35 vis[i*prime[j]]=1; 36 if(i%prime[j]) 37 mu[i*prime[j]]=-mu[i]; 38 else 39 { 40 mu[i*prime[j]]=0; 41 break; 42 } 43 } 44 } 45 } 46 47 long long power(long long a,long long n,long long m) 48 { 49 long long ans=1,tmp=a%m; 50 while(n) 51 { 52 if(n&1) 53 ans=ans*tmp%m; 54 tmp=tmp*tmp%m; 55 n=n/2; 56 } 57 return ans; 58 } 59 60 long long C(long long n,long long m) 61 { 62 if(n==0) 63 return 1; 64 return jie[n]*ni[m]%mod*ni[n-m]%mod; 65 } 66 67 int main() 68 { 69 init(); 70 ni[0]=1; 71 jie[0]=1; 72 for(int i=1;i<maxn;i++) 73 { 74 jie[i]=jie[i-1]*i%mod; 75 ni[i]=power(jie[i],mod-2,mod); 76 } 77 while(~scanf("%d%d%d",&n,&m,&k)) 78 { 79 for(int i=1;i<=n;i++) 80 scanf("%d",&a[i]); 81 memset(num,0,sizeof(num)); 82 for(int i=1;i<=n;i++) 83 num[a[i]]++; 84 for(int i=1;i<=m;i++) 85 { 86 long long p=0; 87 for(int j=1;j*i<=m;j++) 88 p+=num[i*j]; 89 if(k-n+p<0) 90 F[i]=0; 91 else 92 F[i]=C(p,k-n+p)*power(m/i-1,k-n+p,mod)%mod*power(m/i,n-p,mod)%mod; 93 } 94 for(int i=1;i<=m;i++) 95 { 96 long long ans=0; 97 for(int j=1;i*j<=m;j++) 98 { 99 ans+=mu[j]*F[i*j]; 100 ans=(ans%mod+mod)%mod; 101 } 102 if(i==1) 103 printf("%lld",ans); 104 else 105 printf(" %lld",ans); 106 } 107 printf("\n"); 108 } 109 return 0; 110 }