dune-grid-glue  2.7.0
ringcomm.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 /* IMPLEMENTATION OF CLASS G R I D G L U E */
4 
6 #define CheckMPIStatus(A,B) {}
7 
8 #include <mpi.h>
9 #include <functional>
10 #include <utility>
11 
12 #include <dune/common/fvector.hh>
13 #include <dune/common/hybridutilities.hh>
14 
15 #include <dune/geometry/type.hh>
16 
17 namespace Dune {
18 namespace Parallel {
19 
20  namespace Impl {
21 
23  template<typename T>
24  struct MPITypeInfo {};
25 
26  template<>
27  struct MPITypeInfo< int >
28  {
29  static const unsigned int size = 1;
30  static inline MPI_Datatype getType()
31  {
32  return MPI_INT;
33  }
34  };
35 
36  template<typename K, int N>
37  struct MPITypeInfo< Dune::FieldVector<K,N> >
38  {
39  static const unsigned int size = N;
40  static inline MPI_Datatype getType()
41  {
42  return Dune::MPITraits<K>::getType();
43  }
44  };
45 
46  template<>
47  struct MPITypeInfo< unsigned int >
48  {
49  static const unsigned int size = 1;
50  static inline MPI_Datatype getType()
51  {
52  return MPI_UNSIGNED;
53  }
54  };
55 
56  template<>
57  struct MPITypeInfo< Dune::GeometryType >
58  {
59  static const unsigned int size = 1;
60  static inline MPI_Datatype getType()
61  {
62  return Dune::MPITraits< Dune::GeometryType >::getType();
63  }
64  };
65 
66  template<typename T>
67  void MPI_SetVectorSize(
68  std::vector<T> & data,
69  MPI_Status & status)
70  {
71  typedef MPITypeInfo<T> Info;
72  int sz;
73  MPI_Get_count(&status, Info::getType(), &sz);
74  assert(sz%Info::size == 0);
75  data.resize(sz/Info::size);
76  }
77 
87  template<typename T>
88  void MPI_SendVectorInRing(
89  std::vector<T> & data,
90  std::vector<T> & next,
91  int tag,
92  int rightrank,
93  int leftrank,
94  MPI_Comm comm,
95  MPI_Request& r_send,
96  MPI_Request& r_recv
97  )
98  {
99  // mpi status stuff
100  int result DUNE_UNUSED;
101  result = 0;
102  typedef MPITypeInfo<T> Info;
103  // resize next buffer to maximum size
104  next.resize(next.capacity());
105  // send data (explicitly send data.size elements)
106  result =
107  MPI_Isend(
108  &(data[0]), Info::size*data.size(), Info::getType(), rightrank, tag,
109  comm, &r_send);
110  // receive up to maximum size. The acutal size is stored in the status
111  result =
112  MPI_Irecv(
113  &(next[0]), Info::size*next.size(), Info::getType(), leftrank, tag,
114  comm, &r_recv);
115  // // check result
116  // MPI_Status status;
117  // CheckMPIStatus(result, status);
118  }
119 
120  template<typename T>
121  using ptr_t = T*;
122 
123  /* these helper structs are needed as long as we still support
124  C++11, as we can't use variadic lambdas */
125  template<typename... Args>
126  struct call_MPI_SendVectorInRing
127  {
128  std::tuple<Args...> & remotedata;
129  std::tuple<Args...> & nextdata;
130  int & tag;
131  int & rightrank;
132  int & leftrank;
133  MPI_Comm & mpicomm;
134  std::array<MPI_Request,sizeof...(Args)> & requests_send;
135  std::array<MPI_Request,sizeof...(Args)> & requests_recv;
136 
137  template<typename I>
138  void operator()(I i)
139  {
140  MPI_SendVectorInRing(
141  std::get<i>(remotedata),
142  std::get<i>(nextdata),
143  tag+i,
144  rightrank, leftrank, mpicomm,
145  requests_send[i],
146  requests_recv[i]);
147  }
148  };
149  template<typename... Args>
150  struct call_MPI_SetVectorSize
151  {
152  std::tuple<Args...> & nextdata;
153  std::array<MPI_Status,sizeof...(Args)> & status_recv;
154 
155  template<typename I>
156  void operator()(I i)
157  {
158  MPI_SetVectorSize(std::get<i>(nextdata),status_recv[i]);
159  }
160  };
161 
162  template<typename OP, std::size_t... Indices, typename... Args>
163  void MPI_AllApply_impl(MPI_Comm mpicomm,
164  OP && op,
165  std::index_sequence<Indices...> indices,
166  const Args&... data)
167  {
168  constexpr std::size_t N = sizeof...(Args);
169  int myrank = 0;
170  int commsize = 0;
171 #if HAVE_MPI
172  MPI_Comm_rank(mpicomm, &myrank);
173  MPI_Comm_size(mpicomm, &commsize);
174 #endif // HAVE_MPI
175 
176  if (commsize > 1)
177  {
178 #ifdef DEBUG_GRIDGLUE_PARALLELMERGE
179  std::cout << myrank << " Start Communication, size " << commsize << std::endl;
180 #endif
181 
182  // get data sizes
183  std::array<unsigned int, N> size({ ((unsigned int)data.size())... });
184 
185  // communicate max data size
186  std::array<unsigned int, N> maxSize;
187  MPI_Allreduce(&size, &maxSize,
188  size.size(), MPI_UNSIGNED, MPI_MAX, mpicomm);
189 #ifdef DEBUG_GRIDGLUE_PARALLELMERGE
190  std::cout << myrank << " maxSize " << "done... " << std::endl;
191 #endif
192 
193  // allocate receiving buffers with maxsize to ensure sufficient buffer size for communication
194  std::tuple<Args...> remotedata { Args(maxSize[Indices])... };
195 
196  // copy local data to receiving buffer
197  remotedata = std::tie(data...);
198 
199  // allocate second set of receiving buffers necessary for async communication
200  std::tuple<Args...> nextdata { Args(maxSize[Indices])... };
201 
202  // communicate data in the ring
203  int rightrank = (myrank + 1 + commsize) % commsize;
204  int leftrank = (myrank - 1 + commsize) % commsize;
205 
206  std::cout << myrank << ": size = " << commsize << std::endl;
207  std::cout << myrank << ": left = " << leftrank
208  << " right = " << rightrank << std::endl;
209 
210  // currently the remote data is our own data
211  int remoterank = myrank;
212 
213  for (int i=1; i<commsize; i++)
214  {
215  // in this iteration we will receive data from nextrank
216  int nextrank = (myrank - i + commsize) % commsize;
217 
218  std::cout << myrank << ": next = " << nextrank << std::endl;
219 
220  // send remote data to right neighbor and receive from left neighbor
221  std::array<MPI_Request,N> requests_send;
222  std::array<MPI_Request,N> requests_recv;
223 
224  int tag = 12345678;
225  Dune::Hybrid::forEach(indices,
226  // [&](auto i){
227  // MPI_SendVectorInRing(
228  // std::get<i>(remotedata),
229  // std::get<i>(nextdata),
230  // tag+i,
231  // rightrank, leftrank, mpicomm,
232  // requests_send[i],
233  // requests_recv[i]);
234  // });
235  call_MPI_SendVectorInRing<Args...>({
236  remotedata,
237  nextdata,
238  tag,
239  rightrank, leftrank, mpicomm,
240  requests_send,
241  requests_recv
242  }));
243 
244  // apply operator
245  op(remoterank,std::get<Indices>(remotedata)...);
246 
247  // wait for communication to finalize
248  std::array<MPI_Status,N> status_send;
249  std::array<MPI_Status,N> status_recv;
250  MPI_Waitall(N,&requests_recv[0],&status_recv[0]);
251 
252  // we finished receiving from nextrank and thus remoterank = nextrank
253  remoterank = nextrank;
254 
255  // get current data sizes
256  // and resize vectors
257  Dune::Hybrid::forEach(indices,
258  // [&](auto i){
259  // MPI_SetVectorSize(std::get<i>(nextdata),status_recv[i]);
260  // });
261  call_MPI_SetVectorSize<Args...>({
262  nextdata, status_recv
263  }));
264 
265  MPI_Waitall(N,&requests_send[0],&status_send[0]);
266 
267  // swap the communication buffers
268  std::swap(remotedata,nextdata);
269  }
270 
271  // last apply (or the only one in the case of sequential application)
272  op(remoterank,std::get<Indices>(remotedata)...);
273  }
274  else // sequential
275  {
276  op(myrank,data...);
277  }
278  }
279 
280  } // end namespace Impl
281 
295 template<typename OP, typename... Args>
296 void MPI_AllApply(MPI_Comm mpicomm,
297  OP && op,
298  const Args& ... data)
299 {
300  Impl::MPI_AllApply_impl(
301  mpicomm,
302  std::forward<OP>(op),
303  std::make_index_sequence<sizeof...(Args)>(),
304  data...
305  );
306 }
307 
308 } // end namespace Parallel
309 } // end namespace Dune
Dune
Definition: gridglue.hh:35
Dune::Parallel::MPI_AllApply
void MPI_AllApply(MPI_Comm mpicomm, OP &&op, const Args &... data)
apply an operator locally to a difstributed data set
Definition: ringcomm.hh:296