spoj 4191. Sky Code(容斥)

http://www.spoj.com/problems/MSKYCODE/

给出一个集合，含有n个元素，每次任意从中取出4个使得他们的gcd是1，问有多少种取法。

可以先考虑问题的反面，就是取出的4个数的gcd为d，d > 1的方案数。总的方案数c(n,4)，减去d > 1的就是问题所求。

对于任意一个d(d > 1)，总能分解成若干素数的乘积。那么就要考虑到把d和素数联系在一起，那就是容斥原理啦.

以下d的范围：d > 1

先定义一个函数f(d)，表示能被d整除的四元组的个数。那么我们得到原题的逆问题的解

deg(d)表示d的素因子的个数。

注意这里的d不是任意一个大于1的整数，还要满足它的质因子没有重复这一条件。例如当d=3时，我们求出了所有能被3整除的四元组，再当d=9时，我们求出了能被9整除的四元组，若把9这个四元组加上的话就会加了重复的，因为在计算d = 3的时候已经计算过了。

#include <stdio.h>
#include <iostream>
#include <map>
#include <set>
#include <bitset>
#include <list>
#include <stack>
#include <vector>
#include <math.h>
#include <string.h>
#include <queue>
#include <string>
#include <stdlib.h>
#include <algorithm>
//#define LL __int64
//#define LL long long
#define eps 1e-9
#define PI acos(-1.0)
typedef long long LL;
using namespace std;

const int maxn = 10010;

int p_num[maxn];
bool flag[maxn];
int sef[maxn];

LL num[maxn];

LL cal(LL n)
{
	long long tmp = n;
	return tmp*(tmp-1LL)*(tmp-2LL)*(tmp-3LL)/(24LL);
}

void init()
{
	memset(flag,false,sizeof(flag));
	memset(p_num,0,sizeof(p_num));

	for(int i = 1; i < maxn; i++)
		sef[i] = 1;

	for(int i = 2; i < maxn; i++)
	{
		if(flag[i] == false)
		{
			for(int j = i; j < maxn; j += i)
			{
				sef[j] *= i;
				p_num[j]++;
				if(j != i)
					flag[j] = true;
			}
		}
	}
}

int main()
{
	LL n;
	LL a[maxn];
	LL num[maxn];
	init();
	while(~scanf("%lld",&n))
	{
		memset(num,0,sizeof(num));
		LL Max = 0;
		for(int i = 0; i < n; i++)
		{
			scanf("%lld",&a[i]);
			num[a[i]] ++;
			Max = max(Max,a[i]);
		}
		if(n < 4)
		{
			printf("0\n");
			continue;
		}
		for(int i = 2; i <= Max; i++)
		{
			for(int j = i+i; j <= Max; j += i)
				num[i] += num[j];
		}

		LL ans = 0;

		for(LL i = 2; i <= Max; i++)
		{
			if(sef[i] == i)
			{
				if(p_num[i] & 1)
					ans += cal(num[i]);
				else
					ans -= cal(num[i]);
			}
		}
		ans = cal(n) - ans;
		printf("%lld\n",ans);
	}
	return 0;
}