static void gather(const communicator& comm, const std::vector<T>& in, std::vector< std::vector<T> >& out, int root)
    {
      std::vector<int>  counts(comm.size());
      Collectives<int,void*>::gather(comm, (int) in.size(), counts, root);

      std::vector<int>  offsets(comm.size(), 0);
      for (unsigned i = 1; i < offsets.size(); ++i)
        offsets[i] = offsets[i-1] + counts[i-1];

      std::vector<T> buffer(offsets.back() + counts.back());
      MPI_Gatherv(Datatype::address(const_cast<T&>(in[0])),
                  in.size(),
                  Datatype::datatype(),
                  Datatype::address(buffer[0]),
                  &counts[0],
                  &offsets[0],
                  Datatype::datatype(),
                  root, comm);

      out.resize(comm.size());
      size_t cur = 0;
      for (unsigned i = 0; i < (unsigned)comm.size(); ++i)
      {
          out[i].reserve(counts[i]);
          for (unsigned j = 0; j < (unsigned)counts[i]; ++j)
              out[i].push_back(buffer[cur++]);
      }
    }

    void operator()(const communicator& comm, const config& cfg,
                    T alpha, MatrixA& A, MatrixB& B, T beta, MatrixC& C)
    {
        using namespace matrix_constants;

        len_type m = (Mat == MAT_A ? A.length(0) : Mat == MAT_B ? B.length(0) : C.length(0));
        len_type n = (Mat == MAT_A ? A.length(1) : Mat == MAT_B ? B.length(1) : C.length(1));

        if (!rscat)
        {
            if (comm.master())
            {
                scat_buffer = Pool.allocate<stride_type>(2*m + 2*n);
                rscat = scat_buffer.get<stride_type>();
            }

            comm.broadcast(rscat);

            cscat = rscat+m;
            rbs = cscat+n;
            cbs = rbs+m;
        }

        matrify_and_run<Mat>(*this, comm, cfg, alpha, A, B, beta, C);
    }

void
broadcast_test(const communicator& comm, const T& bc_value,
               const char* kind, int root = -1)
{
  if (root == -1) {
    for (root = 0; root < comm.size(); ++root)
      broadcast_test(comm, bc_value, kind, root);
  } else {
    using boost::mpi::broadcast;

    T value;
    if (comm.rank() == root) {
      value = bc_value;
      std::cout << "Broadcasting " << kind << " from root " << root << "...";
      std::cout.flush();
    }

    broadcast(comm, value, root);
    BOOST_CHECK(value == bc_value);
    if (comm.rank() == root && value == bc_value)
      std::cout << "OK." << std::endl;
  }

  (comm.barrier)();
}

void
gather_test(const communicator& comm, Generator generator,
            const char* kind, int root = -1)
{
  typedef typename Generator::result_type value_type;
  value_type value = generator(comm.rank());

  if (root == -1) {
    for (root = 0; root < comm.size(); ++root)
      gather_test(comm, generator, kind, root);
  } else {
    using boost::mpi::gather;

    std::vector<value_type> values;
    if (comm.rank() == root) {
      std::cout << "Gathering " << kind << " from root "
                << root << "..." << std::endl;
    }

    gather(comm, value, values, root);

    if (comm.rank() == root) {
      std::vector<value_type> expected_values;
      for (int p = 0; p < comm.size(); ++p)
        expected_values.push_back(generator(p));
      BOOST_CHECK(values == expected_values);
    } else {
      BOOST_CHECK(values.empty());
    }
  }

  (comm.barrier)();
}

void block_scatter(const communicator& comm, tensor_matrix<T>& A,
                   stride_type* rscat, len_type MB, stride_type* rbs,
                   stride_type* cscat, len_type NB, stride_type* cbs)
{
    len_type m = A.length(0);
    len_type n = A.length(1);

    len_type first, last;
    std::tie(first, last, std::ignore) = comm.distribute_over_threads(m, MB);

    A.length(0, last-first);
    A.shift(0, first);
    A.fill_block_scatter(0, rscat+first, MB, rbs+first/MB);
    A.shift(0, -first);
    A.length(0, m);

    std::tie(first, last, std::ignore) = comm.distribute_over_threads(n, NB);

    A.length(1, last-first);
    A.shift(1, first);
    A.fill_block_scatter(1, cscat+first, NB, cbs+first/NB);
    A.shift(1, -first);
    A.length(1, n);

    comm.barrier();
}

void
scatter_test(const communicator& comm, Generator generator,
            const char* kind, int root = -1)
{
  typedef typename Generator::result_type value_type;

  if (root == -1) {
    for (root = 0; root < comm.size(); ++root)
      scatter_test(comm, generator, kind, root);
  } else {
    using boost::mpi::scatter;

    value_type value;

    if (comm.rank() == root) {
      std::vector<value_type> values;

      for (int p = 0; p < comm.size(); ++p)
        values.push_back(generator(p));

      if (comm.rank() == root) {
        std::cout << "Scattering " << kind << " from root " << root << "...";
        std::cout.flush();
      }

      scatter(comm, values, value, root);
    } else {
      scatter(comm, value, root);
    }

    BOOST_CHECK(value == generator(comm.rank()));
  }

  (comm.barrier)();
}

 //---------
 static void broadcast(communicator c, A &a, int root) {
     check_is_contiguous(a);
     auto sh = a.shape();
     MPI_Bcast(&sh[0], sh.size(), mpi_datatype<typename decltype(sh)::value_type>::invoke(), root, c.get());
     if (c.rank() != root) a.resize(sh);
     MPI_Bcast(a.data_start(), a.domain().number_of_elements(), D(), root, c.get());
 }

 std::vector<T> mpi_reduce(std::vector<T> const &a, communicator c, int root, bool all, MPI_Op op, std::true_type) {
  std::vector<T> b(a.size());
  if (!all)
   MPI_Reduce((void *)a.data(), b.data(), a.size(), mpi_datatype<T>(), op, root, c.get());
  else
   MPI_Allreduce((void *)a.data(), b.data(), a.size(), mpi_datatype<T>(), op, c.get());
  return b;
 }

 template <typename T> REQUIRES_IS_BASIC(T, T) mpi_reduce(T a, communicator c = {}, int root = 0, bool all = false, MPI_Op op = MPI_SUM) {
   T b;
   auto d = datatype<T>();
   if (!all)
     MPI_Reduce(&a, &b, 1, d, op, root, c.get());
   else
     MPI_Allreduce(&a, &b, 1, d, op, c.get());
   return b;
 }

void
gather(const communicator& comm, const T* in_values, int n, 
       std::vector<T>& out_values, int root)
{
  if (comm.rank() == root) {
    out_values.resize(comm.size() * n);
    ::boost::mpi::gather(comm, in_values, n, &out_values[0], root);
  } 
  else
    ::boost::mpi::gather(comm, in_values, n, root);
}

  void
  gather_impl(const communicator& comm, const T* in_values, int n, 
              T* out_values, int root, mpl::false_)
  {
    int tag = environment::collectives_tag();
    int size = comm.size();

    for (int src = 0; src < size; ++src) {
      if (src == root)
        std::copy(in_values, in_values + n, out_values + n * src);
      else
        comm.recv(src, tag, out_values + n * src, n);
    }
  }

object gather(const communicator& comm, object value, int root)
{
  if (comm.rank() == root) {
    std::vector<object> values;
    boost::mpi::gather(comm, value, values, root);

    boost::python::list l;
    for (int i = 0; i < comm.size(); ++i)
      l.append(values[i]);
    return boost::python::tuple(l);
  } else {
    boost::mpi::gather(comm, value, root);
    return object();
  }
}

object scatter(const communicator& comm, object values, int root)
{
  object result;

  if (comm.rank() == root) {
    std::vector<object> values_vec(comm.size());
    object iterator = object(handle<>(PyObject_GetIter(values.ptr())));
    for (int i = 0; i < comm.size(); ++i)
      values_vec[i] = object(handle<>(PyIter_Next(iterator.ptr())));

    boost::mpi::scatter(comm, values_vec, result, root);
  } else {
    boost::mpi::scatter(comm, result, root);
  }
  return result;
}

int espresso::esutil::Collectives::locateItem(bool here, int controller, communicator world) {
  int node = here ? world.rank() : UniqueReduce::NotHere;

  if (world.rank() != controller) {
    reduce(world, node, UniqueReduce(), controller);
    return None;
  }
  else {
    int owner;
    reduce(world, node, owner, UniqueReduce(), controller);
    if (owner == UniqueReduce::Duplicate) {
      throw DuplicateError();
    }
    return owner;
  }
}

object all_to_all(const communicator& comm, object in_values)
{
  // Build input values
  std::vector<object> in_values_vec(comm.size());
  object iterator = object(handle<>(PyObject_GetIter(in_values.ptr())));
  for (int i = 0; i < comm.size(); ++i)
    in_values_vec[i] = object(handle<>(PyIter_Next(iterator.ptr())));

  std::vector<object> out_values_vec(comm.size());
  boost::mpi::all_to_all(comm, in_values_vec, out_values_vec);

  boost::python::list l;
  for (int i = 0; i < comm.size(); ++i)
    l.append(out_values_vec[i]);
  return boost::python::tuple(l);
}

void
scatterv_impl(const communicator& comm, T* out_values, int n, int root, 
              mpl::false_ isnt_mpi_type)
{
  assert(root != comm.rank());
  scatterv_impl(comm, (T const*)0, out_values, n, (int const*)0, (int const*)0, root, isnt_mpi_type);
}

 void all_gather(const communicator& comm, const T& in_val, std::vector<T>& out_vals) {
     out_vals.resize(comm.size());
     MPI_Allgather((void*)&in_val, 1, detail::mpi_type<T>(),
                   &out_vals.front(), 1, detail::mpi_type<T>(),
                   comm);
     // throw std::logic_error(std::string("all_gather() is not implemented, called for type T=")
     //                        +typeid(T).name());
 }

object
communicator_iprobe(const communicator& comm, int source, int tag)
{
  if (pdalboost::optional<status> result = comm.iprobe(source, tag))
    return object(*result);
  else
    return object();
}

request_with_value 
communicator_irecv(const communicator& comm, int source, int tag)
{
  pdalboost::shared_ptr<object> result(new object());
  request_with_value req(comm.irecv(source, tag, *result));
  req.m_internal_value = result;
  return req;
}

void
ring_array_test(const communicator& comm, const T* pass_values,
                int n, const char* kind, int root = 0)
{
  T* transferred_values = new T[n];
  int rank = comm.rank();
  int size = comm.size();

  if (rank == root) {

    std::cout << "Passing " << kind << " array around a ring from root "
              << root  << "...";
    comm.send((rank + 1) % size, 0, pass_values, n);
    comm.recv((rank + size - 1) % size, 0, transferred_values, n);
    bool okay = std::equal(pass_values, pass_values + n,
                           transferred_values);
    BOOST_CHECK(okay);
    if (okay) std::cout << " OK." << std::endl;
  } else {
    status stat = comm.probe(boost::mpi::any_source, 0);
    boost::optional<int> num_values = stat.template count<T>();
    if (boost::mpi::is_mpi_datatype<T>())
      BOOST_CHECK(num_values && *num_values == n);
    else
      BOOST_CHECK(!num_values || *num_values == n);     
    comm.recv(stat.source(), 0, transferred_values, n);
    BOOST_CHECK(std::equal(pass_values, pass_values + n,
                           transferred_values));
    comm.send((rank + 1) % size, 0, transferred_values, n);
  }
  (comm.barrier)();
  delete [] transferred_values;
}