Note: the next iteration at Natural merge sort - follow-up
I have compiled natural merge sort algorithm from Java to C++.
Natural merge sort sacrifices at most \$2N\$ amount of work in order to recognize the runs in the input range. A run is a contiguous subsequence which is ascending or strictly descending. Every descending run is reversed. We don't allow non-strictly descending runs, since reversing them would rearrange possible equal elements in it and make the entire algorithm non-stable.
After it has figured the runs, it puts them all into a queue. Then while there is runs to merge, it pops two runs, merges them, and appends the new merged run to the tail of the queue. It continues in this manner until there is only one run in the queue left, which corresponds to the sorted range.
So, all in all, the time complexity of this sort is \$\Omega(N) \cap \mathcal{O}(N \log N)\$. \$\Theta(N)\$ space complexity.
Is that a good way of writing (potentially) reusable code? Are there some C++-idioms I should have adhered to?
natural_merge_sort.h:
#ifndef NATURAL_MERGE_SORT_H
#define NATURAL_MERGE_SORT_H
#include <algorithm>
#include <iostream>
#include <iterator>
/*******************************************************************************
* Implements a simple, array-based queue of integers. All three operations run *
* in constant time. This queue, however, does not check for under-/overflow of *
* underlying buffer because of performance considerations. *
*******************************************************************************/
class UnsafeIntQueue {
private:
const size_t MINIMUM_CAPACITY = 256;
size_t m_head;
size_t m_tail;
size_t m_size;
size_t m_mask;
size_t* m_buffer;
/***************************************************************************
* Makes sure a capacity is at least 'MINIMUM_CAPACITY' and is a power of *
* two. *
***************************************************************************/
size_t fixCapacity(size_t capacity)
{
capacity = std::max(capacity, MINIMUM_CAPACITY);
size_t s = 1;
while (s < capacity)
{
s <<= 1;
}
return s;
}
public:
/***************************************************************************
* Constructs a new integer queue, which can accommodate 'capacit' amount *
* integers. *
***************************************************************************/
UnsafeIntQueue(size_t capacity) :
m_head{0},
m_tail{0},
m_size{0}
{
capacity = fixCapacity(capacity);
m_mask = capacity - 1;
m_buffer = new size_t[capacity];
}
/***************************************************************************
* Destroys this queue, which releases the underlying buffer. *
***************************************************************************/
~UnsafeIntQueue()
{
delete[] m_buffer;
}
/***************************************************************************
* Appends the input integer to the tail of this queue. *
***************************************************************************/
void enqueue(const size_t element)
{
m_buffer[m_tail & m_mask] = element;
m_tail = (m_tail + 1) & m_mask;
m_size++;
}
/***************************************************************************
* Removes and returns the integer at the head of this queue. *
***************************************************************************/
size_t dequeue()
{
const size_t ret = m_buffer[m_head];
m_head = (m_head + 1) & m_mask;
m_size--;
return ret;
}
/***************************************************************************
* Returns the amount of integers in this queue. *
***************************************************************************/
size_t size() const
{
return m_size;
}
};
/*******************************************************************************
* Scans the range [first, lst) and returns the queue containing sizes of each *
* run in the order they appear while scanning from left to right. *
*******************************************************************************/
template<class RandomIt, class Cmp>
std::unique_ptr<UnsafeIntQueue> build_run_size_queue(RandomIt first,
RandomIt lst,
Cmp cmp)
{
const size_t length = std::distance(first, lst);
UnsafeIntQueue* p_q = new UnsafeIntQueue(length / 2 + 1);
RandomIt head;
RandomIt left = first;
RandomIt right = left + 1;
const RandomIt last = lst - 1;
while (left < last)
{
head = left;
if (cmp(*right++, *left++))
{
// Reading a strictly descending run.
while (left < last && cmp(*right, *left))
{
++left;
++right;
}
p_q->enqueue(right - head);
std::reverse(head, right);
}
else
{
// Reading a ascending run.
while (left < last && !cmp(*right, *left))
{
++left;
++right;
}
p_q->enqueue(left - head + 1);
}
++left;
++right;
}
if (left == last)
{
// Handle the case of an orphan element at the end of the range.
p_q->enqueue(1);
}
return std::unique_ptr<UnsafeIntQueue>(p_q);
}
/*******************************************************************************
* Returns the amount of leading zeros in 'num'. *
*******************************************************************************/
size_t leading_zeros(const size_t num)
{
size_t count = 0;
for (size_t t = (size_t) 1 << (8 * sizeof(t) - 1); t; t >>= 1, ++count)
{
if ((t & num))
{
return count;
}
}
return count;
}
/*******************************************************************************
* Returns the amount of merge passes needed to sort a range with 'run_amount' *
* runs. *
*******************************************************************************/
size_t get_pass_amount(size_t run_amount)
{
return 8 * sizeof(run_amount) - leading_zeros(run_amount - 1);
}
/*******************************************************************************
* Implements the merging routine. Runs in in time O(n + m), where 'n' is the *
* the length of the left chunk and 'm' is the length of the right chunk. *
*******************************************************************************/
template<class RandomIt, class Cmp>
void merge(RandomIt source,
RandomIt target,
const size_t offset,
const size_t left_run_length,
const size_t right_run_length,
Cmp cmp)
{
RandomIt left = source + offset;
RandomIt right = left + left_run_length;
const RandomIt left_bound = right;
const RandomIt right_bound = right + right_run_length;
RandomIt target_iter = target + offset;
while (left < left_bound && right < right_bound)
{
*target_iter++ = cmp(*right, *left) ?
*right++:
*left++;
}
std::copy(left, left_bound, target_iter);
std::copy(right, right_bound, target_iter);
}
/*******************************************************************************
* Implements the natural merge sort, which sacrifices one pass over the input *
* range in order to establish an implicit queue of runs. A run is the longest *
* consecutive subsequence, in which all elements are ascending or strictly *
* descending. Every descending run is reversed to ascending run. We cannot *
* consider non-strictly descending runs, since that would sacrifice the stabi- *
* lity of the algorithm. After the run queue is establish, the algorithm re- *
* moves two runs from the head of the queue, merges them into one run, which *
* is then appended to the tail of the run queue. Merging continues until the *
* queue contains only one run, which denotes that the entire input range is *
* sorted. *
* *
* The best-case complexity is O(N), the average and worst-case complexity is *
* O(N log N). Space complexity is O(N). *
*******************************************************************************/
template<class RandomIt, class Cmp>
void natural_merge_sort(RandomIt first, RandomIt last, Cmp cmp)
{
const size_t length = std::distance(first, last);
if (length < 2)
{
// Trivially sorted.
return;
}
typedef typename std::iterator_traits<RandomIt>::value_type value_type;
// Scan the runs.
std::unique_ptr<UnsafeIntQueue> p_queue =
build_run_size_queue(first, last, cmp);
// Request a buffer.
RandomIt buffer = new value_type[length];
std::copy(first, last, buffer);
// Count the amount of merge passes over the array required to bring order.
const size_t merge_passes = get_pass_amount(p_queue->size());
RandomIt source;
RandomIt target;
// Make sure that after the last merge pass, all data ends up in the input
// container.
if ((merge_passes & 1) == 1)
{
source = buffer;
target = first;
}
else
{
source = first;
target = buffer;
}
size_t runs_left = p_queue->size();
size_t offset = 0;
// While there is runs to merge, do...
while (p_queue->size() > 1)
{
// Remove two runs from the head of the run queue.
const size_t left_run_length = p_queue->dequeue();
const size_t right_run_length = p_queue->dequeue();
// Do the merge.
merge(source,
target,
offset,
left_run_length,
right_run_length,
cmp);
// Append the merged run to the tail of the queue.
p_queue->enqueue(left_run_length + right_run_length);
runs_left -= 2;
offset += left_run_length + right_run_length;
// The current pass over the array is almost complete.
switch (runs_left)
{
case 1:
{
const size_t single_length = p_queue->dequeue();
std::copy(source + offset,
source + offset + single_length,
target + offset);
p_queue->enqueue(single_length);
}
// FALL THROUGH!
case 0:
{
runs_left = p_queue->size();
offset = 0;
RandomIt tmp = source;
source = target;
target = tmp;
break;
}
}
}
delete[] buffer;
}
#endif
main.cpp:
#include <algorithm>
#include <chrono>
#include <iostream>
#include <random>
#include "natural_merge_sort.h"
/*******************************************************************************
* Prints an array. *
*******************************************************************************/
template<class T>
static void print_int_array(const T* begin, const T* last)
{
while (begin < last)
{
std::cout << *begin++ << " ";
}
std::cout << std::endl;
}
/*******************************************************************************
* Checks that the input range is sorted (is in ascending order). *
*******************************************************************************/
template<class T, class Cmp>
bool is_sorted(T* begin, T* end, Cmp cmp)
{
while (begin < end - 1)
{
if (cmp(*(begin + 1), *begin))
{
return false;
}
++begin;
}
return true;
}
/*******************************************************************************
* Creates a random integer array of length 'length', minimum integer *
* 'minimum', maximum integer 'maximum', using seed 'seed'. *
*******************************************************************************/
static int* get_random_int_array(const size_t length,
const int minimum,
const int maximum,
const unsigned int seed)
{
std::default_random_engine generator(seed);
std::uniform_int_distribution<int> distribution(minimum, maximum);
int* array = new int[length];
for (size_t i = 0; i < length; ++i)
{
array[i] = distribution(generator);
}
return array;
}
/*******************************************************************************
* Create an array of pointers to integers. *
*******************************************************************************/
static int** get_random_int_pointer_array(const size_t length,
const int minimum,
const int maximum,
const unsigned seed)
{
std::default_random_engine generator(seed);
std::uniform_int_distribution<int> distribution(minimum, maximum);
int** array = new int*[length];
for (size_t i = 0; i < length; ++i)
{
array[i] = new int(distribution(generator));
}
return array;
}
/*******************************************************************************
* Returns a strongly presorted array of integers. *
*******************************************************************************/
static int* get_presorted_int_array(const size_t length)
{
int* array = new int[length];
int num = 0;
for (size_t i = 0; i < length / 2; ++i)
{
array[i] = num++;
}
for (size_t i = length / 2; i < length; ++i)
{
array[i] = num--;
}
return array;
}
/*******************************************************************************
* Returns the milliseconds since the Unix epoch. *
*******************************************************************************/
static unsigned long long get_milliseconds()
{
return std::chrono::duration_cast<std::chrono::milliseconds>(
std::chrono::system_clock::now().time_since_epoch()).count();
}
/*******************************************************************************
* Checks that the two ranges are of the same length and content. *
*******************************************************************************/
template <class T>
bool are_equal(const T* begin1,
const T* end1,
const T* begin2,
const T* end2)
{
if (std::distance(begin1, end1) != std::distance(begin2, end2))
{
return false;
}
while (begin1 < end1)
{
if (*begin1++ != *begin2++)
{
return false;
}
}
return true;
}
/*******************************************************************************
* Compares two integer pointers by the values they point to. *
*******************************************************************************/
static bool compare_dereference(const int* a, const int* b)
{
return *a < *b;
}
/*******************************************************************************
* Compares two integers. *
*******************************************************************************/
static bool compare_int(const int a, const int b)
{
return a < b;
}
/*******************************************************************************
* Profiles the 'std::stable_sort' agains the range ['begin', 'end') using the *
* comparator 'cmp'. *
*******************************************************************************/
template<class T, class Cmp>
static void profile_stable_sort(T* begin, T* end, Cmp cmp)
{
unsigned long long ta = get_milliseconds();
std::stable_sort(begin, end, cmp);
unsigned long long tb = get_milliseconds();
std::cout << "std::stable_sort in "
<< (tb - ta)
<< " milliseconds. "
<< "Sorted: "
<< is_sorted(begin, end, cmp)
<< std::endl;
}
/*******************************************************************************
* Profiles the 'natural_merge_sort' agains the range ['begin', 'end') using *
* the comparator 'cmp'. *
*******************************************************************************/
template<class T, class Cmp>
static void profile_natural_merge_sort(T* begin, T* end, Cmp cmp)
{
unsigned long long ta = get_milliseconds();
natural_merge_sort(begin, end, cmp);
unsigned long long tb = get_milliseconds();
std::cout << "natural_merge_sort in "
<< (tb - ta)
<< " milliseconds. "
<< "Sorted: "
<< is_sorted(begin, end, cmp)
<< std::endl;
}
/*******************************************************************************
* Profiles the sorting algorithms on a random integer array. *
*******************************************************************************/
static void profile_on_random_array(const size_t sz,
const int minimum,
const int maximum,
const unsigned seed)
{
int* array1 = get_random_int_array(sz, minimum, maximum, seed);
int* array2 = new int[sz];
std::copy(array1, array1 + sz, array2);
std::cout << "--- PROFILING ON RANDOM ARRAY OF LENGTH "
<< sz
<< " ---"
<< std::endl;
profile_stable_sort(array1, array1 + sz, compare_int);
profile_natural_merge_sort(array2, array2 + sz, compare_int);
std::cout << "Same contents: "
<< are_equal(array1, array1 + sz, array2, array2 + sz)
<< std::endl
<< std::endl;
}
/*******************************************************************************
* Profiles the sorting algorithms on an array of pointers to random integers. *
*******************************************************************************/
static void profile_on_integer_pointer_array(const size_t sz,
const int minimum,
const int maximum,
const unsigned seed)
{
std::cout << "--- PROFILING ON RANDOM POINTER ARRAY OF LENGTH "
<< sz
<< " ---"
<< std::endl;
int** array1 = get_random_int_pointer_array(sz,
minimum,
maximum,
seed);
int** array2 = new int*[sz];
std::copy(array1, array1 + sz, array2);
profile_stable_sort(array1, array1 + sz, compare_dereference);
profile_natural_merge_sort(array2, array2 + sz, compare_dereference);
std::cout << "Same contents: "
<< are_equal(array1, array1 + sz, array2, array2 + sz)
<< std::endl
<< std::endl;
}
/*******************************************************************************
* Profiles the sorting algorithms on a presorted array. *
*******************************************************************************/
static void profile_on_presorted_array(const size_t sz)
{
std::cout << "--- PROFILING ON PRESORTED ARRAY OF LENGTH "
<< sz
<< " ---"
<< std::endl;
int* array1 = get_presorted_int_array(sz);
int* array2 = new int[sz];
std::copy(array1, array1 + sz, array2);
profile_stable_sort(array1, array1 + sz, compare_int);
profile_natural_merge_sort(array2, array2 + sz, compare_int);
std::cout << "Same contents: "
<< are_equal(array1, array1 + sz, array2, array2 + sz)
<< std::endl
<< std::endl;
}
/*******************************************************************************
* The entry point to a demo program. *
*******************************************************************************/
int main(int argc, const char * argv[]) {
unsigned long long seed = get_milliseconds();
////
std::cout << "Seed: "
<< seed
<< std::endl
<< std::endl;
////
const size_t length = 5000000;
const int min_int = -100;
const int max_int = 300;
std::cout << std::boolalpha;
profile_on_random_array(length, min_int, max_int, seed);
profile_on_integer_pointer_array(length, min_int, max_int, seed);
profile_on_presorted_array(length);
return 0;
}