poj3904 Sky Code（容斥原理）

关键词：容斥原理、莫比乌斯函数
题意：n个数中gcd为1的4数对(a1,a2,a3,a4)的组数
解法：$A_{p}$：最大公约数是p的倍数的4数对组数
ans=$|∩(!A_{pi})|$，pi是n个数中的所有质因子。意思是gcd不被pi中的任意一个整除，即gcd==1
该式可由莫比乌斯函数化简

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>
#include <vector>
#include <algorithm>
#include <queue>
#define mem(a , b) memset(a , b , sizeof(a))
#define pb(a) push_back(a)
#define mp(a,b) make_pair(a,b)
using namespace std;
typedef long long ll;
const int maxn = 10000+10;

bool not_prime[maxn];
vector<ll> prime;
ll mu[maxn];

void make_prime(){
    mu[1]=1;
    for(ll i=2;i<maxn-5;i++){
        if(!not_prime[i]) { prime.push_back(i); mu[i]=-1; }
        for(ll j=0;j<prime.size()&&i*prime[j]<maxn-5;j++){
            not_prime[i*prime[j]]=1; mu[i*prime[j]]=-mu[i];
            if(!(i%prime[j])) { mu[i*prime[j]]=0; break; }
        }
    }
}

ll n,a[maxn];
ll beinum[maxn];
ll maxx,ans;
ll c4[maxn];

void add(ll x){
    for(ll i=1;i*i<=x;i++){
        if(x%i==0){
            beinum[i]++;
            if(i!=(x/i)) beinum[x/i]++;
        }
    }
}

void init(){
    c4[4]=1;
    for(ll i=5;i<maxn-5;i++){
        c4[i]=(c4[i-1]*i)/(i-4);
    }
}

int main(){
    make_prime();
    init();
    //freopen("a.txt","r",stdin);
    while(scanf("%lld",&n)!=EOF){
        memset(beinum,0,sizeof(beinum));
        maxx=0;
        for(ll i=1;i<=n;i++){
            scanf("%lld",&a[i]);
            add(a[i]);
            maxx=max(maxx,a[i]);
        }
        ans=0;
        for(ll i=1;i<=maxx;i++) ans+=(mu[i]*c4[beinum[i]]);
        printf("%lld\n",ans);
    }
    return 0;
}