/******************************************************************
*
* Uebungen zu "Parallele und verteilte Algorithmen"
* SS 2007
*
* Zu Aufgabe 5.3: Odd-Even-Transposition Sort 
*             
*******************************************************************/


#include "mpi.h"
#include <limits.h> // INT_MAX
#include <stdlib.h>
#include <stdio.h>
#include <string.h>

#define LIMIT 1023

// Message-Tags
#define FROMLEFT  1
#define FROMRIGHT 2

// debugging:
#ifdef OUTPUT
# define DEBUG(s, arg) printf("%d" s "\n",rank, arg)
#else
# define DEBUG(s,arg) /* Nothing */
#endif

// Globale Variablen: 
int size; // MPI
int rank; // MPI

// shutdown in case of errors
void error_exit(char* message) {
  fprintf(stderr,"[%x] Fehler: %s\n", rank, message);
  MPI_Finalize(); // MPI verlassen
  exit(-1);   // mit Fehler terminieren
}

//  Exchange the contents of  * X  and  * Y .
void  Swap  (int* X, int* Y)
  {
   int Save;
   Save = * X;
   * X = * Y;
   * Y = Save;
   return;
  } // Swap

//  Sort array  A [0 .. (N-1)]  with bubblesort
void  Sort  (int* A, int N)
  {
   int  i, j;
   for  (i = 0;  i < N - 1;  i ++)
     for  (j = 0;  j < N - i - 1;  j ++)
       if  (A [j] > A [j + 1])           // wrong order
           Swap (& A [j], & A [j + 1]);
   return;
  } // Sort


// merge sorted arrays A and B into result (must be allocated!)
void Merge(int* A, int a_size, int* B, int b_size, int* result) {
  int ai=0,bi=0,ri=0;
  
  while (ai < a_size && bi < b_size) {
    if (A[ai] < B[bi]) 
      result[ri++] = A[ai++];
    else 
      result[ri++] = B[bi++];
  }
  if (ai==a_size)
    memcpy(result+ri, B+bi, sizeof(int)*(b_size-bi));
  else
    memcpy(result+ri, A+ai, sizeof(int)*(a_size-ai));

  return;
}

// fill array with n random elements between 0 and LIMIT
void fillRnd(int* array, int n) {

  srand((unsigned) array); // :-)

  while (n) {
    array[--n] = rand() & LIMIT;
  }
  return;
}

void  ArrayOutput(int* array, int size) {
  int j;
  printf("array:[");
  for ( j = 0; j < size; j++) {
    printf("%d%s", array[j], j==size-1?"]\n":", ");
  }
}

int main(int argc, char** argv) {

  int nElements;
  int addElements;
  int* numbers, *numbers2;
  int iteration;

  MPI_Init(&argc, &argv);
  MPI_Comm_rank(MPI_COMM_WORLD, &rank);
  MPI_Comm_size(MPI_COMM_WORLD, &size);

  DEBUG("Process starting, size is %d",size);

  nElements = argc>1?
                atoi(argv[1]): // arg. 1 gives no. of elements
                size*20;       // or use default

  if (!rank) {
    printf("Sorting %d random numbers with OET-Sort\n", 
	   nElements);
  }

  // ensure that nElements can be divided without problems
  addElements = nElements % size;
  if (addElements) // need correction?
    addElements = size - addElements;
  nElements += addElements;
  // additional elements could be initialized LIMIT+1, and removed
  // from last process after sorting (here: random values anyway)

  // now start
  nElements /= size; // local no. of elements

  // allocate additional space for exchange
  numbers = (int*) malloc(sizeof(int)*2*nElements);
  // fill with values (only first half)
  fillRnd(numbers, nElements);

  DEBUG("array created", NULL);

#ifdef OUTPUT
  if (!rank) {
    numbers2 = (int*) malloc(sizeof(int)*size*nElements);
  }
  MPI_Gather(numbers, nElements, MPI_INT,
	     numbers2, nElements, MPI_INT,
	     0, MPI_COMM_WORLD);
  if (!rank) {
    ArrayOutput(numbers2, size*nElements);
    free(numbers2);
  }
#else 
  // sync
  MPI_Barrier(MPI_COMM_WORLD);
#endif

  // local sort
  Sort(numbers,nElements);

  DEBUG("array sorted", NULL);

  // allocate space to receive and merge
  numbers2 = (int*) malloc(sizeof(int)*2*nElements); 

  // iterate OET
  iteration = 0;
  while (iteration < size) {
    DEBUG("Iteration %d",iteration);

    // send to left, receive from left
    if (((iteration++) + rank) & 1) {
      // odd iteration, or odd rank, not  both

      if (rank-1 >=0) { // neighbour exists
	DEBUG("Sending/Receiving array", NULL);
	// send, then receive
	MPI_Send(numbers,nElements, MPI_INT,
		 rank-1, FROMLEFT, MPI_COMM_WORLD);
	MPI_Recv(numbers, nElements, MPI_INT,
		 rank-1, FROMRIGHT, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
      }
    } else {
      if (rank+1 < size) { // neighbour exists
	int* pointer;

	DEBUG("Receiving/Merging/Sending array", NULL);
	// receive (smaller values from neighbour)
	MPI_Recv(numbers+nElements, nElements, MPI_INT,
		 rank+1, FROMLEFT, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
	// merge into own list
	Merge(numbers,nElements,numbers+nElements,nElements, numbers2);
	// send back higher half
	MPI_Send(numbers2+nElements, nElements, MPI_INT,
		 rank+1, FROMRIGHT, MPI_COMM_WORLD);
	
	// swap array roles
	pointer = numbers2;
	numbers2 = numbers;
	numbers = pointer;
      }
    }
  }

  printf("%d: done\n", rank);
#ifdef OUTPUT
  if (!rank) {
    realloc(numbers2, sizeof(int)*nElements*size);
  }
  MPI_Gather(numbers,nElements, MPI_INT,
	     numbers2, nElements, MPI_INT,
	     0, MPI_COMM_WORLD);
  if (!rank) {
    ArrayOutput(numbers2, size*nElements);
  }
#endif

  // free allocated memory
  free(numbers);
  free(numbers2);

  MPI_Finalize();
}