HDU5072 容斥原理

G - Coprime

Time Limit:1000MS     Memory Limit:262144KB     64bit IO Format:%I64d & %I64u

Submit Status Practice HDU 5072

Description

There are n people standing in a line. Each of them has a unique id number. 

Now the Ragnarok is coming. We should choose 3 people to defend the evil. As a group, the 3 people should be able to communicate. They are able to communicate if and only if their id numbers are pairwise coprime or pairwise not coprime. In other words, if their id numbers are a, b, c, then they can communicate if and only if [(a, b) = (b, c) = (a, c) = 1] or [(a, b) ≠ 1 and (a, c) ≠ 1 and (b, c) ≠ 1], where (x, y) denotes the greatest common divisor of x and y. 

We want to know how many 3-people-groups can be chosen from the n people.

 

Input

The first line contains an integer T (T ≤ 5), denoting the number of the test cases. 

For each test case, the first line contains an integer n(3 ≤ n ≤ 10 ⁵), denoting the number of people. The next line contains n distinct integers a ₁, a ₂, . . . , a _n(1 ≤ a _i ≤ 10 ⁵) separated by a single space, where a _i stands for the id number of the i-th person.

 

Output

For each test case, output the answer in a line.

 

Sample Input

1 5 1 3 9 10 2

 

Sample Output

4

 

题意：给出n(3 ≤ n ≤ 10⁵)个数字，每个数ai满足1 ≤ a_i ≤ 10⁵，求有多少对(a,b,c)

　　　满足[(a, b) = (b, c) = (a, c) = 1] or [(a, b) ≠ 1 and (a, c) ≠ 1 and (b, c) ≠ 1]，

　　　都互素或都不互素。

思路：由于数据范围不大10^5以内，总组合数C(n,3) longlong不会爆。

　　　abc两两互质和两两不互质，就对应着两个互质另两个不互质，这两个集合构成了全集U。

　　　不妨把前者称为集合A,后者称为集合B，那么A并B等于U，且A交B为空。U的大小为C(n,3)。

 　　  如果a，b，c不符合条件，必然有一对互质，一对不互质，不妨设a，b互质，b，c不互质，

　　　于是我们可以枚举b来统计所有的三元组：如果a，c互质那么这样的三元组中b，c可以互换位置；

　　　如果a，c不互质，那么a，b可以互换位置。每个答案被算了两遍。

　　   所以只要枚举每个b，统计出k个和它不互质的，那么剩下n-1-k个就是和它互质的，

　　　那么三元组就有k*(n-1-k)/2种。

　　　对于b不超过10^5，质因子的个数不超过6个（2*3*5*7*11*13 *17>10^5）。

　　　用状压搜索质因子组成的每个因数，如果某数是该因数的倍数，

　　　那么就说明该数和b是不互质的。利用容斥原理统计出与b不互质的数的综述。

　　　由于数据范围不超过10^5，预处理筛除出每个质数和每个质因子，复杂度为nlogn。

　　　对于具体的n个数，再筛出在n个数中以他们为倍数的数的个数也是nlogn。(代码中用cntExtend[]记录)

#include <cstdio>
#include <cstring>
#include <iostream>
#include <cstring>
using namespace std;
typedef long long LL;
const int maxn = 111111;
int prim[maxn],isprim[maxn];
int primnum=0;
void initprim(){
    memset(isprim,-1,sizeof isprim);
    isprim[0]=isprim[1]=0;
    prim[primnum++] = 2;
    for(int i=4;i<maxn;i+=2){
        isprim[i]=0;
    }
    for(int i=3;i<maxn;i+=2){
        if(isprim[i]){
            prim[primnum++]=i;
            for(int j=i+i;j<maxn;j+=i){
                isprim[j]=0;
            }
        }
    }
}
int has[maxn];
int factor[maxn][6];
int factornum[maxn];
void getfactor(){
    for(int num=2;num<maxn;num++){
        int n=num,cnt = 0;
        for(int i=0;i<primnum;i++){
            if(isprim[n]) {
                factor[num][cnt++]=n;
                break;
            }
            if(n%prim[i]==0){
                factor[num][cnt++]=prim[i];
                while(n%prim[i]==0){
                    n/=prim[i];
                }
            }
        }
        factornum[num]=cnt;
    }
}
int num[maxn],cntExtend[maxn];
void factorExtend(int len){
    memset(cntExtend,0,sizeof cntExtend);
    for(int i=1; i<maxn; i++){
        for(int j=i; j<maxn; j+=i){
            if(has[j])
                cntExtend[i]++;
        }
    }
}
LL solve(int len){
    LL re = 0;
    for(int i=0; i<len; i++){
        int n = num[i];
        if(n==1) continue;
        int facnum = factornum[n];
        LL sum=0;
        for(int k=(1<<facnum)-1; k>0; k--){
            int mul=1,b=0;
            for(int j=0; j<facnum; j++){
                if((1<<j) & k) {
                    mul*=factor[n][j];
                    b^=1;
                }
            }
            if(b){
                sum+=cntExtend[mul]-1;
            }else{
                sum-=cntExtend[mul]-1;
            }
        }
        re+=(sum)*(len-1-sum);
    }
    return re;
}

int main(){
    int T,n,x;
    initprim();
    getfactor();
    scanf("%d",&T);
    while(T--){
        scanf("%d",&n);
        memset(has,0,sizeof has);
        for(int i=0;i<n;i++){
            scanf("%d",&x);
            has[x]++;
            num[i]=x;
        }
        factorExtend(n);
        LL ans = (LL)n*(n-1)*(n-2)/6 - solve(n)/2;
        printf("%I64d
",ans);
    }

}