dune-grid-glue 2.8.0
Loading...
Searching...
No Matches
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 [[maybe_unused]] int result = 0;
101 typedef MPITypeInfo<T> Info;
102 // resize next buffer to maximum size
103 next.resize(next.capacity());
104 // send data (explicitly send data.size elements)
105 result =
106 MPI_Isend(
107 &(data[0]), Info::size*data.size(), Info::getType(), rightrank, tag,
108 comm, &r_send);
109 // receive up to maximum size. The acutal size is stored in the status
110 result =
111 MPI_Irecv(
112 &(next[0]), Info::size*next.size(), Info::getType(), leftrank, tag,
113 comm, &r_recv);
114 // // check result
115 // MPI_Status status;
116 // CheckMPIStatus(result, status);
117 }
118
119 template<typename T>
120 using ptr_t = T*;
121
122 /* these helper structs are needed as long as we still support
123 C++11, as we can't use variadic lambdas */
124 template<typename... Args>
125 struct call_MPI_SendVectorInRing
126 {
127 std::tuple<Args...> & remotedata;
128 std::tuple<Args...> & nextdata;
129 int & tag;
130 int & rightrank;
131 int & leftrank;
132 MPI_Comm & mpicomm;
133 std::array<MPI_Request,sizeof...(Args)> & requests_send;
134 std::array<MPI_Request,sizeof...(Args)> & requests_recv;
135
136 template<typename I>
137 void operator()(I i)
138 {
139 MPI_SendVectorInRing(
140 std::get<i>(remotedata),
141 std::get<i>(nextdata),
142 tag+i,
143 rightrank, leftrank, mpicomm,
144 requests_send[i],
145 requests_recv[i]);
146 }
147 };
148 template<typename... Args>
149 struct call_MPI_SetVectorSize
150 {
151 std::tuple<Args...> & nextdata;
152 std::array<MPI_Status,sizeof...(Args)> & status_recv;
153
154 template<typename I>
155 void operator()(I i)
156 {
157 MPI_SetVectorSize(std::get<i>(nextdata),status_recv[i]);
158 }
159 };
160
161 template<typename OP, std::size_t... Indices, typename... Args>
162 void MPI_AllApply_impl(MPI_Comm mpicomm,
163 OP && op,
164 std::index_sequence<Indices...> indices,
165 const Args&... data)
166 {
167 constexpr std::size_t N = sizeof...(Args);
168 int myrank = 0;
169 int commsize = 0;
170#if HAVE_MPI
171 MPI_Comm_rank(mpicomm, &myrank);
172 MPI_Comm_size(mpicomm, &commsize);
173#endif // HAVE_MPI
174
175 if (commsize > 1)
176 {
177#ifdef DEBUG_GRIDGLUE_PARALLELMERGE
178 std::cout << myrank << " Start Communication, size " << commsize << std::endl;
179#endif
180
181 // get data sizes
182 std::array<unsigned int, N> size({ ((unsigned int)data.size())... });
183
184 // communicate max data size
185 std::array<unsigned int, N> maxSize;
186 MPI_Allreduce(&size, &maxSize,
187 size.size(), MPI_UNSIGNED, MPI_MAX, mpicomm);
188#ifdef DEBUG_GRIDGLUE_PARALLELMERGE
189 std::cout << myrank << " maxSize " << "done... " << std::endl;
190#endif
191
192 // allocate receiving buffers with maxsize to ensure sufficient buffer size for communication
193 std::tuple<Args...> remotedata { Args(maxSize[Indices])... };
194
195 // copy local data to receiving buffer
196 remotedata = std::tie(data...);
197
198 // allocate second set of receiving buffers necessary for async communication
199 std::tuple<Args...> nextdata { Args(maxSize[Indices])... };
200
201 // communicate data in the ring
202 int rightrank = (myrank + 1 + commsize) % commsize;
203 int leftrank = (myrank - 1 + commsize) % commsize;
204
205 std::cout << myrank << ": size = " << commsize << std::endl;
206 std::cout << myrank << ": left = " << leftrank
207 << " right = " << rightrank << std::endl;
208
209 // currently the remote data is our own data
210 int remoterank = myrank;
211
212 for (int i=1; i<commsize; i++)
213 {
214 // in this iteration we will receive data from nextrank
215 int nextrank = (myrank - i + commsize) % commsize;
216
217 std::cout << myrank << ": next = " << nextrank << std::endl;
218
219 // send remote data to right neighbor and receive from left neighbor
220 std::array<MPI_Request,N> requests_send;
221 std::array<MPI_Request,N> requests_recv;
222
223 int tag = 0;
224 Dune::Hybrid::forEach(indices,
225 // [&](auto i){
226 // MPI_SendVectorInRing(
227 // std::get<i>(remotedata),
228 // std::get<i>(nextdata),
229 // tag+i,
230 // rightrank, leftrank, mpicomm,
231 // requests_send[i],
232 // requests_recv[i]);
233 // });
234 call_MPI_SendVectorInRing<Args...>({
235 remotedata,
236 nextdata,
237 tag,
238 rightrank, leftrank, mpicomm,
239 requests_send,
240 requests_recv
241 }));
242
243 // apply operator
244 op(remoterank,std::get<Indices>(remotedata)...);
245
246 // wait for communication to finalize
247 std::array<MPI_Status,N> status_send;
248 std::array<MPI_Status,N> status_recv;
249 MPI_Waitall(N,&requests_recv[0],&status_recv[0]);
250
251 // we finished receiving from nextrank and thus remoterank = nextrank
252 remoterank = nextrank;
253
254 // get current data sizes
255 // and resize vectors
256 Dune::Hybrid::forEach(indices,
257 // [&](auto i){
258 // MPI_SetVectorSize(std::get<i>(nextdata),status_recv[i]);
259 // });
260 call_MPI_SetVectorSize<Args...>({
261 nextdata, status_recv
262 }));
263
264 MPI_Waitall(N,&requests_send[0],&status_send[0]);
265
266 // swap the communication buffers
267 std::swap(remotedata,nextdata);
268 }
269
270 // last apply (or the only one in the case of sequential application)
271 op(remoterank,std::get<Indices>(remotedata)...);
272 }
273 else // sequential
274 {
275 op(myrank,data...);
276 }
277 }
278
279 } // end namespace Impl
280
294template<typename OP, typename... Args>
295void MPI_AllApply(MPI_Comm mpicomm,
296 OP && op,
297 const Args& ... data)
298{
299 Impl::MPI_AllApply_impl(
300 mpicomm,
301 std::forward<OP>(op),
302 std::make_index_sequence<sizeof...(Args)>(),
303 data...
304 );
305}
306
307} // end namespace Parallel
308} // 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:295