Is partial_sort_copy the fastest C++ partial sort?

Question

Consider the following function, median:

  real_t median(const std::initializer_list<real_t> vars) {
    real_t tmp[15];
    const unsigned x = vars.size() / 2;
    if (x & 1) {
      std::partial_sort_copy(vars.begin(), vars.end(), &tmp[0], &tmp[x]);
      return tmp[x];
    }
    const unsigned y = x + 1;
    std::partial_sort_copy(vars.begin(), vars.end(), &tmp[0], &tmp[y]);
    return (tmp[x] + tmp[y]) / 2;
  }

I am using a partial sort to decrease complexity, as I only need to sort half of the list.

Further, I have assumed that std::partial_sort_copy is faster than std::partial_sort or std::nth_element because there is no shuffling required in the sort algorithm (It1 != It2). Is my assumption correct?

---

NB: Assume real_t could be a double, so please don't criticise the use of division.

NBB: I'm using -pedantic and vars is known to not be longer than 15 elements.

Answer 1

Using the following code

#include <chrono>
#include <iostream>
#include <thread>
#include <string>
#include <array>
#include <algorithm>

volatile int answer;

const int size = 15;

std::array<std::array<int, size>, 0x100> fresh_data;
std::array<std::array<int, size>, 0x100> data;

void naive(int n) {
    auto & a = data[n];
    std::sort(a.begin(), a.end());
    answer = a[size / 2];
}

void fancy(int n) {
    auto & a = data[n];
    std::partial_sort(a.begin(), a.begin() + (size / 2 + 1), a.end());
    answer = a[size / 2 ];
}

void ghoul(int n) {
    auto & a = data[n];
    std::array<int, size / 2 + 1> temp;
    std::partial_sort_copy(a.begin(), a.end(), temp.begin(), temp.end());
    answer = temp[size / 2];

}

void nthel(int n) {
    auto & a = data[n];
    std::nth_element(a.begin(), a.begin() + size / 2, a.end());
    answer = a[size / 2];
}

void gen_data() {
    for (auto & a : fresh_data)
    for (auto & b : a)
        b = rand();
}

void regen_data() {
    data = fresh_data;
}


template <typename T>
void test(T f, std::string n) {
    regen_data();
    auto a = std::chrono::high_resolution_clock::now();
    for (auto i = 0; i < 10000; ++i)
    for (auto i = 0; i < 0x100; ++i)
        f(i);
    auto b = std::chrono::high_resolution_clock::now();
    std::cout << n << ": " << std::chrono::duration_cast<std::chrono::milliseconds>(b - a).count() << std::endl;
}

int main() {
    gen_data();
    test(naive, "             std::sort");
    test(fancy, "     std::partial_sort");
    test(ghoul, "std::partial_sort_copy");
    test(nthel, "      std::nth_element");
}

I get the following results:

             std::sort: 141
     std::partial_sort: 359
std::partial_sort_copy: 831
      std::nth_element: 149

Tested using Visual Studio 2013 in 64bit release mode on an AMD Phenom II x4 2.5GHz.

Answer 2

If I could choose, I would go for a Partial Quicksort.

Info on Partial Quicksort

But if you have to compare these two only... then Partial sort is better vs partial sort copy. Here you have more info about these two methods:

Info on Partial Sort

Info on Partial Sort Copy

Here you also find an algorithm code example for Partial Quicksort - it was implemented in C and matlab:

Example - Partial Quicksort

Answer 3

Have you tested your code?

std::partial_sort_copy(vars.begin(), vars.end(), &tmp[0], &tmp[x]); will not copy anything to tmp[x] since &tmp[x] is considered the end of the half-open range (ie., it's just past the last valid element). So your return statements access indeterminate or default constructed array elements.

Try the following:

real_t median(const std::initializer_list<real_t> vars) 
{
    real_t tmp[15];

    size_t siz = vars.size();
    if ((siz == 0) || (15 < siz)) return 0;     // or throw some sort of exception or ???

    const unsigned x = vars.size() / 2;
    std::partial_sort_copy(vars.begin(), vars.end(), &tmp[0], &tmp[x+1]);

    if (siz % 2 == 0) {
        return (tmp[x-1] + tmp[x]) / 2;
    }

    return tmp[x];
}

Note that if given an initializer_list as the source of data, a modify-in-place algorithm such as nth_element or partial_sort won't work since initializer lists can't be modified (whether the parameter is marked const or not - the iterators into the initializer_list are const qualified). So a copy must be done to find the median using the standard algorithm functions, either by copying the list before calling the algorithm, or by using an algorithm variant that performs the copy as part of its work, such as partial_sort_copy().