Cod sursa(job #619781)

Utilizator claudiumihailClaudiu Mihail claudiumihail Data 16 octombrie 2011 04:40:02
Problema Principiul includerii si excluderii Scor 100
Compilator cpp Status done
Runda Arhiva educationala Marime 2.48 kb
#include <fstream>
#include <iostream>
#include <bitset>
#include <vector>
#include <cmath>

#define MAX_NUMS 1001001

using namespace std;

typedef unsigned long uint32;
typedef long long int64;
typedef unsigned long long uint64;

template<unsigned long numbits>
struct PrimeSieve
{
	uint64 num;
	bitset<numbits> bits;

	void GeneratePrimes(const uint64 n)
	{
		num=1;
		for (uint32 i=1; ((i<<1)+1)*((i<<1)+1)<n; ++i)
		{
			if (!bits[i])
			{
				for (uint64 j=i; ((i<<1)+1)*((j<<1)+1)<n; ++j)
				{
					const long long index=((i<<1)+1)*((j<<1)+1);
					//cout<<"("<<((i<<1)+1)*((j<<1)+1)<<" "<<(i+j)+1<<")"<<" ";
					bits[index>>1]=1;
				}
			}
		}
	}
	
	const bitset<numbits>& GetPrimes() const
	{
		return bits;
	}
};

uint32 divisors[MAX_NUMS];
PrimeSieve<16000002> primeSieve;

uint64 InclusionExclusion
(
	const uint64 A,
	const uint32 *divisors,
	const int numDivisors
)
{
	int idx=1;
	int64 sum=0;
	
	while (idx < (1<<numDivisors))
	{
		int aux = idx;
		int nr = -1;
		int64 prod = 1;

		while (aux)
		{
			const int index = __builtin_ffs(aux);
			//cout << divisors[index-1] << " ";
			prod *= divisors[index-1];
			nr++;
			aux = aux & (~(1<<(index-1)));
		}
		//cout << "\n";
		idx++;
		
		if (nr%2)
			nr = -1;
		else
			nr = 1;
		
		sum = sum + 1LL * nr * (uint64)(A/prod);
	}
	
	return sum;
}

int main()
{
	int m;
	uint64 a=50,b=30;
	vector<uint32> primes;
	fstream fin("pinex.in", fstream::in);
	fstream fout("pinex.out", fstream::out);

	primeSieve.GeneratePrimes(MAX_NUMS);
	//cout << endl;

	primes.push_back(2);
	for(int i=1; (i<<1)+1<MAX_NUMS; ++i)
	{
		if(! (primeSieve.GetPrimes())[i] )
		{
			//cout<<(i<<1)+1<<" ";
			primes.push_back((i<<1)+1);
		}
	}
	//cout << endl;
	
	fin >> m;
	//cout << m << endl;
	for (int i=0; i<m; ++i)
	{
		fin >> a >> b;
		//cout << a << " " << b << endl;

		uint32 k=0, numDivs=0;
		uint32 sq = (int)ceil(sqrt(b));
		
		//cout << sq << endl;
		
		while (k<primes.size() && b>1)
		{
			if ((b%primes[k]) == 0)
			{
				divisors[numDivs] = primes[k];
				numDivs++;
				
				while ((b%primes[k]) == 0)
					b /= primes[k];
			}
			k++;
			
			if (primes[k] > sq && b>1)
			{
				divisors[numDivs] = b;
				numDivs++;
				break;
			}
		}
		
		//for (int j=0; j<numDivs; ++j)
		//	cout << divisors[j] << " ";
		//cout << endl;
		fout << a - InclusionExclusion(a, divisors, numDivs) << "\n";
	}
	
	fin.close();
	fout.close();
	return 0;
}