Tpetra_Details_shortSort.hpp

Go to the documentation of this file.

// @HEADER
// *****************************************************************************
//          Tpetra: Templated Linear Algebra Services Package
//
// Copyright 2008 NTESS and the Tpetra contributors.
// SPDX-License-Identifier: BSD-3-Clause
// *****************************************************************************
// @HEADER

#ifndef TPETRA_DETAILS_SHORTSORT_HPP
#define TPETRA_DETAILS_SHORTSORT_HPP


#include "TpetraCore_config.h"
#include "Kokkos_Macros.hpp"
#include <type_traits>

namespace Tpetra {
namespace Details {

// Make sure that the macro defined below wasn't defined somewhere else.
#ifdef TPETRA_DETAILS_SWAP_KEYSANDVALUES
#error "The TPETRA_DETAILS_SWAP_KEYSANDVALUES macro is already defined."
#endif  // TPETRA_DETAILS_SWAP_KEYSANDVALUES

#define TPETRA_DETAILS_SWAP_KEYSANDVALUES(i, j) \
  if (keys[i] > keys[j]) {                      \
    const KeyType tmpKey(keys[i]);              \
    keys[i] = keys[j];                          \
    keys[j] = tmpKey;                           \
    const ValueType tmpVal(values[i]);          \
    values[i] = values[j];                      \
    values[j] = tmpVal;                         \
  }

// Make sure that the macro defined below wasn't defined somewhere else.
#ifdef TPETRA_DETAILS_SWAP_KEYS
#error "The TPETRA_DETAILS_SWAP_KEYS macro is already defined."
#endif  // TPETRA_DETAILS_SWAP_KEYSANDVALUES

#define TPETRA_DETAILS_SWAP_KEYS(i, j) \
  if (keys[i] > keys[j]) {             \
    const KeyType tmpKey(keys[i]);     \
    keys[i] = keys[j];                 \
    keys[j] = tmpKey;                  \
  }

template <class KeyType, class ValueType>
KOKKOS_FUNCTION void
shortSortKeysAndValues_2(KeyType keys[2], ValueType values[2]) {
  // Since this function takes a constant number of entries, I use a
  // sorting network here.  For 2 entries, the sorting network is
  // nearly trivial.
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(0, 1);
}

template <class KeyType>
KOKKOS_FUNCTION void
shortSortKeys_2(KeyType keys[2]) {
  // Since this function takes a constant number of entries, I use a
  // sorting network here.  For 2 entries, the sorting network is
  // nearly trivial.
  TPETRA_DETAILS_SWAP_KEYS(0, 1);
}

template <class KeyType, class ValueType>
KOKKOS_FUNCTION void
shortSortKeysAndValues_3(KeyType keys[3], ValueType values[3]) {
  // Since this function takes a constant number of entries, I use a
  // sorting network here.  To make the network, I used the generator
  // at
  //
  // http://pages.ripco.net/~jgamble/nw.html
  //
  // "Best" algorithm in this case is the Bose-Nelson Algorithm.
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(1, 2);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(0, 2);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(0, 1);
}

template <class KeyType>
KOKKOS_FUNCTION void
shortSortKeys_3(KeyType keys[3]) {
  // Since this function takes a constant number of entries, I use a
  // sorting network here.  To make the network, I used the generator
  // at
  //
  // http://pages.ripco.net/~jgamble/nw.html
  //
  // "Best" algorithm in this case is the Bose-Nelson Algorithm.
  TPETRA_DETAILS_SWAP_KEYS(1, 2);
  TPETRA_DETAILS_SWAP_KEYS(0, 2);
  TPETRA_DETAILS_SWAP_KEYS(0, 1);
}

template <class KeyType, class ValueType>
KOKKOS_FUNCTION void
shortSortKeysAndValues_4(KeyType keys[4], ValueType values[4]) {
  // Since this function takes a constant number of entries, I use a
  // sorting network here.  To make the network, I used the generator
  // at
  //
  // http://pages.ripco.net/~jgamble/nw.html
  //
  // "Best" algorithm in this case is the Bose-Nelson Algorithm.
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(0, 1);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(2, 3);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(0, 2);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(1, 3);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(1, 2);
}

template <class KeyType>
KOKKOS_FUNCTION void
shortSortKeys_4(KeyType keys[4]) {
  // Since this function takes a constant number of entries, I use a
  // sorting network here.  To make the network, I used the generator
  // at
  //
  // http://pages.ripco.net/~jgamble/nw.html
  //
  // "Best" algorithm in this case is the Bose-Nelson Algorithm.
  TPETRA_DETAILS_SWAP_KEYS(0, 1);
  TPETRA_DETAILS_SWAP_KEYS(2, 3);
  TPETRA_DETAILS_SWAP_KEYS(0, 2);
  TPETRA_DETAILS_SWAP_KEYS(1, 3);
  TPETRA_DETAILS_SWAP_KEYS(1, 2);
}

template <class KeyType, class ValueType>
KOKKOS_FUNCTION void
shortSortKeysAndValues_8(KeyType keys[8], ValueType values[8]) {
  // Since this function takes a constant number of entries, I use a
  // sorting network here.  To make the network, I used the generator
  // at
  //
  // http://pages.ripco.net/~jgamble/nw.html
  //
  // "Best" algorithm in this case is the Bose-Nelson Algorithm.
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(0, 1);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(2, 3);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(0, 2);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(1, 3);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(1, 2);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(4, 5);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(6, 7);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(4, 6);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(5, 7);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(5, 6);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(0, 4);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(1, 5);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(1, 4);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(2, 6);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(3, 7);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(3, 6);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(2, 4);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(3, 5);
  TPETRA_DETAILS_SWAP_KEYSANDVALUES(3, 4);
}

template <class KeyType>
KOKKOS_FUNCTION void
shortSortKeys_8(KeyType keys[8]) {
  // Since this function takes a constant number of entries, I use a
  // sorting network here.  To make the network, I used the generator
  // at
  //
  // http://pages.ripco.net/~jgamble/nw.html
  //
  // "Best" algorithm in this case is the Bose-Nelson Algorithm.
  TPETRA_DETAILS_SWAP_KEYS(0, 1);
  TPETRA_DETAILS_SWAP_KEYS(2, 3);
  TPETRA_DETAILS_SWAP_KEYS(0, 2);
  TPETRA_DETAILS_SWAP_KEYS(1, 3);
  TPETRA_DETAILS_SWAP_KEYS(1, 2);
  TPETRA_DETAILS_SWAP_KEYS(4, 5);
  TPETRA_DETAILS_SWAP_KEYS(6, 7);
  TPETRA_DETAILS_SWAP_KEYS(4, 6);
  TPETRA_DETAILS_SWAP_KEYS(5, 7);
  TPETRA_DETAILS_SWAP_KEYS(5, 6);
  TPETRA_DETAILS_SWAP_KEYS(0, 4);
  TPETRA_DETAILS_SWAP_KEYS(1, 5);
  TPETRA_DETAILS_SWAP_KEYS(1, 4);
  TPETRA_DETAILS_SWAP_KEYS(2, 6);
  TPETRA_DETAILS_SWAP_KEYS(3, 7);
  TPETRA_DETAILS_SWAP_KEYS(3, 6);
  TPETRA_DETAILS_SWAP_KEYS(2, 4);
  TPETRA_DETAILS_SWAP_KEYS(3, 5);
  TPETRA_DETAILS_SWAP_KEYS(3, 4);
}

template <class KeyType, class ValueType, class IndexType>
KOKKOS_FUNCTION void
shellSortKeysAndValues(KeyType keys[],
                       ValueType values[],
                       const IndexType n) {
  static_assert(std::is_integral<IndexType>::value,
                "IndexType must be a signed integer type.");
  static_assert(std::is_signed<IndexType>::value,
                "IndexType must be a signed integer type.  "
                "This implementation does a count-down loop, "
                "and may thus loop forever "
                "if one attempts to use it with unsigned IndexType.");
  constexpr IndexType ZERO = 0;
  IndexType midpoint       = n / static_cast<IndexType>(2);

  while (midpoint > ZERO) {
    // Avoid names like "max" in case they collide with macros.
    const IndexType theMax = n - midpoint;
    for (IndexType j = 0; j < theMax; ++j) {
      // IndexType is signed, so it's legit to compare >= 0.
      for (IndexType k = j; k >= 0; k -= midpoint) {
        if (keys[k + midpoint] >= keys[k]) {
          break;
        }
        const KeyType tmpKey   = keys[k + midpoint];
        keys[k + midpoint]     = keys[k];
        keys[k]                = tmpKey;
        const ValueType tmpVal = values[k + midpoint];
        values[k + midpoint]   = values[k];
        values[k]              = tmpVal;
      }
    }
    midpoint = midpoint / 2;
  }
}

template <class KeyType, class IndexType>
KOKKOS_FUNCTION void
shellSortKeys(KeyType keys[], const IndexType n) {
  static_assert(std::is_integral<IndexType>::value,
                "IndexType must be a signed integer type.");
  static_assert(std::is_signed<IndexType>::value,
                "IndexType must be a signed integer type.  "
                "This implementation does a count-down loop, "
                "and may thus loop forever "
                "if one attempts to use it with unsigned IndexType.");
  constexpr IndexType ZERO = 0;
  IndexType midpoint       = n / static_cast<IndexType>(2);

  while (midpoint > ZERO) {
    // Avoid names like "max" in case they collide with macros.
    const IndexType theMax = n - midpoint;
    for (int j = 0; j < theMax; ++j) {
      // IndexType must be signed, so it's legit to compare >= 0.
      for (int k = j; k >= 0; k -= midpoint) {
        if (keys[k + midpoint] >= keys[k]) {
          break;
        }
        const KeyType tmpKey = keys[k + midpoint];
        keys[k + midpoint]   = keys[k];
        keys[k]              = tmpKey;
      }
    }
    midpoint = midpoint / 2;
  }
}

template <typename KeyType, typename ValueType, typename IndexType>
KOKKOS_FUNCTION void
shellSortKeysAndValues2(KeyType* keys, ValueType* values, IndexType n) {
  IndexType ngaps = 10;
  // Use Ciura's gap choices
  IndexType gaps[] = {3548, 1577, 701, 301, 132, 57, 23, 10, 4, 1};
  for (IndexType gapIndex = 0; gapIndex < ngaps; gapIndex++) {
    auto gap = gaps[gapIndex];
    if (n < gap)
      continue;
    // insertion sort the array {keys[0*gap], keys[1*gap], keys[2*gap], ...}
    for (IndexType gapOffset = 0; gapOffset < gap; gapOffset++) {
      for (IndexType i = gap + gapOffset; i < n; i += gap) {
        // avoid extra swaps: scan for the final position of keys[i]
        if (keys[i - gap] > keys[i]) {
          IndexType finalPos = i - gap;
          while (finalPos - gap >= 0 && keys[finalPos - gap] > keys[i]) {
            finalPos -= gap;
          }
          // save keys/values [i], then shift up all keys/values between finalPos and i-gap (inclusive)
          auto tempKey = keys[i];
          auto tempVal = values[i];
          for (IndexType j = i - gap; j >= finalPos; j -= gap) {
            keys[j + gap]   = keys[j];
            values[j + gap] = values[j];
          }
          keys[finalPos]   = tempKey;
          values[finalPos] = tempVal;
        }
      }
    }
  }
#undef SHELL_SWAP
}

}  // namespace Details
}  // namespace Tpetra

#endif  // TPETRA_DETAILS_SHORTSORT_HPP