欢迎您访问程序员文章站本站旨在为大家提供分享程序员计算机编程知识!
您现在的位置是: 首页

2018.06.29 NOIP模拟 Gcd(容斥原理)

程序员文章站 2022-05-08 19:01:19
...

Gcd
题目背景
SOURCE:NOIP2015-SHY-2
题目描述
给出n个正整数,放入数组 a 里。
问有多少组方案,使得我从 n 个数里取出一个子集,这个子集的 gcd 不为 1 ,然后我再从剩下的数中取出一个数,把他放进刚刚取出的子集里,使得 gcd 为 1 。
输出方案数 mod (10^9 + 7)。
输入格式
第一行一个数 n 。
第二行 n 个数,表示 a 数组。
输出格式
输出一个数表示答案。
样例数据 1
输入
3
2 3 2
输出
5
备注
【样例说明】
set :2 number:3 count:2
set :3 number:2 count:2
set :2 2 number:3 count:1
total = 2 + 2 + 1 = 5
【数据范围】
对 30% 的输入数据 :1≤n≤10;
对 100% 的输入数据 :1≤n≤500000;2≤a中元素≤10000000 。

这是我做过的最可(dudu)做(liuliu)的NOIPNOIP模拟级别的数论题。考场上打算写随机算法,想了想3030暴力更稳,于是交了暴力,结果1010分滚粗了。

这题让我们联想到正难则反的思想,题目上要我们求出所有可能解的方案数,那么我们这样想,我们先把总方案数求出来,显然是n(2n11)n∗(2n−1−1)

然后我们将不合法的情况去掉,我们设加入的数为numbernumber,那么我们要去掉的就是在加入numbernumber之前就已经使得gcdgcd11的集合对应的numbernumber的选择方案数之和,以及在加入numbernumber之后gcdgcd仍然不为11的方案数。这东西怎么算啊?变形代数式看看吧。

首先解决第二种情况,即在加入numbernumber之后gcdgcd仍然不为11的方案数。这种情况直接统计所有集合中使得gcdgcd不为一的方案数即可,即用总方案数减去所有集合中使得gcdgcd11的方案数。
另外,我们还知道,情况二其实是让我们求这个东西:|S|,gcd(S)!=1|S|∑|S|,gcd(S)!=1|S|

这样的话,我们似乎将第二种情况转化成了第一种情况的子问题。怎么求出所有集合中使得gcdgcd11的方案数呢?

考虑容斥原理,我们先将gcdgcd11的倍数的集合数求出,减去gcdgcd22的倍数的集合数……一直这样迭代下去,然后有个显然的结论,每个gcdgcd的容斥系数就对应着它们的莫比乌斯函数值。于是只需要提前用线性筛筛出所有数的莫比乌斯函数值,每次将情况二的答案加上mu[i](2sum[i]1)mu[i]∗(2sum[i]−1)就行了,其中sum[i]sum[i]表示所有数列中所有是ii的倍数的数的个数。

情况一?考虑用情况二反推情况一,再次将式子变形:我们要求的是在加入numbernumber之前就已经使得gcdgcd11的集合对应的numbernumber的选择方案数之和,即S,gcd(S)==1(n|S|)∑S,gcd(S)==1(n−|S|),显然这个式子可以展开成S,(gcd)==11S,gcd(S)==1|S|∑S,(gcd)==11−∑S,gcd(S)==1|S|,哇,后面一坨不是可以跟情况二合并吗?好像还可以抵掉一坨东西。最后发现只需要容斥计算乱搞一下就行了啊。

代码如下:

#include<bits/stdc++.h>
#define inf 10000005
#define mod 1000000007
#define ll long long
using namespace std;
int n,tot=0,mu[inf],prime[inf],cnt[inf],sum[inf];
ll mul[inf],a=0,b=0;
bool pr[inf];
inline ll read(){
    ll ans=0;
    char ch=getchar();
    while(!isdigit(ch))ch=getchar();
    while(isdigit(ch))ans=(ans<<3)+(ans<<1)+ch-'0',ch=getchar();
    return ans;
}
int main(){
    n=read();
    for(int i=1;i<=n;++i){
        int x=read();
        ++cnt[x];
    }
    mul[0]=1;
    for(int i=1;i<=n;++i)mul[i]=mul[i-1]<<1,mul[i]%=mod;
    mu[1]=1;
    for(int i=2;i<=inf;++i){
        if(!pr[i])prime[++tot]=i,mu[i]=-1;
        for(int j=1;j<=tot&&i*prime[j]<=inf;++j){
            int k=i*prime[j];
            pr[k]=1;
            if(i%prime[j]==0){
                mu[k]=0;
                break;
            }
            mu[k]=-mu[i];
        }
    }
    for(int i=1;i<=inf;++i)
        for(int j=i;j<=inf;j+=i)sum[i]+=cnt[j];
    for(int i=1;i<=inf;++i)
        if(mu[i]&&sum[i]){
            a=(a+mu[i]*mul[sum[i]]-mu[i])%mod;
            b=(b+mu[i]*mul[sum[i]-1]*sum[i])%mod;
        }
    a=a*n%mod;
    b=(b<<1)%mod;
    a=(b-a+mod)%mod;
    printf("%lld",a);
    return 0;
}