Muster
 All Classes Namespaces Files Functions Variables Typedefs Macros
gather.h
Go to the documentation of this file.
1 #ifndef MUSTER_GATHER_H
2 #define MUSTER_GATHER_H
3 
4 #include <mpi.h>
5 #include <numeric>
6 #include <algorithm>
7 #include "mpi_bindings.h"
8 #include "mpi_utils.h"
9 #include "binomial.h"
10 
11 namespace cluster {
12 
13  ///
14  /// Packs and gathers a buffer full of packed representation of src's. All packed
15  /// src's are stored in the destination buffer in binomial traversal order on completion.
16  /// That is, element i in dest is the value for rank binomial.reverse_relative_rank(i).
17  ///
18  /// This allows you to use native MPI operations like bcast on the packed buffer
19  /// once it's gathered.
20  ///
21  /// @see gather() for a version of this that will unpack the gathered data for you.
22  ///
23  template <class T>
24  void gather_packed(const T& src, std::vector<char>& dest, const binomial_embedding binomial, MPI_Comm comm) {
25  int rank, size;
26  CMPI_Comm_rank(comm, &rank);
27  CMPI_Comm_size(comm, &size);
28 
29  int parent = binomial.parent(rank);
30  std::vector<int> children = binomial.children(rank);
31 
32  std::vector<int> sizes; // sizes of buffers to receive.
33  sizes.push_back(src.packed_size(comm)); // size of local packed data
34 
35  for (size_t i=0; i < children.size(); i++) {
36  // Receive sizes from all children
37  int cur_size;
38  CMPI_Recv(&cur_size, 1, MPI_INT, children[i], 0, comm, MPI_STATUS_IGNORE);
39  sizes.push_back(cur_size);
40  }
41 
42  // construct offsets to receive buffers into.
43  std::vector<int> offsets;
44  offsets.push_back(0);
45  std::partial_sum(sizes.begin(), sizes.end(), back_inserter(offsets));
46 
47  // create a buffer to receive into, then to send to parent
48  std::vector<char> sendbuf(accumulate(sizes.begin(), sizes.end(), 0));
49 
50  // pack local object before its children
51  int pos = 0;
52  src.pack(&sendbuf[0], sendbuf.size(), &pos, comm);
53 
54  // receive from all children
55  for (size_t i = 0; i < children.size(); i++) {
56  CMPI_Recv(&sendbuf[offsets[i+1]], sizes[i+1], MPI_PACKED, children[i], 0, comm, MPI_STATUS_IGNORE);
57  }
58 
59  if (parent != -1) {
60  // Done receiving. Send to parent.
61  int size = sendbuf.size();
62  CMPI_Send(&size, 1, MPI_INT, parent, 0, comm);
63  CMPI_Send(&sendbuf[0], sendbuf.size(), MPI_PACKED, parent, 0, comm);
64 
65  } else {
66  // put packed data in the destination.
67  dest.swap(sendbuf);
68  }
69  }
70 
71 
72  ///
73  /// Unpacks a packed vector in binomial order into objects in rank order in the destination vector.
74  ///
75  template <class T>
76  void unpack_binomial(const std::vector<char>& src, std::vector<T>& dest, const binomial_embedding binomial,
77  MPI_Comm comm) {
78  int pos = 0;
79  dest.resize(binomial.size());
80  for (size_t i=0; i < binomial.size(); i++) {
81  dest[binomial.reverse_relative_rank(i)] = T::unpack(const_cast<char*>(&src[0]), src.size(), &pos, comm);
82  }
83  }
84 
85 
86  ///
87  /// Binomial gather of char buffers into a single agglomerated clump of buffers
88  ///
89  template <class T>
90  void gather(const T& src, std::vector<T>& dest, MPI_Comm comm, int root = 0) {
91  int rank, size;
92  CMPI_Comm_rank(comm, &rank);
93  CMPI_Comm_size(comm, &size);
94 
95  // gather everything to a packed buffer at the root.
96  binomial_embedding binomial(size, root);
97  std::vector<char> packed;
98  gather_packed(src, packed, binomial, comm);
99 
100  // now unpack everything.
101  if (rank == root) {
102  unpack_binomial(packed, dest, binomial, comm);
103  }
104  }
105 
106 
107  ///
108  /// Allgather for variable-length data.
109  ///
110  template <class T>
111  void allgather(const T& src, std::vector<T>& dest, MPI_Comm comm, int root = 0) {
112  int rank, size;
113  CMPI_Comm_rank(comm, &rank);
114  CMPI_Comm_size(comm, &size);
115 
116  // gather everything to a packed buffer at the root.
117  binomial_embedding binomial(size, root);
118  std::vector<char> packed;
119  gather_packed(src, packed, binomial, comm);
120 
121  size_t packed_size = packed.size();
122  CMPI_Bcast(&packed_size, 1, MPI_SIZE_T, root, comm);
123 
124  packed.resize(packed_size);
125  CMPI_Bcast(const_cast<char*>(&packed[0]), packed.size(), MPI_PACKED, root, comm);
126 
127  unpack_binomial(packed, dest, binomial, comm);
128  }
129 
130 } // namespace cluster
131 
132 #endif // MUSTER_GATHER_H
#define MPI_SIZE_T
Definition: mpi_utils.h:39
void unpack_binomial(const std::vector< char > &src, std::vector< T > &dest, const binomial_embedding binomial, MPI_Comm comm)
Unpacks a packed vector in binomial order into objects in rank order in the destination vector...
Definition: gather.h:76
void gather_packed(const T &src, std::vector< char > &dest, const binomial_embedding binomial, MPI_Comm comm)
Packs and gathers a buffer full of packed representation of src&#39;s.
Definition: gather.h:24
void allgather(const T &src, std::vector< T > &dest, MPI_Comm comm, int root=0)
Allgather for variable-length data.
Definition: gather.h:111
#define CMPI_Comm_rank
Definition: mpi_bindings.h:81
int reverse_relative_rank(int rank) const
Reverse rank permutation.
Definition: binomial.cpp:17
std::vector< int > children(int rank) const
Same as get_children, but returns vector.
Definition: binomial.cpp:21
int parent(int rank) const
Get the parent of a particular rank.
Definition: binomial.cpp:27
Overloaded utility functions to convert between arbitrary C/C++ types and MPI types, custom typedefs for cstdlib types like size_t, and a wrapper for MPI_Pack_Size.
#define CMPI_Send
Definition: mpi_bindings.h:86
void gather(const T &src, std::vector< T > &dest, MPI_Comm comm, int root=0)
Binomial gather of char buffers into a single agglomerated clump of buffers.
Definition: gather.h:90
#define CMPI_Bcast
Definition: mpi_bindings.h:80
#defines for switching between MPI and PMPI bindings.
#define CMPI_Recv
Definition: mpi_bindings.h:85
#define CMPI_Comm_size
Definition: mpi_bindings.h:82
Muster. Copyright © 2010, Lawrence Livermore National Laboratory, LLNL-CODE-433662.
Distribution of Muster and its documentation is subject to terms of the Muster LICENSE.
Generated on Thu Sep 1 2016 using Doxygen 1.8.5