dune-grid-glue 2.8-git
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
17namespace Dune {
18namespace 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 = 0;
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
295template<typename OP, typename... Args>
296void 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
Definition: gridglue.hh:35
void MPI_AllApply(MPI_Comm mpicomm, OP &&op, const Args &... data)
apply an operator locally to a difstributed data set
Definition: ringcomm.hh:296