MinMax.h

#ifndef DASH__ALGORITHM__MIN_MAX_H__
#define DASH__ALGORITHM__MIN_MAX_H__

#include <dash/internal/Config.h>

#include <dash/Allocator.h>

#include <dash/algorithm/LocalRange.h>

#include <dash/util/Config.h>
#include <dash/util/Trace.h>
#include <dash/util/UnitLocality.h>

#include <dash/internal/Logging.h>
#include <dash/iterator/GlobIter.h>

#include <algorithm>
#include <memory>

#ifdef DASH_ENABLE_OPENMP
#include <omp.h>
#endif


namespace dash {

template <
  class ElementType,
  class Compare = std::less<const ElementType &> >
const ElementType * min_element(
  const ElementType * l_range_begin,
  const ElementType * l_range_end,
  Compare             compare
    = std::less<const ElementType &>())
{
#ifdef DASH_ENABLE_OPENMP
#define ROUNDUP(T, A) ((T + A - 1) & ~(A - 1))

  typedef typename std::decay<ElementType>::type      value_t;
  dash::util::UnitLocality uloc;
  //TODO: ask dyloc for available memory spaces
  auto n_threads = uloc.num_domain_threads();
  DASH_LOG_DEBUG("dash::min_element", "thread capacity:",  n_threads);

  // TODO: Should also restrict on elements/units > ~10240.
  //       Find a model for the minimum work laod.
  if (n_threads > 1) {
    auto          l_size     = l_range_end - l_range_begin;
    int           min_idx_l  = 0;
    ElementType   min_val_l  = *l_range_begin;

    typedef struct min_pos_t { value_t val; size_t idx; } min_pos;

    DASH_LOG_TRACE("dash::min_element", "sizeof(min_pos):", sizeof(min_pos));

    dash::HostSpace hostSpace;

    int       align_bytes         = uloc.cache_line_size(0);
    int       single_element_sz   = ROUNDUP(sizeof(min_pos), align_bytes);

    size_t    min_vals_t_bytes    = single_element_sz * n_threads;

    auto min_vals_t_raw = static_cast<uint8_t *>(
        hostSpace.allocate(min_vals_t_bytes, align_bytes));

    DASH_LOG_TRACE("dash::min_element", "min * alloc:",   min_vals_t_raw);
    //DASH_LOG_TRACE("dash::min_element", "min * aligned:", min_vals_t);
    DASH_LOG_TRACE("dash::min_element", "min * size:",    min_vals_t_bytes);
    DASH_ASSERT_GE(min_vals_t_bytes, n_threads * sizeof(min_pos),
                   "Aligned buffer of min_pos has insufficient size");
    DASH_ASSERT_MSG(nullptr != min_vals_t_raw,
                    "Aligned allocation of min_pos returned nullptr");

    // Cannot use user-defined reduction (OpenMP 4.0) as the compare
    // predicate cannot be used in `omp declare reduction`.
    // Avoid omp for + omp critical section by using array of
    // thread-local minimum values, aligned to prevent false sharing:
    int t_id;
    #pragma omp parallel num_threads(n_threads) private(t_id)
    {
      // Documentation of Intel MIC intrinsics, see:
      // https://software.intel.com/de-de/node/523533
      // https://software.intel.com/de-de/node/523387
      t_id = omp_get_thread_num();
      DASH_LOG_TRACE("dash::min_element", "starting thread", t_id);
      auto & min_val_t = *(reinterpret_cast<min_pos *>(min_vals_t_raw + t_id * single_element_sz));
      min_val_t.idx = min_idx_l;
      min_val_t.val = min_val_l;
      // Cannot use explicit private(min_val_t) as ElementType might
      // not be default-constructible:
      #pragma omp for schedule(static)
      for (int i = 0; i < l_size; i++) {
        const ElementType & val_t = *(l_range_begin + i);
        if (compare(val_t, min_val_t.val)) {
          min_val_t.val = val_t;
          min_val_t.idx = i;
        }
      }
      DASH_LOG_TRACE("dash::min_element", "local minimum at thread", t_id,
                     "idx:", min_val_t.idx,
                     "val:", min_val_t.val);
    }

    min_pos min_pos_l = * (reinterpret_cast<min_pos *>(min_vals_t_raw));

    for (int t = 1; t < n_threads; t++) {
      const min_pos & mpt = *(reinterpret_cast<min_pos *>(min_vals_t_raw + t * single_element_sz));
      if (compare(mpt.val, min_pos_l.val)) {
        min_pos_l = mpt;
      }
    }

    hostSpace.deallocate(min_vals_t_raw, min_vals_t_bytes, align_bytes);

    return (l_range_begin + min_pos_l.idx);
  }
#endif // DASH_ENABLE_OPENMP
  return ::std::min_element(l_range_begin, l_range_end, compare);
}

template <
    typename GlobInputIt,
    class Compare = std::less<
        const typename dash::iterator_traits<GlobInputIt>::value_type &> >
GlobInputIt min_element(
    const typename std::enable_if<
        dash::iterator_traits<GlobInputIt>::is_global_iterator::value,
        GlobInputIt>::type &first,
    const GlobInputIt &last,
    Compare compare = Compare())
{
  typedef typename GlobInputIt::pattern_type     pattern_t;
  typedef typename pattern_t::index_type         index_t;
  typedef typename std::decay<
      typename dash::iterator_traits<GlobInputIt>::value_type>::type value_t;

  // return last for empty array
  if (first == last) {
    DASH_LOG_DEBUG("dash::min_element >",
                   "empty range, returning last", last);
    return last;
  }

  dash::util::Trace trace("min_element");

  auto & pattern = first.pattern();
  auto & team    = pattern.team();
  DASH_LOG_DEBUG("dash::min_element()",
                 "allocate minarr, size", team.size());
  // Global position of end element in range:
  auto    gi_last            = last.gpos();
  // Find the local min. element in parallel
  // Get local address range between global iterators:
  auto    local_idx_range    = dash::local_index_range(first, last);
  // Pointer to local minimum element:
  const   value_t * lmin = nullptr;
  // Local offset of local minimum element, or -1 if no element found:
  index_t l_idx_lmin         = -1;
  if (local_idx_range.begin == local_idx_range.end) {
    // local range is empty
    DASH_LOG_DEBUG("dash::min_element", "local range empty");
  } else {
    trace.enter_state("local");

    // Pointer to first element in local memory:
    auto *lbegin = dash::local_begin(
        static_cast<typename GlobInputIt::const_pointer>(first), team.myid());

    // Pointers to first / final element in local range:
    const auto * l_range_begin = lbegin + local_idx_range.begin;
    const auto * l_range_end   = lbegin + local_idx_range.end;

    lmin = dash::min_element(l_range_begin, l_range_end, compare);

    if (lmin != l_range_end) {
      DASH_LOG_TRACE_VAR("dash::min_element", *lmin);
      // Offset of local minimum in local memory:
      l_idx_lmin = lmin - lbegin;
    }

    trace.exit_state("local");
  }
  DASH_LOG_TRACE("dash::min_element",
                 "local index of local minimum:", l_idx_lmin);
  DASH_LOG_TRACE("dash::min_element",
                 "waiting for local min of other units");

  trace.enter_state("barrier");
  team.barrier();
  trace.exit_state("barrier");

  typedef struct {
    value_t  value;
    index_t  g_index;
  } local_min_t;

  std::vector<local_min_t> local_min_values(team.size());

  // Set global index of local minimum to -1 if no local minimum has been
  // found:
  local_min_t local_min;
  local_min.value   = l_idx_lmin < 0
                      ? value_t()
                      : *lmin;
  local_min.g_index = l_idx_lmin < 0
                      ? -1
                      : pattern.global(l_idx_lmin);

  DASH_LOG_TRACE("dash::min_element", "sending local minimum: {",
                 "value:",   local_min.value,
                 "g.index:", local_min.g_index, "}");

  DASH_LOG_TRACE("dash::min_element", "dart_allgather()");
  trace.enter_state("allgather");
  DASH_ASSERT_RETURNS(
    dart_allgather(
      &local_min,
      local_min_values.data(),
      sizeof(local_min_t),
      DART_TYPE_BYTE,
      team.dart_id()),
    DART_OK);
  trace.exit_state("allgather");

#ifdef DASH_ENABLE_LOGGING
  for (int lmin_u = 0; lmin_u < local_min_values.size(); lmin_u++) {
    auto lmin_entry = local_min_values[lmin_u];
    DASH_LOG_TRACE("dash::min_element", "dart_allgather >",
                   "unit:",    lmin_u,
                   "value:",   lmin_entry.value,
                   "g_index:", lmin_entry.g_index);
  }
#endif

  auto gmin_elem_it  = ::std::min_element(
                           local_min_values.begin(),
                           local_min_values.end(),
                           [&](const local_min_t & a,
                               const local_min_t & b) {
                             // Ignore elements with global index -1 (no
                             // element found):
                             return (b.g_index < 0 ||
                                     (a.g_index > 0 &&
                                      compare(a.value, b.value)));
                           });

  if (gmin_elem_it == local_min_values.end()) {
    DASH_LOG_DEBUG_VAR("dash::min_element >", last);
    return last;
  }

  auto gi_minimum    = gmin_elem_it->g_index;

  DASH_LOG_TRACE("dash::min_element",
                 "min. value:", gmin_elem_it->value,
                 "at unit:",    (gmin_elem_it - local_min_values.begin()),
                 "global idx:", gi_minimum);

  DASH_LOG_TRACE_VAR("dash::min_element", gi_minimum);
  if (gi_minimum < 0 || gi_minimum == gi_last) {
    DASH_LOG_DEBUG_VAR("dash::min_element >", last);
    return last;
  }
  // iterator 'first' is relative to start of input range, convert to start
  // of its referenced container (= container.begin()), then apply global
  // offset of minimum element:
  auto minimum = (first - first.gpos()) + gi_minimum;
  DASH_LOG_DEBUG("dash::min_element >", minimum,
                 "=", static_cast<value_t>(*minimum));

  return minimum;
}

template <
    class GlobIter,
    class Compare = std::greater<const typename GlobIter::value_type &> >
GlobIter max_element(
    const GlobIter &first,
    const GlobIter &last,
    Compare compare = Compare())
{
  // Same as min_element with different compare function
  return dash::min_element(first, last, compare);
}

template <class ElementType, class Compare = std::greater<ElementType &> >
const ElementType *max_element(
    const ElementType *first,
    const ElementType *last,
    Compare compare = Compare())
{
  // Same as min_element with different compare function
  return dash::min_element(first, last, compare);
}

} // namespace dash

#endif // DASH__ALGORITHM__MIN_MAX_H__