Binary heap performance optimization

Question

I have implemented a binary heap as a exercise in data structures. I realize that I'm missing some useful functions (increase/decrease key for example) but I don't consider them interesting for the exercise.

I'm interested in learning if there is anything that I could do to improve performance of the pop operation which by my benchmarks is the slowest operation. I have already taken it down from the initial 6.5s to 2.5s on my machine.

Compiling and running the below on my machine yields:

g++ -O3 --std=c++11 -Wall -Wextra -pedantic test_heap.cc

/a.exe

Push time: 177ms Pop time: 2504ms Sum: -2014260032

I have used gprof to profile but the compiler inlines everything into the main function so gprof unfortunately doesn't give anything useful in this case. And I don't feel like the results are useful if disabling optimizations to measure performance...

Note that I'm printing the sum only to avoid the compiler from removing the benchmark code because it isn't used.

Here is the implementation heap.hpp:

#include <vector>
#include <functional> // std::less
#include <utility> // std::swap
#include <stdexcept>

namespace elads{

  // A binary min-heap implemented using the given comparator.
  template<typename T, typename Comparator = std::less<T> >
  class heap{
  public:
    void push(const T& value){
      auto idx = m_data.size();
      m_data.emplace_back(value);
      sift_up(idx);
    }

    void pop(){
      if(m_data.empty()){
        throw std::runtime_error("Stack is empty!");
      }

      if(m_data.size() == 1){
        m_data.clear();
        return;
      }

      swap(0, m_data.size() - 1);
      m_data.erase(m_data.end()-1);
      sift_down(0);
    }

    const T& top() const {
      if(m_data.empty()){
        throw std::runtime_error("Stack is empty!");
      }
      return m_data.front();
    }

    bool empty() const{
      return m_data.empty();
    }

  private:
    std::vector<T> m_data;
    Comparator m_cmp;

    bool less(size_t lhs, size_t rhs) const{
      return m_cmp(m_data[lhs], m_data[rhs]);
    }

    void swap(size_t lhs, size_t rhs){
      using std::swap;
      swap(m_data[lhs], m_data[rhs]);
    }

    void sift_up(size_t child){
      if(root(child))
        return;
      auto parent = parent_of(child);
      if(!less(parent, child)){
        swap(parent, child);
        sift_up(parent);
      }
    }

    void sift_down(size_t parent){
      auto left = left_of(parent);
      auto right = right_of(parent);

      // Note: "right <= parent" means that the index calculation of the right child overflowed.
      // Assume: 0 <= parent < 2^k  (i.e. parent is valid index) for some positive k depending on the index type.
      // The index of right is calculated as: 
      //
      //     right = 2*parent + 2 (0)
      //
      // Then an overflow that the above test would return false for can only occurr if:
      //     
      //     (2*parent + 2 =) (parent + x) % (2^k) > parent  (1)
      //
      // Where x = parent + 2. Inequality (1) can only be true if:
      // 
      //     x = n*2^k + y; where "n >= 1" and "0 < y < 2^k - parent" (2).
      //
      // Now we know that "parent < 2^k" thus "n >= 1" can only occur if "2^k - 2 <= parent < 2^k"
      // so we can test the two cases:
      // 
      // Case: parent = 2^k - 2
      // From the above we know that n=1. Inserting "parent = 2^k - 2"  into (2) yields: 0 < y < 2 
      // and calculating y yields y = 0, thus this case contradicts (1).
      // 
      // Case: parent = 2^k - 1
      // From the above we know that n=1. Inserting "parent = 2^k - 1" into (2) yields: 0 < y < 1
      // and calculating y yields y = 1. This case too contradicts (1). 
      //
      // Thus we have proven that if "right <= parent" is sufficienct and necessary conditions for
      // detecting if an overflow has occurred. For the left index replace "2*parent + 2" with "2*parent + 1"
      // and repeat the same steps. The result will be that the first case cannot occurr.

      auto has_right = right < m_data.size() && right > parent; // right > parent == !(right <= parent) - overflow
      auto has_left = left < m_data.size() && left > parent;

      if(has_right && less(right, parent)){
        // By design, left < right so if has_right is true then has_left is also true.
        if(less(left, right)){
          swap(left, parent);
          sift_down(left);
        }
        else{
          swap(right, parent);
          sift_down(right);
        }
      }
      else if(has_left && less(left, parent)){
        swap(left, parent);
        sift_down(left);
      }
      // Else, well we're either on the last level or the heap property holds
    }

    size_t parent_of(size_t index) const{
        return (index - 1)/2;
      }

    size_t left_of(size_t index) const{
        return 2*index + 1;
      }

    size_t right_of(size_t index) const{
        return 2*index + 2;
      }

    bool root(size_t index) const{
        return index == 0;
    } 
  };

}

Here is the benchmark/test harness test_heap.cpp:

#include "heap.hpp"

#include <vector>
#include <numeric>
#include <algorithm>
#include <random>
#include <chrono>
#include <ctime>
#include <cassert>
#include <iostream>

template<typename TimeT = std::chrono::milliseconds>
struct measure
{
     template<typename F, typename ...Args>
     static typename TimeT::rep execution(
     F func, Args&&... args)
     {
         auto start = std::chrono::system_clock::now();
         func(std::forward<Args>(args)...);
         auto duration = std::chrono::duration_cast<
            TimeT>(
             std::chrono::system_clock::now() - start);
         return duration.count();
     }
};

int main(){
  // Setup test data
  const auto n = 10000000;
  auto rng = std::mt19937(std::clock());
  std::vector<int> dataset(n);
  std::iota(dataset.begin(), dataset.end(), 0);
  std::shuffle(dataset.begin(), dataset.end(), rng);

  // Begin test
  elads::heap<int> heap;

  auto push_time_ms = measure<>::execution([&](){
      for(auto v : dataset){
        heap.push(v);
      }
    });

  int sum = 0;
  auto pop_time_ms = measure<>::execution([&](){
      int expected = 0;
      while(!heap.empty()){
        assert(expected++ == heap.top());
        sum += heap.top();
        heap.pop();
      }
    });

  std::cout<<"Push time: "<<push_time_ms<<"ms Pop time: "<<pop_time_ms<<"ms Sum: "<<sum<<std::endl;

  return 0;
}

Answer 1

Drop the checks

Right now, your pop() and top() both check for emptiness and throw. assert() would be better - you don't want these functions throwing. It's likely that your user will check for empty() first anyway - and if they don't, make sure they understand that they have to. With that, our top() becomes trivial:

const T& top() const {
    return m_data.front();
}

And our pop() is a few lines shorter.

Swapping

Two unrelated comments to make about swapping. First, in your pop(), you swap the first and last elements and then erase the last one. But you don't actually need to do that. We don't care that the last element is the first one when we erase it. Just move the last one into the first one:

m_data.front() = std::move(m_data.back());
m_data.pop_back();

That'll save some unnecessary operations.

The second comment is about your member function swap(). It's a little confusing in that it doesn't swap the arguments passed in - it swaps those indices. So prefer to name it swap_idx() or something along those lines just for clarity.

Special cases

In your pop(), you start with the special case that size() == 1. But there's nothing special about that case for the purpose of your algorithm. sift_down() will just do nothing. So you don't need it. That reduces your pop() down to:

void pop() {
    m_data.front() = std::move(m_data.back());
    m_data.pop_back();
    sift_down(0);
}

Overflow

You have this enormous comment block and logic in sift_down() describing overflow cases for indices. You are really not going to have a heap large enough to have potential for overflow. I wouldn't worry about it. It just makes your logic more complicated than it needs to be. Drop it.

Push it. Push it real good

If you have push(const T&), you may as well add push(T&& ), and then while we're at it, also add emplace(Args&&...). You could have them all forward to a common implementation:

void push(const T& value) { emplace(value); }
void push(T&& value) { emplace(std::move(value); }

template <class... Args>
void emplace(Args&&... args) {
    auto idx = m_data.size();
    m_data.emplace_back(std::forward<Args>(args)...);
    sift_up(idx);
}

Arbitrary Comparators

If you allow for a templated Comparator - you need to allow the user to pass one in via a constructor. As-is, it really only works with default-constructible types. Furthermore, I'd suggest inheriting from the comparator to take advantage of the empty base optimization - std::less<> for instance is empty, so can save that space.