【数论】hdu5072 Coprime （容斥求互质+组合计数）

题目链接：

http://acm.split.hdu.edu.cn/showproblem.php?pid=5072

题目大意：

给N个数，要求选出其中3个，使他们都互质，或都不互质，问有多少种选法。

解题思路：

3个数互质或不互质，即为0组或3组互质，首先可以题目转化为求3个数中1组互质或2组互质的情况。（如下图）

问题转化为，3个数中，至少1对互质且至少1对不互质的情况。

这个很好求，设A=所有与a[i]互质的数，B=所有与a[i]不互质的数

对某个数a[i]，ans=A*B；且B=N-A-1；

由于N=10^5;需要用容斥原理来求。

容斥原理求互质数的原理如下：

假设n=12,m=30.

第一步：求出m的质因子：2,3,5；

第二步：(1,n)中是n的因子的倍数当然就不互质了(2,4,6,8,10)->n/2  6个,(3,6,9,12)->n/3  4个,(5,10)->n/5  2个。

如果是粗心的同学就把它们全部加起来就是：6+4+2=12个了，那你就大错特错了,里面明显出现了重复的,我们现在要处理的就是如何去掉那些重复的了！

第三步：这里就需要用到容斥原理了，公式就是：n/2+n/3+n/5-n/(2*3)-n/(2*5)-n/(3*5)+n/(2*3*5).

(奇数个因子相加，偶数个因子相减)

由于题目不是1-n的连续区间，稍微转化一下即可。

代码如下：

#include <cstdio>
#include <cstring>
#include <string.h>
#include <iostream>
#include <algorithm>
#include <cmath>
using namespace std;
#define N 100010
#define M 2000
typedef long long LL;


bool isprime[M];
int prime[M];

LL tot=0,n;
int factor[50];

bool mark[N];
int a[N];
int s[N];

void getprime() {
    for(int i=2; i<=M; i++) {
        if(isprime[i]) continue;
        prime[tot++]=i;
        for(int j=i*i; j<=M; j+=i)
            isprime[i]=true;
    }
}

LL solve() {
    LL ans=0;
    for(int i=1; i<=n; i++) {
        int m=a[i],num=0;
        LL sum=0;
        //求a[i]的质因子
        for(int j=0; j<tot&&prime[j]*prime[j]<=m; j++) {
            if(m%prime[j]==0) {
                factor[num++]=prime[j];
                while(m%prime[j]==0)
                    m/=prime[j];
            }
        }
        if(m!=1) factor[num++]=m;

//容斥法求互质
        for(int j=1; j<(1<<num); j++) {
            int op=0,tmp=1;
            for(int k=0; k<num; k++) {
                if((1<<k)&j) {
                    tmp*=factor[k];
                    op++;
                }
            }
            if(op&1) sum+=s[tmp];
            else sum-=s[tmp];
        }

        if(sum==0) continue;
        ans+=(sum-1)*(n-sum);  //sum不互质的数包括了a[i]自身故减1
    }
    return ans/2;

}
int main() {
    // freopen("in.txt","r",stdin);
    int t;
    scanf("%d",&t);
    getprime();
    while(t--) {
        scanf("%d",&n);
        memset(mark,0,sizeof(mark));
        memset(s,0,sizeof(s));
        for(int i=1; i<=n; i++) {
            scanf("%d",&a[i]);
            mark[a[i]]=true;
        }

        //若a[i]为质因数i的倍数，则s[i]++;(a[i]中能整除i的个数）
        //不用质数表遍历，因为要考虑质数倍数的情况
        for(int i=2; i<N; i++)
            for(int j=i; j<N; j+=i)
                if(mark[j]) s[i]++;

        LL ans=n*(n-1)*(n-2)/6;
        ans-=solve();
        cout<<ans<<endl;
    }
    return 0;
}