#include <iostream>
#include <complex>
//#include <limits>
#include <vector>
#include <deque>
#include <list>
#include <algorithm>

#include <ctime>

#include <pthread.h>
//#include <machine/cpufunc.h>
#include <sched.h>

typedef u_int64_t integer_t;
//typedef unsigned __int32 integer_t;
typedef std::vector<integer_t> num_list_t;
typedef std::deque<integer_t> cand_list_t;
typedef std::vector<cand_list_t> cand_list_list_t;

namespace
{
	const integer_t MAX_TARGET_VALUE = 4294967296;

	//bool divided;

	typedef struct
	{
		cand_list_t* nl;
		integer_t ord;
		integer_t dv;
	} thread_arg_t;


	class prscrmsk
	{
		cpu_set_t i_msk_;
	public:
		prscrmsk(u_int msk)
		{
			sched_getaffinity(0, sizeof(cpu_set_t), &i_msk_);

			cpu_set_t s;
			CPU_ZERO(&s);
			CPU_SET(msk, &s);

			if(sched_setaffinity(0, sizeof(cpu_set_t), &s))
			{
				throw 0;
			}
		}
		~prscrmsk()
		{
			sched_setaffinity(0, sizeof(cpu_set_t), &i_msk_);
		}
	};

	__inline__ u_int64_t rdtsc()
	{
		u_int64_t x;
		// x86_64 だと、なぜか次の 1 行がうまくいかない
		//__asm__ volatile (".byte 0x0f, 0x31" : "=A" (x));

		u_int32_t a, d;
		__asm__ volatile ("rdtsc" : "=d" (d), "=a" (a));
		x = static_cast<u_int64_t>(d) << 32;
		x += a;

		return x;
	}

	u_int64_t tsc()
	{
		prscrmsk m(0);

		return rdtsc();
	}

	u_int64_t time()
	{
		std::time_t t;
		std::time(&t);

		return t;
	}

	u_int det_phy_cores()
	{
		cpu_set_t t;
		sched_getaffinity(0, sizeof(cpu_set_t), &t);
		return CPU_COUNT(&t);
	}

	void init(cand_list_list_t& cll, integer_t start, integer_t end, u_int32_t cores)
	{
		cand_list_t tmp;

		for(u_int32_t i = 0; i < cores; ++i)
		{
			cll.push_back(tmp);
		}

		integer_t unit = end / cores;

		int unum = 0;
		for(integer_t n = start; n < end; ++n)
		{
			cll[unum].push_back(n);

			if(!(n % unit))
			{
				++unum;
			}
		}
	}

	void print(integer_t i)
	{
		std::cout << i << " ";
	}

	template <typename Integer> class Divi
	{
		Integer N_;

	public:
		Divi(Integer N) : N_(N) {};
		bool operator()(Integer i) const
		{
			return !(i % N_);
		};
		void setn(Integer N) { N_ = N; };
	};

	void* thread(void *args)
	{
		thread_arg_t* tmp = static_cast<thread_arg_t*>(args);

		Divi<integer_t> d(tmp->dv);

		tmp->nl->erase(std::remove_if(tmp->nl->begin(), tmp->nl->end(), d), tmp->nl->end());

		return 0;
	}
}

int main()
{
	using std::cout;
	using std::endl;

	u_int64_t r = tsc();
	u_int64_t t = time();

	u_int cores = det_phy_cores();

	cout << "初期化中" << endl;

	num_list_t prime;
	prime.push_back(2);

	std::vector<cand_list_t> candid;
	init(candid, 3, MAX_TARGET_VALUE/(2 * 2 * 2 * 2 * 2 * 2), cores);

	std::vector<thread_arg_t> ta;
	for(unsigned int i = 0; i < cores; ++i)
	{
		thread_arg_t tmp;
		tmp.nl = &candid[i];
		tmp.ord = MAX_TARGET_VALUE / cores * i;
		tmp.dv = 2;
		ta.push_back(tmp);
	}

	cout << cores << " スレッドで計算をスタートします。" << endl;

	std::vector<pthread_t> ptv;
	// int の最大値まで計算する場合はループの終了条件に注意
	for(Divi<integer_t> d(2); prime.back() * prime.back() < candid[cores -1].back(); d.setn(prime.back()))
	{
		ptv.reserve(cores);
		for(u_int i = 0; i < cores; ++i)
		{
			pthread_t pt;
			pthread_create(&pt, 0, thread, static_cast<void*>(&ta[i]));

			ptv.push_back(pt);
		}

		for(u_int i = 0; i < cores; ++i)
		{
			pthread_join(ptv[i], NULL);
		}

		prime.push_back(candid[0].front());
		candid[0].pop_front();

		for(u_int i = 0; i < cores; ++i)
		{
			ta[i].dv = prime.back();
		}
		ptv.clear();
	}

	for(u_int32_t i = 0; i < cores; ++i)
	{
		prime.insert(prime.end(), candid[i].begin(), candid[i].end());
	}

	u_int64_t s = tsc();
	s -= r;
	cout << s << " クロックかそこら ";
	std::time_t u = time();
	u -= t;
	cout << u << " 秒で計算完了 " << endl;
	cout << s / u << " Hz くらいの CPU だと思われます。" << endl;
//	std::for_each(prime.begin(), prime.end(), print);
	cout << "計算した範囲で最大の素数は " << prime.back() << " でした。" << endl;

	return 0;
}