poj3904 容斥原理或莫比乌斯反演

传送门

题意：从给出的n个数中选取4个最大公因数为1的四元组的选取数

思路：感觉容斥问题都有一个很相似的地方，就是要从问题的逆方向来考虑，所谓正难则反，这也是思维中一个很重要的地方。我们来考虑这个逆问题，便是选取4个最大公因数不为1的四元组，总选取方式为C(n,4)，那么也就是这四个数的最大公因数里至少有一个素数，换句话说，至少有一个相同的素数能够将这四个数整除，如果这n个数中有m个数是这个素数的倍数，那么C(m,4)便是这个素数下互相不互素的选取方式数，举个例子：

给出10个数：12，14，15，16，17，18，20，21，24，30我们从素数2开始看，则2的倍数有7个，选取的组合数为C(7,4)，意味着从这五个数中选取4个数一定是不互素的，我们再看下一个素数3，则3的倍数有6个，选取的组合数为C(6,4)，这些组合数中有重复的部分，6的倍数有4个，所以还需要减去C(4,4)，这里就很明显看出来是一个容斥原理了

所以问题的关键是求素数组合的倍数的个数（好绕啊！！！）。由于数据范围不超1e4，所以对每个数进行素因子分解，然后再求每个数下的素数组合，每求得一个素数组合（这个素数组合不会大于1e4）便在当前记录下加一，具体看代码吧……

最后不要忘记用C(n,4)减去得到的结果……

有错误请指正……

完整代码：

#include <cstdio>
#include <cstring>
#include <iostream>
#include <vector>
using namespace std;
typedef long long LL;
const int N=10005;
int data[N];
int p[26]={2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79,83,89,97};
LL num[N];
LL Count[N];
LL cn;
LL cnt;
vector<int> temp;
int main()
{
    LL n;
    while(scanf("%lld",&n)!=-1)
    {
        memset(Count,0,sizeof(Count));
        for(LL i=0;i<n;i++)
            scanf("%d",&data[i]);
        for(LL i=0;i<n;i++)
        {
            temp.clear();
            cnt=0;
            int a=data[i];
            for(int j=0;j<25&&p[j]*p[j]<=a;j++)
            {
                if(a%p[j]==0)
                {
                    temp.push_back(p[j]);
                    while(a%p[j]==0)
                        a/=p[j];
                }
            }
            if(a>1)
                temp.push_back(a);
            for(LL j=1;j<(1<<temp.size());j++)
            {
                LL mul=1;
                LL cnt=0;
                for(LL k=0;k<temp.size();k++)
                {
                    if(j&(1<<k))
                    {
                        mul=mul*temp[k];
                        if(mul>N) break;
                        cnt++;
                    }
                }
                if(mul>N) continue;
                Count[mul]++;
                num[mul]=cnt;
            }
        }
        LL S=0;
        for(LL i=0;i<N;i++)
        {
            if(Count[i]>=4)
            {
                if(num[i]&1)
                    S=S+Count[i]*(Count[i]-1)*(Count[i]-2)*(Count[i]-3)/24;
                else
                    S=S-Count[i]*(Count[i]-1)*(Count[i]-2)*(Count[i]-3)/24;
            }
        }
        LL ans=n*(n-1)*(n-2)*(n-3)/24;
        ans=ans-S;
        printf("%lld\n",ans);
    }
    return 0;
}

这个题也可以用莫比乌斯反演的思想做，其实最后得到的结果是一致的

我们设f(k)=四个数最大公因数为k的方法组合数

       设F(k)=四个数的最大公因数为k的倍数的方法组合数

那么我们要求的f(1)=∑mu(d)F(d)=mu(1)F(1)+mu(2)F(2)+...+mu(max)F(max)

那么关键就是F(k)怎么求了，仔细想F(k)为最大公因数是k的倍数的方法数，那就从是k的倍数的数中任意取4个就可以了，那就跟上面的思想一样了，对每个数素因子分解，记录每个素因子组合个数就可以了

莫比乌斯反演代码：

#include<cstdio>
#include<iostream>
#include<cstring>
#include<vector>
using namespace std;
typedef long long LL;
const int N=10005;
int mu[N];
bool check[N+5];
int prime[N+5];
LL data[N];
LL num[N];
void miu()
{
    memset(check,false,sizeof(check));
    mu[1]=1;
    int tot=0;
    for(int i=2;i<=N;i++)
    {
        if(!check[i])
        {
            prime[tot++]=i;
            mu[i]=-1;
        }
        for(int j=0;j<tot;j++)
        {
            if(i*prime[j]>N) break;
            check[i*prime[j]]=true;
            if(i%prime[j]==0)
            {
                mu[i*prime[j]]=0;
                break;
            }
            else
                mu[i*prime[j]]=-mu[i];
        }
    }
}
LL c(LL n)
{
    return n*(n-1)*(n-2)*(n-3)/24;
}
vector<LL> p;
int main()
{
    LL n;
    miu();
    while(scanf("%lld",&n)!=-1)
    {
        LL cnt=0;
        memset(num,0,sizeof(num));
        for(LL i=0;i<n;i++)
        {
            scanf("%lld",&data[i]);
            p.clear();
            LL t=data[i];
            for(LL j=2;j*j<=t;j++)
            {
                if(t%j==0)
                {
                    p.push_back(j);
                    while(t%j==0)
                        t/=j;
                }
            }
            if(t>1)
                p.push_back(t);
            for(LL j=1;j<(1<<p.size());j++)
            {
                LL mul=1;
                for(LL k=0;k<p.size();k++)
                {
                    if(j&(1<<k))
                    {
                        mul=mul*p[k];
                        if(mul>N) break;
                    }
                }
                if(mul>N) continue;
                num[mul]++;
            }
        }
        LL ans=c(n);
        for(LL i=1;i<N;i++)
        {
            if(num[i])
                ans+=c(num[i])*mu[i];
        }
        printf("%lld\n",ans);
    }
    return 0;
}