path: root/src/quicksort.c

#include <assert.h>
#include <stdio.h>
#include <stdlib.h>
#include <time.h>

static int
bigrand(void)
{
	int x = (rand() << 24) | (rand() << 16) | (rand() << 8) | rand();
	if ( x < 0 )
		return -x;
	return x;
}

static int
randint(int l, int u)
{
	return l + bigrand() % (u - l + 1);
}

typedef int T;

static const int CUTOFF = 128;

static inline void
swap(T a[], int i, int j)
{
	T t  = a[i];
	a[i] = a[j];
	a[j] = t;
}

static void
insertsort(T a[], int n)
{
	int i, j;
	for ( i = 1; i < n; ++i ) {
		T t = a[i];
		for ( j = i; j > 0 && a[j - 1] > t; --j )
			a[j] = a[j - 1];
		a[j] = t;
	}
}

static void
quicksort(T a[], int n)
{
	int i, last;

	if ( n <= CUTOFF )
		return;

	swap(a, 0, bigrand() % n);

	last = 0;
	for ( i = 1; i < n; ++i )
		if ( a[i] < a[0] )
			swap(a, ++last, i);

	swap(a, 0, last);

	quicksort(a, last);
	quicksort(a + last + 1, n - last - 1);
}

void
my_quicksort(T a[], int n)
{
	quicksort(a, n);
	insertsort(a, n);
}

static int
cmp(const void *a, const void *b)
{
	const int *pa = (const T *) a;
	const int *pb = (const T *) b;

	if ( *pa < *pb )
		return -1;
	if ( *pa > *pb )
		return 1;
	return 0;
}

void
c_quicksort(T a[], int n)
{
	assert(n >= 0);
	qsort(a, (size_t) n, sizeof(T), cmp);
}

/* This function takes last element as pivot, places
   the pivot element at its correct position in sorted
    array, and places all smaller (smaller than pivot)
   to left of pivot and all greater elements to right
   of pivot */
int
partition(int arr[], int low, int high)
{
	int pivot = arr[high]; // pivot
	int i     = (low - 1); // Index of smaller element

	for ( int j = low; j <= high - 1; j++ ) {
		// If current element is smaller than or
		// equal to pivot
		if ( arr[j] <= pivot ) {
			i++; // increment index of smaller element
			swap(arr, i, j);
		}
	}
	swap(arr, i + 1, high);
	return (i + 1);
}

/* The main function that implements QuickSort
 arr[] --> Array to be sorted,
  low  --> Starting index,
  high  --> Ending index */
void
quickSort(int arr[], int low, int high)
{
	while ( low < high ) {
		/* pi is partitioning index, arr[p] is now
		   at right place */
		int pi = partition(arr, low, high);

		if ( pi - low < high - pi ) {
			quickSort(arr, low, pi - 1);
			low = pi + 1;
		}
		else {
			quickSort(arr, pi + 1, high);
			high = pi - 1;
		}
	}
}

void
g4g_quicksort(T a[], int n)
{
	quickSort(a, 0, n - 1);
}

static void
three_way_quicksort(int a[], int l, int r)
{
	int k;
	T   v = a[r];

	if ( r <= l )
		return;

	int i = l - 1, j = r, p = l - 1, q = r;

	for ( ;; ) {
		while ( a[++i] < v )
			;
		while ( v < a[--j] )
			if ( j == l )
				break;
		if ( i >= j )
			break;

		swap(a, i, j);

		if ( a[i] == v )
			swap(a, ++p, i);
		if ( v == a[j] )
			swap(a, --q, j);
	}
	swap(a, i, r);
	j = i - 1;
	i = i + 1;

	for ( k = l; k <= p; ++k, --j )
		swap(a, k, j);
	for ( k = r - 1; k >= q; --k, ++i )
		swap(a, k, i);

	three_way_quicksort(a, l, j);
	three_way_quicksort(a, i, r);
}

void
sed_quicksort(int a[], int n)
{
	three_way_quicksort(a, 0, n - 1);
}

/* ====================== HEAPSORT ======================= */

#define LEFT(idx) (idx * 2 + 1)
#define RIGHT(idx) (idx * 2 + 2)
#define PARENT(idx) ((idx - 1) / 2)

static void
fixdown(T heap[], int i, int n)
{
	for ( ;; ) {
		const int l = LEFT(i);
		const int r = RIGHT(i);
		int       m = i;

		if ( l < n && heap[m] < heap[l] )
			m = l;

		if ( r < n && heap[m] < heap[r] )
			m = r;

		if ( m == i )
			break;

		swap(heap, m, i);

		i = m;
	}
}

static void
heapify(T heap[], int n)
{
	for ( int i = n / 2 - 1; i >= 0; --i )
		fixdown(heap, i, n);
}

void
my_heapsort(T a[], int n)
{
	heapify(a, n);

	for ( int i = n - 1; i >= 0; --i ) {
		swap(a, 0, i);
		fixdown(a, 0, i);
	}
}

void
heapsort_bu(T *data, int n) // zu sortierendes Feld und seine Länge
{
	T   val;
	int parent, child;
	int root  = n >> 1; // erstes Blatt im Baum
	int count = 0;      // Zähler für Anzahl der Vergleiche

	for ( ;; ) {
		if ( root ) { // Teil 1: Konstruktion des Heaps
			parent = --root;
			val    = data[root]; // zu versickernder Wert
		}
		else if ( --n ) {      // Teil 2: eigentliche Sortierung
			val     = data[n]; // zu versickernder Wert vom Heap-Ende
			data[n] = data[0]; // Spitze des Heaps hinter den Heap in
			parent  = 0;       //   den sortierten Bereich verschieben
		}
		else // Heap ist leer; Sortierung beendet
			break;

		while ( (child = (parent + 1) << 1) < n )         // zweites (!) Kind;
		{                                                 // Abbruch am Ende des Heaps
			if ( ++count, data[child - 1] > data[child] ) // größeres Kind wählen
				--child;

			data[parent] = data[child]; // größeres Kind nach oben rücken
			parent       = child;       // in der Ebene darunter weitersuchen
		}

		if ( child == n )                          // ein einzelnes Kind am Heap-Ende
		{                                          //   ist übersprungen worden
			if ( ++count, data[--child] >= val ) { // größer als der zu versick-
				data[parent] = data[child];        //   ernde Wert, also noch nach oben
				data[child]  = val;                // versickerten Wert eintragen
				continue;
			}

			child = parent; // 1 Ebene nach oben zurück
		}
		else {
			if ( ++count, data[parent] >= val ) { // das Blatt ist größer als der
				data[parent] = val;               //   zu versickernde Wert, der damit
				continue;                         //   direkt eingetragen werden kann
			}

			child = (parent - 1) >> 1; // 2 Ebenen nach oben zurück
		}

		while ( child != root ) // maximal zum Ausgangspunkt zurück
		{
			parent = (child - 1) >> 1; // den Vergleichswert haben wir bereits
			//   nach oben verschoben
			if ( ++count, data[parent] >= val ) // größer als der zu versickernde
				break;                          //   Wert, also Position gefunden

			data[child] = data[parent]; // Rückverschiebung nötig
			child       = parent;       // 1 Ebene nach oben zurück
		}

		data[child] = val; // versickerten Wert eintragen
	}
	printf("%d\n", count);
}

/*----------------------------------------------------------------------*/
/*                         BOTTOM-UP HEAPSORT                           */
/* Written by J. Teuhola <teuhola@cs.utu.fi>; the original idea is      */
/* probably due to R.W. Floyd. Thereafter it has been used by many      */
/* authors, among others S. Carlsson and I. Wegener. Building the heap  */
/* bottom-up is also due to R. W. Floyd: Treesort 3 (Algorithm 245),    */
/* Communications of the ACM 7, p. 701, 1964.                           */
/*----------------------------------------------------------------------*/

#define element float

/*-----------------------------------------------------------------------*/
/* The sift-up procedure sinks a hole from v[i] to leaf and then sifts   */
/* the original v[i] element from the leaf level up. This is the main    */
/* idea of bottom-up heapsort.                                           */
/*-----------------------------------------------------------------------*/
static void
siftup(T v[], int i, int n)
{
	int j, start;
	T   x;

	start = i;
	x     = v[i];
	j     = i << 1;
	while ( j <= n ) {
		if ( j < n )
			if ( v[j] < v[j + 1] )
				j++;
		v[i] = v[j];
		i    = j;
		j    = i << 1;
	}
	j = i >> 1;
	while ( j >= start ) {
		if ( v[j] < x ) {
			v[i] = v[j];
			i    = j;
			j    = i >> 1;
		}
		else
			break;
	}
	v[i] = x;
} /* End of siftup */

/*----------------------------------------------------------------------*/
/* The heapsort procedure; the original array is r[0..n-1], but here    */
/* it is shifted to vector v[1..n], for convenience.                     */
/*----------------------------------------------------------------------*/
void
bottom_up_heapsort(T r[], int n)
{
	int k;
	T   x;
	T  *v;

	v = r - 1; /* The address shift */

	/* Build the heap bottom-up, using siftup. */
	for ( k = n >> 1; k > 1; k-- )
		siftup(v, k, n);

	/* The main loop of sorting follows. The root is swapped with the last  */
	/* leaf after each sift-up. */
	for ( k = n; k > 1; k-- ) {
		siftup(v, 1, k);
		x    = v[k];
		v[k] = v[1];
		v[1] = x;
	}
} /* End of bottom_up_heapsort */

/* ====================== HEAPSORT ======================= */

/* ====================== INTROSORT ======================= */

static void
introsort(T a[], int n, int h)
{
	int i, last;

	if ( n <= CUTOFF )
		return;

	if ( --h == 1 ) {
		my_heapsort(a, n);
		return;
	}

	swap(a, 0, bigrand() % n);

	last = 0;
	for ( i = 1; i < n; ++i )
		if ( a[i] < a[0] )
			swap(a, ++last, i);

	swap(a, 0, last);

	introsort(a, last, h);
	introsort(a + last + 1, n - last - 1, h);
}

void
my_introsort(T a[], int n)
{
	int h = 1;

	for ( int nn = 1; nn < n; nn <<= 1 )
		++h;

	introsort(a, n, h);
	insertsort(a, n);
}

static void
pp_quicksort_impl(T a[], int l, int u)
{
	if ( u - l < CUTOFF )
		return;

	swap(a, l, randint(l, u));

	T   t = a[l];
	int i = l;
	int j = u + 1;
	for ( ;; ) {
		do
			i++;
		while ( /*i <= u &&*/ a[i] < t );
		do
			j--;
		while ( a[j] > t );
		if ( i > j )
			break;
		swap(a, i, j);
	}
	swap(a, l, j);
	pp_quicksort_impl(a, l, j - 1);
	pp_quicksort_impl(a, j + 1, u);
}

void
pp_quicksort(T a[], int n)
{
	pp_quicksort_impl(a, 0, n - 1);
	insertsort(a, n);
}

static void
pp_quicksort_impl_it(T a[], int l, int u, int h)
{
	if ( --h == 1 ) {
		my_heapsort(a + l, u - l + 1);
		return;
	}

	while ( u - l >= CUTOFF ) {
		swap(a, l, randint(l, u));

		T   t = a[l];
		int i = l;
		int j = u + 1;

		for ( ;; ) {
			do
				i++;
			while ( /*i <= u && */ a[i] < t );
			do
				j--;
			while ( a[j] > t );
			if ( i > j )
				break;
			swap(a, i, j);
		}
		swap(a, l, j);

		if ( j - l < u - j ) {
			pp_quicksort_impl_it(a, l, j - 1, h);
			l = j + 1;
		}
		else {
			pp_quicksort_impl_it(a, j + 1, u, h);
			u = j - 1;
		}
	}
}

void
pp_quicksort_it(T a[], int n)
{
	int h = 1;

	for ( int nn = 1; nn < n; nn <<= 1 )
		++h;

	pp_quicksort_impl_it(a, 0, n - 1, h);
	insertsort(a, n);
}

int
binary_search(T x, T v[], int n)
{
	int low, high;

	low  = 0;
	high = n - 1;
	while ( low <= high ) {
		int mid = low + ((high - low) / 2);
		if ( x > v[mid] )
			low = mid + 1;
		else if ( x < v[mid] )
			high = mid - 1;
		else
			return mid;
	}
	return -1;
}

int
lower_bound(T x, T v[], int n)
{
	int low, high;

	low  = 0;
	high = n;
	while ( low < high ) {
		int mid = low + ((high - low) / 2);

		if ( x > v[mid] )
			low = mid + 1;
		else
			high = mid;
	}
	return x == v[low] ? low : -1;
}

int
upper_bound(T x, T v[], int n)
{
	int low, high;

	low  = 0;
	high = n;
	while ( low < high ) {
		int mid = low + ((high - low) / 2);

		if ( x >= v[mid] )
			low = mid + 1;
		else
			high = mid;
	}
	return (low > 0 && x == v[low - 1]) ? low - 1 : -1;
}

/* TESTTREIBER */

void
gen_testset_random(T a[], int n)
{
	for ( int i = 0; i < n; ++i )
		a[i] = bigrand();
}

void
gen_testset_random2(T a[], int n)
{
	for ( int i = 0; i < n; ++i )
		a[i] = bigrand() % 100;
}

void
gen_testset_asc(T a[], int n)
{
	for ( int i = 0; i < n; ++i )
		a[i] = i;
}

void
gen_testset_desc(T a[], int n)
{
	for ( int i = n; i >= 0; --i )
		a[i] = i;
}

void
gen_testset_unique(T a[], int n)
{
	for ( int i = 0; i < n; ++i )
		a[i] = 1;
}

void
test_sorting(T a[], int n)
{
	for ( int i = 1; i < n; ++i )
		if ( a[i - 1] > a[i] ) {
			puts("Fehler!");
			exit(EXIT_SUCCESS);
		}
}

void
do_one_test(int n, void (*do_sort)(T a[], int n), void (*gen_testset)(T a[], int n))
{
	T      *a;
	clock_t start, ende;

	a = malloc(sizeof(T) * (size_t) n);

	(*gen_testset)(a, n);

	start = clock();
	(*do_sort)(a, n);
	ende = clock();

	test_sorting(a, n);

	free(a);

	printf("%10.3lf sec", ((double) ende - start) / CLOCKS_PER_SEC);
	fflush(stdout);
}

void
do_all_tests(const char *msg, int n, void (*do_sort)(T a[], int n))
{
	printf("%-15.15s n=%-10d: ", msg, n);
	do_one_test(n, do_sort, gen_testset_random);
	do_one_test(n, do_sort, gen_testset_random2);
	do_one_test(n, do_sort, gen_testset_asc);
	do_one_test(n, do_sort, gen_testset_desc);
	do_one_test(n, do_sort, gen_testset_unique);
	putchar('\n');
}

/* ENDE TESTTREIBER */

int
main(void)
{
	const int n = 20000000;

	srand(time(0));
	do_all_tests("my_heapsort", n, my_heapsort);
	do_all_tests("heapsort_bu", n, heapsort_bu);
	do_all_tests("bottom_up_heapsort", n, bottom_up_heapsort);
	do_all_tests("c_quicksort", n, c_quicksort);
	// do_all_tests("g4g_quicksort", n, g4g_quicksort);
	// do_all_tests("sed_quicksort", n, sed_quicksort);
	do_all_tests("my_introsort", n, my_introsort);
	// do_all_tests("my_quicksort", n, my_quicksort);
	do_all_tests("pp_quicksort", n, pp_quicksort);
	do_all_tests("pp_quicksort_it", n, pp_quicksort_it);

	return 0;
}

#if 0
int _main(void)
{
	T array[20];

	for ( int i = 0; i != NELEM(array); ++i )
		array[i] = rand() % 10;

	sort(array, NELEM(array));
	print(array, NELEM(array));

	int value;
	while ( scanf("%d", &value) == 1 && value != -1 ) {
		int bs = binary_search(value, array, NELEM(array));
		int lb = lower_bound(value, array, NELEM(array));
		int ub = upper_bound(value, array, NELEM(array));

		printf("bs: %d - lb: %d - ub: %d\n", bs, lb, ub);
	}

	return EXIT_SUCCESS;
}
#endif