/* Copyright 1995 Jonathan C. Hardwick */

#include <stdlib.h>
#include <limits.h>
#include <assert.h>
#include <mpi.h>
#include "globals.h"
#include "main.h"
#include "serial.h"
#include "utils.h"

/* Vectorized inplace quicksort of asrc, using atmp as scratch space. */

void
inplace_sort (double asrc[], double atmp[], int neltHere, int oddeven)
{
#ifdef STEAL
  MPI_Status status;
  int rnkSend;
#endif
#ifdef PIVOT
  double pivot0, pivot1, pivot2;
#endif
  double pivot, elt, eltFirst, eltSave;
  int neltLess, neltGreater;
  int celtLess, celtGreater;
  int i;

  /* Terminate if we have 0 or 1 elements. */
  if (neltHere < 2) {
    return;
  }

  /* Terminate if all elements are equal. */

  /* Insert a sentinel to force termination. */
  eltFirst = asrc [0];
  eltSave = asrc [neltHere - 1];
  asrc [neltHere - 1] = (eltFirst * 2) + 1;

  /* Loop until we find a nonequal element. */
  for (i = 1; asrc [i] == eltFirst; i++) {
    /* NULL */
  }

  /* Remove the sentinel again. */
  asrc [neltHere - 1] = eltSave;

  /* Did we finish the loop? */
  if ((i == neltHere - 1) && (eltSave == eltFirst)) {
    return;
  }

#ifdef PIVOT
  if (neltHere < 8) {
    /* Don't bother taking the branches necessary to pick a good pivot
     * if we only have a few elements. */
    pivot = asrc [neltHere / 2];
  } else {
    pivot0 = asrc [neltHere / 4];
    pivot1 = asrc [neltHere / 2];
    pivot2 = asrc [3 * neltHere / 4];

    if (pivot0 < pivot1) { 	/* 0 < 1 */
      if (pivot1 < pivot2) {	/* 0 < 1 < 2 */
	pivot = pivot1;
      } else {			/* 0 < 1, 2 < 1 */
	if (pivot0 < pivot2) {	/* 0 < 2 < 1 */
	  pivot = pivot2;
	} else {		/* 2 < 0 < 1 */
	  pivot = pivot0;
	}
      }
    } else {			/* 1 < 0 */
      if (pivot0 < pivot2) {	/* 1 < 0 < 2 */
	pivot = pivot0;
      } else {			/* 1 < 0, 2 < 0 */
	if (pivot1 < pivot2) {	/* 1 < 2 < 0 */
	  pivot = pivot2;
	} else {		/* 2 < 1 < 0 */
	  pivot = pivot1;
	}
      }
    }
  }

#else  /* !PIVOT */
  pivot = asrc [neltHere / 2];
#endif

  /* Copy asrc to atmp, counting elements above/below pivot. */
  neltLess = 0;
  if (oddeven) {
    for (i = 0; i < neltHere; i++) {
      elt = asrc [i];
      if (elt <= pivot) neltLess++;
      atmp [i] = elt;
    }
  } else {
    for (i = 0; i < neltHere; i++) {
      elt = asrc [i];
      if (elt < pivot) neltLess++;
      atmp [i] = elt;
    }
  }

  neltGreater = neltHere - neltLess;

  /* Now partition from atmp back into asrc. */
  celtLess = 0;
  celtGreater = neltLess;

  if (oddeven) {
    for (i = 0; i < neltHere; i++) {
      elt = atmp [i];
      if (elt <= pivot) {
	asrc [celtLess++] = elt;
      } else {
	asrc [celtGreater++] = elt;
      }
    }
  } else {
    for (i = 0; i < neltHere; i++) {
      elt = atmp [i];
      if (elt < pivot) {
	asrc [celtLess++] = elt;
      } else {
	asrc [celtGreater++] = elt;
      }
    }
  }
  assert (celtLess == neltLess);
  assert (celtGreater == neltHere);

  oddeven = (oddeven == 0);

  /* Recurse on Less and Greater. */

#ifdef STEAL
  if (neltLess > gcutoff) {
    MPE_Log_event (ASSIGN_START, 0, NULL);
    MPI_Send (&neltLess, 1, MPI_INT, 0, WKR_MGR_HELP_TAG, gcomWorld);
    MPI_Recv (&rnkSend, 1, MPI_INT, 0, MPI_ANY_TAG, gcomWorld,
	      &status);
    MPE_Log_event (ASSIGN_END, 0, NULL);

    if (status.MPI_TAG == MGR_WKR_TO_TAG) {
      MPE_Log_event (SEND_START, 0, NULL);
      MPI_Send (asrc, neltLess, MPI_DOUBLE, rnkSend, WKR_WKR_DATA_TAG,
		gcomWorld);
      MPE_Log_event (SEND_END, 0, NULL);
    } else {
      assert (status.MPI_TAG == MGR_WKR_TOUGH_TAG);
      inplace_sort (asrc, atmp, neltLess, oddeven);
    }
  } else {
    inplace_sort (asrc, atmp, neltLess, oddeven);
  }       

  if (neltGreater > gcutoff) {
    MPE_Log_event (ASSIGN_START, 0, NULL);
    MPI_Send (&neltGreater, 1, MPI_INT, 0, WKR_MGR_HELP_TAG,
	      gcomWorld);
    MPI_Recv (&rnkSend, 1, MPI_INT, 0, MPI_ANY_TAG, gcomWorld,
	      &status);
    MPE_Log_event (ASSIGN_END, 0, NULL);

    if (status.MPI_TAG == MGR_WKR_TO_TAG) {
      MPE_Log_event (SEND_START, 0, NULL);
      MPI_Send (&asrc [neltLess], neltGreater, MPI_DOUBLE, rnkSend,
		WKR_WKR_DATA_TAG, gcomWorld);
      MPE_Log_event (SEND_END, 0, NULL);
   } else {
      assert (status.MPI_TAG == MGR_WKR_TOUGH_TAG);
      inplace_sort (&asrc [neltLess], &atmp [neltLess], neltGreater,
		    oddeven);
    }
  } else {
    inplace_sort (&asrc [neltLess], &atmp [neltLess], neltGreater,
		  oddeven);
  }

#else  /* !STEAL */
  inplace_sort (asrc, atmp, neltLess, oddeven);
  inplace_sort (&asrc [neltLess], &atmp [neltLess], neltGreater,
		oddeven);
#endif
}


void 
setup_sort (int seed, double asrc[], int nelt)
{
  int i, neltHere;

  srandom ((seed * gnprocWkr) + grnkWkr);
  neltHere = nelt_on_proc (grnkWkr, nelt, gnprocWkr);
  for (i = 0; i < neltHere; i++) {
    asrc [i] = ((double) random()) / (double) LONG_MAX;
  }
}


void 
verify_sorted (double asrc[], int nelt)
{
  int i, neltHere;

  neltHere = nelt_on_proc (grnkWkr, nelt, gnprocWkr);
  for (i = 1; i < neltHere; i++) {
    if (asrc [i-1] > asrc [i]) {
      printf ("%d failed on elements %d [%g] and %d [%g]\n", grnkWkr,
	      i-1, asrc[i-1], i, asrc[i]);
    }
  }
}