The Big Programming Thread - Page 14

#include <fstream>
#include <math.h>
#include "StopWatch.h"

using namespace std;

const long TOTAL_NUMS = 2000000;
const int INITIALPRIMES_SIZE = 1415; // INITIALPRIMES_SIZE = sqrt( TOTAL_NUMS )
const long SIEVECHUNK_SIZE = 500000; // arbitrary but anything much bigger crashes
long long primesTotal;
long primeCount;

void FillSievePrimes( int sievePrimes[], long totalNums );
void InitialSieve( bool initialPrimes[], int sieveNum, int sieveMax );
void InitializeTestPrimes( bool testPrimes[], int duplicateCount );
void SieveChunk( int chunkNumber, bool testPrimes[], int sievePrimes[] );
void NormalSieve( bool testPrimes[], long sieveNum, long sieveMin, long sieveMax );
void OutputResults( double elapsedTime );

void main()
{
	bool testPrimes[SIEVECHUNK_SIZE];
	int sievePrimes[INITIALPRIMES_SIZE / 2];
	StopWatch stopWatch;

	stopWatch.StartTimer();

	// fill sievePrimes with the primes up to sqrt( 2,000,000 ) to be used for the rest of the sieves
	FillSievePrimes( sievePrimes, TOTAL_NUMS );

	// prepare the testPrimes array for use
	InitializeTestPrimes( testPrimes, INITIALPRIMES_SIZE );

	// for each of 4 chunks of 500,000 numbers, use the sievePrimes to determine the primes in that chunk
	for( int chunkNumber = 0; chunkNumber < 4; chunkNumber++ )
		SieveChunk( chunkNumber, testPrimes, sievePrimes );

	stopWatch.StopTimer();

	OutputResults( stopWatch.GetElapsedTime() );
}

void FillSievePrimes( int sievePrimes[], long totalNums )
{
	bool initialPrimes[INITIALPRIMES_SIZE];
	int sieveNum = 2;
	int max = int( sqrt( sqrt( float( totalNums ) ) ) ) + 1; // max = 38

	primeCount = 0;

	for( int i = 0; i < INITIALPRIMES_SIZE; i++ )
		initialPrimes[i] = true;

	while( sieveNum < max )
	{
		sievePrimes[primeCount] = sieveNum;
		primeCount++;
		primesTotal += sieveNum;
		InitialSieve( initialPrimes, sieveNum, INITIALPRIMES_SIZE );
		sieveNum++;
		while( initialPrimes[sieveNum] == false )
			sieveNum++;
	}
	for( int i = sieveNum; i < INITIALPRIMES_SIZE; i++ )
		if( initialPrimes[i] == true )
		{
			sievePrimes[primeCount] = i;
			primeCount++;
			primesTotal += i;
		}
	sievePrimes[primeCount] = totalNums;
}

void InitialSieve( bool initialPrimes[], int sieveNum, int sieveMax )
{
	int currentNum = sieveNum * sieveNum;

	while( currentNum < sieveMax )
	{
		initialPrimes[currentNum] = false;
		currentNum += sieveNum;
	}
}

void InitializeTestPrimes( bool testPrimes[], int duplicateCount )
{
	// set all to true initially
	for( int i = 0; i < SIEVECHUNK_SIZE; i++ )
		testPrimes[i] = true;

	// prevent duplicate prime counting since FillSievePrimes already counted primes up to INITIALPRIMES_SIZE
	for( int i = 0; i < duplicateCount; i++ )
		testPrimes[i] = false;
}

void SieveChunk( int chunkNumber, bool testPrimes[], int sievePrimes[] )
{
	long sieveIndex = 0;
	long sieveNum = sievePrimes[sieveIndex];
	long sieveChunkMin = SIEVECHUNK_SIZE * chunkNumber;
	long sieveChunkMax = SIEVECHUNK_SIZE * ( chunkNumber + 1 );
	int max = int( sqrt( float( sieveChunkMax ) ) );

	testPrimes[0] = false;

	while( sieveNum < max )
	{
		NormalSieve( testPrimes, sieveNum, sieveChunkMin, sieveChunkMax );
		sieveIndex++;
		sieveNum = sievePrimes[sieveIndex];
	}

	for( long i = 1; i < SIEVECHUNK_SIZE; i+=2 )
	{
		if( testPrimes[i] == true )
		{
			primeCount++;
			primesTotal += i + sieveChunkMin;
		}
		else
			testPrimes[i] = true;
	}
}

void NormalSieve( bool testPrimes[], long sieveNum, long sieveMin, long sieveMax )
{
	long currentNum = sieveNum * sieveNum;

	if( currentNum < sieveMin )
		currentNum = sieveMin + sieveNum - ( sieveMin % sieveNum );

	while( currentNum < sieveMax )
	{
		if( currentNum > sieveMin )
			testPrimes[currentNum - sieveMin] = false;
		currentNum += sieveNum;
	}
}

void OutputResults( double elapsedTime )
{
	ofstream outputFileStream;

	outputFileStream.open( "output.txt", ios::out );
	outputFileStream << "Total of primes up to 2,000,000: " << primesTotal << endl << endl
									 << "Elapsed time: " << elapsedTime << endl;
	outputFileStream.close();
}