POJ_3904_Sky_Code 容斥定理

#include <stdio.h>
#include <string.h>
#include <algorithm>
using namespace std;
long long p[10005];
long long prime[10005];
long long num[10005];
long long cout[10005];
void init()
{
	long long i;
	memset(num,0,sizeof(num));
	memset(p,0,sizeof(p));
	for(i=4;i<=10005;i++)
	{
		p[i] = i * (i-1) * (i-2) * (i-3) / 24;
	}
}
void find(long long x)
{
	long long i,j;
	long long top = 0;
	for(i=2;i*i<=x;i++)
	{
		if(x%i==0)
		{
			prime[top++] = i;
			while(x%i==0)
			{
				x /= i;
			}
		}
	}
	if(x>1)
	{
		prime[top++] = x;
	}
	long long t, sum;
	for(i = 1 ;i < (1<<top); i++)
	{
		t = 1;
		sum = 0;
		for(j=0;j<top;j++)
		{
			if(i & (1<<j))
			{
				t *= prime[j];
				sum ++ ;
			}
		}
		cout[t] ++ ;   //目前因子的個數
		num[t] = sum;  //目前因子是由多少個素因子組成 
	}	
}
int main()
{
	init();
	long long i,j,n,w;
	while(scanf("%lld",&n)!=EOF)
	{
		long long M = - 0x3f3f3f3f;
		memset(cout,0,sizeof(cout));
		for(i=0;i<n;i++)
		{
			scanf("%lld",&w);
			M = max(w,M);
			find(w);
		}
		long long res  = 0;  
		for(i=0;i<M;i++)
		{
			if(cout[i])
			{
				if(num[i]%2)
				{
					res += p[cout[i]];  //容斥定理 奇數加偶數減 
				}
				else
				{
					res -= p[cout[i]];
				}
			}
		}
		printf("%lld\n",p[n]-res);
	} 
	return 0;
}