/*
 bitonic.c 

 This file contains two different implementations of the bitonic sort
        recursive  version :  recBitonicSort()
        imperative version :  impBitonicSort() 
 

 The bitonic sort is also known as Batcher Sort. 
 For a reference of the algorithm, see the article titled 
 Sorting networks and their applications by K. E. Batcher in 1968 


 The following codes take references to the codes avaiable at 

 http://www.cag.lcs.mit.edu/streamit/results/bitonic/code/c/bitonic.c

 http://www.tools-of-computing.com/tc/CS/Sorts/bitonic_sort.htm

 http://www.iti.fh-flensburg.de/lang/algorithmen/sortieren/bitonic/bitonicen.htm 
 */

/* 
------- ---------------------- 
   Nikos Pitsianis, Duke CS 
-----------------------------
*/

#include <stdio.h>
#include <stdlib.h>
#include <sys/time.h>
#include <time.h>
#include <pthread.h>

typedef enum { false, true } bool;
struct timeval startwtime, endwtime;
double seq_time;

int threadsNum; //number of threads
int N; // data array size
int *a; // data array to be sorted
unsigned randSeed; //seed array initialisation
bool sortPass; //flag showing whether the test passed or not

const int ASCENDING  = 1;
const int DESCENDING = 0;
//min lenght of an array to be sorted in parallel manner
const int REC_BITONIC_SORT_PARALLEL_MIN = (1 << 22) + 1;

//structs for parallel functions
typedef struct recBitonicSortData{
	int threadID;
	int lo;
	int cnt;
	int dir;
} recBitonicSortData;

typedef struct impBitonicSortThreadData{
	int i;
	int N;
	int j;
	int k;
} impBitonicSortThreadData;

pthread_t *threads; //array that holds pointers to all threads
int runningThreads = 0; //number of threads currently running
pthread_mutex_t runningThreadsMutex = PTHREAD_MUTEX_INITIALIZER;

void getArgs(int argc, char** argv);
void init(void);
void print(void);
void sort(void);
void test(void);
void qSortTest(void);
void exchange(int i, int j);
void compare(int i, int j, int dir);
void bitonicMerge(int lo, int cnt, int dir);
void *recBitonicSort(void *threadArgs);
void *impBitonicSortThread(void * threadArgs);
void impBitonicSort(void);
int qSortAscendingCompFuncWithTest (const void * a, const void * b);
int qSortAscending (const void * a, const void * b);
int qSortDescending (const void * a, const void * b);

/** the main program **/ 
int main(int argc, char **argv) {
	getArgs(argc, argv);
	a = (int *) malloc(N * sizeof(int));
	randSeed = (unsigned) time(NULL);

	//allocates space for an array that holds pointers to all threads
	threads = (pthread_t *)malloc(sizeof(pthread_t)*threadsNum);

	//Sorts using the qSort algorithm
	init();

	gettimeofday (&startwtime, NULL);
	qsort(a, N, sizeof(int), qSortAscending);
	gettimeofday (&endwtime, NULL);

	seq_time = (double)((endwtime.tv_usec - startwtime.tv_usec)/1.0e6
		+ endwtime.tv_sec - startwtime.tv_sec);

	printf("qSort wall clock time = %f\n\n", seq_time);

	//Sorts using the recursive bitonic algorithm
	init();

	gettimeofday (&startwtime, NULL);
	sort();
	gettimeofday (&endwtime, NULL);

	seq_time = (double)((endwtime.tv_usec - startwtime.tv_usec)/1.0e6
		+ endwtime.tv_sec - startwtime.tv_sec);

	printf("\nRecursive wall clock time = %f\n", seq_time);

	test();
	qSortTest();

	free(threads);
	free(a);
}

/** -------------- SUB-PROCEDURES  ----------------- **/ 

void getArgs(int argc, char** argv){
	if (argc != 3) {
		printf("Usage: %s p q\nwhere:\n\tP=2^p is the the number of threads(power of two)\n\tN=2^q is problem size (power of two)\n", 
			argv[0]);
		exit(1);
	}
	threadsNum = 1<<atoi(argv[1]);
	N = 1<<atoi(argv[2]);
}

/** procedure test() : verify sort results **/
void test() {
	int pass = 1;
	int i;
	for (i = 1; i < N; i++) {
		pass &= (a[i-1] <= a[i]);
	}
	printf("\tTEST\t\t%s\n",(pass) ? "PASSed" : "FAILed");
}

/** procedure qSortTest() : verify sort results using qsort method **/
void qSortTest(){
	sortPass = true;
	qsort(a, N, sizeof(int), qSortAscendingCompFuncWithTest);
	printf("\tQSORT TEST\t%s\n",(sortPass) ? "PASSed" : "FAILed");
}

/** procedure init() : initialize array "a" with data **/
void init() {
	srand(randSeed);
	int i;
	for (i = 0; i < N; i++) {
		a[i] = rand() % N; // (N - i);
	}
}

/** procedure  print() : print array elements **/
void print() {
	int i;
	for (i = 0; i < N; i++) {
		printf("%d\n", a[i]);
	}
	printf("\n");
}

/** INLINE procedure exchange() : pair swap **/
inline void exchange(int i, int j) {
	int t;
	t = a[i];
	a[i] = a[j];
	a[j] = t;
}

/** procedure compare() 
   The parameter dir indicates the sorting direction, ASCENDING 
   or DESCENDING; if (a[i] > a[j]) agrees with the direction, 
   then a[i] and a[j] are interchanged.
**/
inline void compare(int i, int j, int dir) {
	if (dir==(a[i]>a[j])) 
		exchange(i,j);
}

/** Procedure bitonicMerge() 
   It recursively sorts a bitonic sequence in ascending order, 
   if dir = ASCENDING, and in descending order otherwise. 
   The sequence to be sorted starts at index position lo,
   the parameter cbt is the number of elements to be sorted. 
 **/
void bitonicMerge(int lo, int cnt, int dir) {
	if (cnt>1) {
		int k=cnt/2;
		int i;
		for (i=lo; i<lo+k; i++){
			compare(i, i+k, dir);
		}

		bitonicMerge(lo, k, dir);
		bitonicMerge(lo+k, k, dir);
	}
}

/** function recBitonicSort() 
    first produces a bitonic sequence by recursively sorting 
    its two halves in opposite sorting orders, and then
    calls bitonicMerge to make them in the same order 
 **/
void *recBitonicSort(void *threadArgs) {
	recBitonicSortData *thisData = (recBitonicSortData *) threadArgs;
	int lo = thisData->lo;
	int cnt = thisData->cnt;
	int dir = thisData->dir;

	if (cnt>1) {
		if (false){ //small sub-arrays best sorted using qsort
			qsort(a+lo, cnt, sizeof(int), (dir == 1 ? qSortAscending : qSortDescending));
		} else {
			//plitts problem in two halfs and tries to create a new thread for each
			//holds each half's index (from threads array)
			int firstHalf = -1, myOtherHalf = -1;
			//attribute for joinable threads
			pthread_attr_t attr;
			pthread_attr_init(&attr);
			pthread_attr_setdetachstate(&attr, PTHREAD_CREATE_JOINABLE);

			int k=cnt/2;
			pthread_mutex_lock(&runningThreadsMutex); //locks running threads variable
			if (runningThreads < threadsNum){
				//arguments for first half
				recBitonicSortData newThreadArgs;
				newThreadArgs.lo = lo;
				newThreadArgs.cnt = k;
				newThreadArgs.dir = ASCENDING;
				//tries to create a new thread with a pointer to it saved inside array "threads"
				if (pthread_create(&threads[runningThreads], &attr, recBitonicSort, (void *) &newThreadArgs)){
					printf("Error creating thread. Aborting.");
					pthread_mutex_unlock(&runningThreadsMutex);
					exit(1);
				} else { //thread was successfully created
					//keeps thread's index and updates running threads counter
					firstHalf = runningThreads;
					++runningThreads;
					pthread_mutex_unlock(&runningThreadsMutex);
				}
			} else { //max threads number reached, does the job sequentially using qsort
				pthread_mutex_unlock(&runningThreadsMutex);
				qsort(a+lo, k, sizeof(int), qSortAscending);
			}

			pthread_mutex_lock(&runningThreadsMutex); //locks running threads variable
			if (runningThreads < threadsNum){
				//arguments for second half
				recBitonicSortData newThreadArgs;
				newThreadArgs.lo = lo+k;
				newThreadArgs.cnt = k;
				newThreadArgs.dir = DESCENDING;
				//tries to create a new thread with a pointer to it saved inside array "threads"
				if (pthread_create(&threads[runningThreads], &attr, recBitonicSort, (void *) &newThreadArgs)){
					printf("Error creating thread. Aborting.");
					pthread_mutex_unlock(&runningThreadsMutex);
					exit(1);
				} else { //thread was successfully created
					//keeps thread's index and updates running threads counter
					myOtherHalf = runningThreads;
					++runningThreads;
					pthread_mutex_unlock(&runningThreadsMutex);
				}
			} else { //max threads number reached, does the job sequentially using qsort
				pthread_mutex_unlock(&runningThreadsMutex);
				qsort(a+lo+k, k, sizeof(int), qSortDescending);
			}

			pthread_attr_destroy(&attr);

			//achieves synchronization
			if (firstHalf != -1){
				pthread_join (threads[firstHalf], NULL);
			}
			if (myOtherHalf != -1){
				pthread_join (threads[myOtherHalf], NULL);
			}
			bitonicMerge(lo, cnt, dir);
		}
	}
}

/** function sort() 
   Caller of recBitonicSort for sorting the entire array of length N 
   in ASCENDING order
 **/
void sort() {
	recBitonicSortData newThreadArgs;
	newThreadArgs.threadID = 1;
	newThreadArgs.lo = 0;
	newThreadArgs.cnt = N;
	newThreadArgs.dir = ASCENDING;
	recBitonicSort((void *) &newThreadArgs);
}

void *impBitonicSortThread(void * threadArgs){
	//get variable values from struct
	impBitonicSortThreadData *thisThreadData = (impBitonicSortThreadData *) threadArgs;
	int i;
	int N = thisThreadData->N;
	int j = thisThreadData->j;
	int k = thisThreadData->k;

	for (i; i<N; i++) {
		int ij=i^j;
		if ((ij)>i) {
			if ((i&k)==0 && a[i] > a[ij]) 
				exchange(i,ij);
			if ((i&k)!=0 && a[i] < a[ij])
				exchange(i,ij);
		}
	}
	pthread_exit(NULL);
}

/*
  imperative version of bitonic sort
*/
void impBitonicSort() {
	int t,j,k;

	//attribute for joinable threads
	pthread_attr_t attr;
	pthread_attr_init(&attr);
	pthread_attr_setdetachstate(&attr, PTHREAD_CREATE_JOINABLE);
	int threadN = N/threadsNum; //splitts work into 2^q sub-problems

	for (k=2; k<=N; k+=k) {
		for (j=k>>1; j>0; j=j>>1) {
			//creates a new thread for each sub-problem
			for (t = 0; t < threadsNum; ++t){
				impBitonicSortThreadData newThreadArgs;
				newThreadArgs.i = t * threadN;
				newThreadArgs.N = (t + 1) * threadN;
				newThreadArgs.j = j;
				newThreadArgs.k = k;
				int rc = pthread_create(&threads[t], &attr, impBitonicSortThread, (void *) &newThreadArgs);
				if (rc){
					printf("ERROR; return code from pthread_create() is %d\n", rc);
					exit(-1);
				}
			}

			//achieves synchronization
			for(t = 0; t < threadsNum; ++t) {
				pthread_join(threads[t], NULL);
			}
		}
	}
	pthread_attr_destroy(&attr);
}

/** function used by qsort for comparing as well as testing **/
int qSortAscendingCompFuncWithTest (const void * a, const void * b) {
	int result = ( *(int*)a - *(int*)b );
	if (result > 0){
		sortPass = false;
	}
	return result;
}

/** function used by qsort for comparing **/
int qSortAscending (const void * a, const void * b) {
	return ( *(int*)a - *(int*)b );
}

/** function used by qsort for comparing **/
int qSortDescending (const void * a, const void * b) {
	return ( *(int*)b - *(int*)a );
}