MeDiPack  1.3.0
A Message Differentiation Package
SciComp TU Kaiserslautern
Loading...
Searching...
No Matches
reverseFunctions.hpp
Go to the documentation of this file.
1/*
2 * MeDiPack, a Message Differentiation Package
3 *
4 * Copyright (C) 2015-2024 Chair for Scientific Computing (SciComp), University of Kaiserslautern-Landau
5 * Homepage: http://scicomp.rptu.de
6 * Contact: Prof. Nicolas R. Gauger (codi@scicomp.uni-kl.de)
7 *
8 * Lead developers: Max Sagebaum (SciComp, University of Kaiserslautern-Landau)
9 *
10 * This file is part of MeDiPack (http://scicomp.rptu.de/software/codi).
11 *
12 * MeDiPack is free software: you can redistribute it and/or
13 * modify it under the terms of the GNU Lesser General Public
14 * License as published by the Free Software Foundation, either
15 * version 3 of the License, or (at your option) any later version.
16 *
17 * MeDiPack is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty of
19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
20 *
21 * See the GNU Lesser General Public License for more details.
22 * You should have received a copy of the GNU
23 * Lesser General Public License along with MeDiPack.
24 * If not, see <http://www.gnu.org/licenses/>.
25 *
26 * Authors: Max Sagebaum, Tim Albring (SciComp, University of Kaiserslautern-Landau)
27 */
28
29#pragma once
30
31#include "ampiMisc.h"
32#include "async.hpp"
33#include "enums.hpp"
34#include "message.hpp"
36
40namespace medi {
41
42#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
43 template<typename DATATYPE>
44 void AMPI_Send_adj(typename DATATYPE::AdjointType* bufAdjoints, int bufSize, int count, DATATYPE* datatype, int dest, int tag, AMPI_Comm comm) {
45 MEDI_UNUSED(count);
46 MPI_Recv(bufAdjoints, bufSize, datatype->getADTool().getAdjointMpiType(), dest, tag, comm, MPI_STATUS_IGNORE);
47 }
48#endif
49
50#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
51 template<typename DATATYPE>
52 void AMPI_Isend_adj(typename DATATYPE::AdjointType* bufAdjoints, int bufSize, int count, DATATYPE* datatype, int dest, int tag, AMPI_Comm comm, AMPI_Request* request) {
53 MEDI_UNUSED(count);
54 MPI_Irecv(bufAdjoints, bufSize, datatype->getADTool().getAdjointMpiType(), dest, tag, comm, &request->request);
55 }
56#endif
57
58#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
59 template<typename DATATYPE>
60 void AMPI_Bsend_adj(typename DATATYPE::AdjointType* bufAdjoints, int bufSize, int count, DATATYPE* datatype, int dest, int tag, AMPI_Comm comm) {
61 MEDI_UNUSED(count);
62 MPI_Recv(bufAdjoints, bufSize, datatype->getADTool().getAdjointMpiType(), dest, tag, comm, MPI_STATUS_IGNORE);
63 }
64#endif
65
66#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
67 template<typename DATATYPE>
68 void AMPI_Ibsend_adj(typename DATATYPE::AdjointType* bufAdjoints, int bufSize, int count, DATATYPE* datatype, int dest, int tag, AMPI_Comm comm, AMPI_Request* request) {
69 MEDI_UNUSED(count);
70 MPI_Irecv(bufAdjoints, bufSize, datatype->getADTool().getAdjointMpiType(), dest, tag, comm, &request->request);
71 }
72#endif
73
74#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
75 template<typename DATATYPE>
76 void AMPI_Ssend_adj(typename DATATYPE::AdjointType* bufAdjoints, int bufSize, int count, DATATYPE* datatype, int dest, int tag, AMPI_Comm comm) {
77 MEDI_UNUSED(count);
78 MPI_Recv(bufAdjoints, bufSize, datatype->getADTool().getAdjointMpiType(), dest, tag, comm, MPI_STATUS_IGNORE);
79 }
80#endif
81
82#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
83 template<typename DATATYPE>
84 void AMPI_Issend_adj(typename DATATYPE::AdjointType* bufAdjoints, int bufSize, int count, DATATYPE* datatype, int dest, int tag, AMPI_Comm comm, AMPI_Request* request) {
85 MEDI_UNUSED(count);
86 MPI_Irecv(bufAdjoints, bufSize, datatype->getADTool().getAdjointMpiType(), dest, tag, comm, &request->request);
87 }
88#endif
89
90#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
91 template<typename DATATYPE>
92 void AMPI_Rsend_adj(typename DATATYPE::AdjointType* bufAdjoints, int bufSize, int count, DATATYPE* datatype, int dest, int tag, AMPI_Comm comm) {
93 MEDI_UNUSED(count);
94 MPI_Recv(bufAdjoints, bufSize, datatype->getADTool().getAdjointMpiType(), dest, tag, comm, MPI_STATUS_IGNORE);
95 }
96#endif
97
98#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
99 template<typename DATATYPE>
100 void AMPI_Irsend_adj(typename DATATYPE::AdjointType* bufAdjoints, int bufSize, int count, DATATYPE* datatype, int dest, int tag, AMPI_Comm comm, AMPI_Request* request) {
101 MEDI_UNUSED(count);
102 MPI_Irecv(bufAdjoints, bufSize, datatype->getADTool().getAdjointMpiType(), dest, tag, comm, &request->request);
103 }
104#endif
105
106#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
107 template<typename DATATYPE>
108 void AMPI_Recv_adj(typename DATATYPE::AdjointType* bufAdjoints, int bufSize, int count, DATATYPE* datatype, int src, int tag, AMPI_Comm comm, AMPI_Status* status, RecvAdjCall reverse_call) {
109 MEDI_UNUSED(count);
110 MEDI_UNUSED(status);
111 if (RecvAdjCall::Send == reverse_call) {
112 MPI_Send(bufAdjoints, bufSize, datatype->getADTool().getAdjointMpiType(), src, tag, comm);
113 } else if (RecvAdjCall::Bsend == reverse_call) {
114 MPI_Bsend(bufAdjoints, bufSize, datatype->getADTool().getAdjointMpiType(), src, tag, comm);
115 } else if (RecvAdjCall::Rsend == reverse_call) {
116 MPI_Rsend(bufAdjoints, bufSize, datatype->getADTool().getAdjointMpiType(), src, tag, comm);
117 } else if (RecvAdjCall::Ssend == reverse_call) {
118 MPI_Ssend(bufAdjoints, bufSize, datatype->getADTool().getAdjointMpiType(), src, tag, comm);
119 } else {
120 MEDI_EXCEPTION("Unimplemented case for RecvAdjCall %d.", (int)reverse_call);
121 }
122 }
123#endif
124
125#if MEDI_MPI_VERSION_3_0 <= MEDI_MPI_TARGET
126 template<typename DATATYPE>
127 void AMPI_Mrecv_adj(typename DATATYPE::AdjointType* bufAdjoints, int bufSize, int count, DATATYPE* datatype, AMPI_Message* message, AMPI_Status* status, RecvAdjCall reverse_call) {
128 AMPI_Recv_adj(bufAdjoints, bufSize, count, datatype, message->src, message->tag, message->comm, status, reverse_call);
129 }
130#endif
131
132#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
133 template<typename DATATYPE>
134 void AMPI_Irecv_adj(typename DATATYPE::AdjointType* bufAdjoints, int bufSize, int count, DATATYPE* datatype, int src, int tag, AMPI_Comm comm, AMPI_Request* request, IrecvAdjCall reverse_call) {
135 MEDI_UNUSED(count);
136 if (IrecvAdjCall::Isend == reverse_call) {
137 MPI_Isend(bufAdjoints, bufSize, datatype->getADTool().getAdjointMpiType(), src, tag, comm, &request->request);
138 } else if (IrecvAdjCall::Ibsend == reverse_call) {
139 MPI_Ibsend(bufAdjoints, bufSize, datatype->getADTool().getAdjointMpiType(), src, tag, comm, &request->request);
140 } else if (IrecvAdjCall::Irsend == reverse_call) {
141 MPI_Irsend(bufAdjoints, bufSize, datatype->getADTool().getAdjointMpiType(), src, tag, comm, &request->request);
142 } else if (IrecvAdjCall::Issend == reverse_call) {
143 MPI_Issend(bufAdjoints, bufSize, datatype->getADTool().getAdjointMpiType(), src, tag, comm, &request->request);
144 } else {
145 MEDI_EXCEPTION("Unimplemented case for IrecvAdjCall %d.", (int)reverse_call);
146 }
147 }
148#endif
149
150#if MEDI_MPI_VERSION_3_0 <= MEDI_MPI_TARGET
151 template<typename DATATYPE>
152 void AMPI_Imrecv_adj(typename DATATYPE::AdjointType* bufAdjoints, int bufSize, int count, DATATYPE* datatype, AMPI_Message* message, AMPI_Request* request, IrecvAdjCall reverse_call) {
153 MEDI_UNUSED(count);
154 AMPI_Irecv_adj(bufAdjoints, bufSize, count, datatype, message->src, message->tag, message->comm, request, reverse_call);
155 }
156#endif
157
158#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
159 template<typename SENDTYPE, typename RECVTYPE>
160 void AMPI_Sendrecv_adj(typename SENDTYPE::AdjointType* sendbuf, int sendbufSize, int sendcount, SENDTYPE* sendtype, int dest, int sendtag,
161#endif
162 typename RECVTYPE::AdjointType* recvbuf, int recvbufSize, int recvcount, RECVTYPE* recvtype, int source, int recvtag, AMPI_Comm comm, AMPI_Status* status) {
163
164 MEDI_UNUSED(sendcount);
165 MEDI_UNUSED(recvcount);
166 MPI_Sendrecv(recvbuf, recvbufSize, recvtype->getADTool().getAdjointMpiType(), source, recvtag, sendbuf, sendbufSize, sendtype->getADTool().getAdjointMpiType(), dest, sendtag, comm, status);
167 }
168
169#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
170 template<typename DATATYPE>
171 void AMPI_Bcast_wrap_adj(typename DATATYPE::AdjointType* &sendbufAdjoints, int sendbufSize, typename DATATYPE::AdjointType* &recvbufAdjoints, int recvbufSize, int count, DATATYPE* datatype, int root, AMPI_Comm comm) {
172 MEDI_UNUSED(count);
173
174 MPI_Gather(recvbufAdjoints, recvbufSize, datatype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufSize, datatype->getADTool().getAdjointMpiType(), root, comm);
175 }
176#endif
177
178#if MEDI_MPI_VERSION_3_0 <= MEDI_MPI_TARGET
179 template<typename DATATYPE>
180 void AMPI_Ibcast_wrap_adj(typename DATATYPE::AdjointType* &sendbufAdjoints, int sendbufSize, typename DATATYPE::AdjointType* &recvbufAdjoints, int recvbufSize, int count, DATATYPE* datatype, int root, AMPI_Comm comm, AMPI_Request* request) {
181 MEDI_UNUSED(count);
182
183 MPI_Igather(recvbufAdjoints, recvbufSize, datatype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufSize, datatype->getADTool().getAdjointMpiType(), root, comm, &request->request);
184 }
185#endif
186
187#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
188 template<typename SENDTYPE, typename RECVTYPE>
189 void AMPI_Scatter_adj(typename SENDTYPE::AdjointType* &sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE* sendtype, typename RECVTYPE::AdjointType* &recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE* recvtype, int root, AMPI_Comm comm) {
190 MEDI_UNUSED(sendcount);
191 MEDI_UNUSED(recvcount);
192
193 MPI_Gather(recvbufAdjoints, recvbufSize, recvtype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufSize, sendtype->getADTool().getAdjointMpiType(), root, comm);
194 }
195#endif
196
197#if MEDI_MPI_VERSION_3_0 <= MEDI_MPI_TARGET
198 template<typename SENDTYPE, typename RECVTYPE>
199 void AMPI_Iscatter_adj(typename SENDTYPE::AdjointType* &sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE* sendtype, typename RECVTYPE::AdjointType* &recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE* recvtype, int root, AMPI_Comm comm, AMPI_Request* request) {
200 MEDI_UNUSED(sendcount);
201 MEDI_UNUSED(recvcount);
202
203 MPI_Igather(recvbufAdjoints, recvbufSize, recvtype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufSize, sendtype->getADTool().getAdjointMpiType(), root, comm, &request->request);
204 }
205#endif
206
207#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
208 template<typename SENDTYPE, typename RECVTYPE>
209 void AMPI_Scatterv_adj(typename SENDTYPE::AdjointType* &sendbufAdjoints, int* sendbufCounts, MEDI_OPTIONAL_CONST int* sendbufDispl, MEDI_OPTIONAL_CONST int* sendcount, MEDI_OPTIONAL_CONST int* displs, SENDTYPE* sendtype, typename RECVTYPE::AdjointType* &recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE* recvtype, int root, AMPI_Comm comm) {
210 MEDI_UNUSED(sendcount);
211 MEDI_UNUSED(recvcount);
212
213 MPI_Gatherv(recvbufAdjoints, recvbufSize, recvtype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufCounts, sendbufDispl, sendtype->getADTool().getAdjointMpiType(), root, comm);
214 }
215#endif
216
217#if MEDI_MPI_VERSION_3_0 <= MEDI_MPI_TARGET
218 template<typename SENDTYPE, typename RECVTYPE>
219 void AMPI_Iscatterv_adj(typename SENDTYPE::AdjointType* &sendbufAdjoints, int* sendbufCounts, MEDI_OPTIONAL_CONST int* sendbufDispl, MEDI_OPTIONAL_CONST int* sendcount, MEDI_OPTIONAL_CONST int* displs, SENDTYPE* sendtype, typename RECVTYPE::AdjointType* &recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE* recvtype, int root, AMPI_Comm comm, AMPI_Request* request) {
220 MEDI_UNUSED(sendcount);
221 MEDI_UNUSED(recvcount);
222
223 MPI_Igatherv(recvbufAdjoints, recvbufSize, recvtype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufCounts, sendbufDispl, sendtype->getADTool().getAdjointMpiType(), root, comm, &request->request);
224 }
225#endif
226
227#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
228 template<typename SENDTYPE, typename RECVTYPE>
229 void AMPI_Gather_adj(typename SENDTYPE::AdjointType* &sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE* sendtype, typename RECVTYPE::AdjointType* &recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE* recvtype, int root, AMPI_Comm comm) {
230 MEDI_UNUSED(sendcount);
231 MEDI_UNUSED(recvcount);
232
233 MPI_Scatter(recvbufAdjoints, recvbufSize, recvtype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufSize, sendtype->getADTool().getAdjointMpiType(), root, comm);
234 }
235#endif
236
237#if MEDI_MPI_VERSION_3_0 <= MEDI_MPI_TARGET
238 template<typename SENDTYPE, typename RECVTYPE>
239 void AMPI_Igather_adj(typename SENDTYPE::AdjointType* &sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE* sendtype, typename RECVTYPE::AdjointType* &recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE* recvtype, int root, AMPI_Comm comm, AMPI_Request* request) {
240 MEDI_UNUSED(sendcount);
241 MEDI_UNUSED(recvcount);
242
243 MPI_Iscatter(recvbufAdjoints, recvbufSize, recvtype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufSize, sendtype->getADTool().getAdjointMpiType(), root, comm, &request->request);
244 }
245#endif
246
247#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
248 template<typename SENDTYPE, typename RECVTYPE>
249 void AMPI_Gatherv_adj(typename SENDTYPE::AdjointType* &sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE* sendtype, typename RECVTYPE::AdjointType* &recvbufAdjoints, int* recvbufCounts, MEDI_OPTIONAL_CONST int* recvbufDispls, MEDI_OPTIONAL_CONST int* recvcounts, MEDI_OPTIONAL_CONST int* displs, RECVTYPE* recvtype, int root, AMPI_Comm comm) {
250 MEDI_UNUSED(sendcount);
251 MEDI_UNUSED(recvcounts);
252 MEDI_UNUSED(displs);
253
254 MPI_Scatterv(recvbufAdjoints, recvbufCounts, recvbufDispls, recvtype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufSize, sendtype->getADTool().getAdjointMpiType(), root, comm);
255 }
256#endif
257
258#if MEDI_MPI_VERSION_3_0 <= MEDI_MPI_TARGET
259 template<typename SENDTYPE, typename RECVTYPE>
260 void AMPI_Igatherv_adj(typename SENDTYPE::AdjointType* &sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE* sendtype, typename RECVTYPE::AdjointType* &recvbufAdjoints, int* recvbufCounts, MEDI_OPTIONAL_CONST int* recvbufDispls, MEDI_OPTIONAL_CONST int* recvcounts, MEDI_OPTIONAL_CONST int* displs, RECVTYPE* recvtype, int root, AMPI_Comm comm, AMPI_Request* request) {
261 MEDI_UNUSED(sendcount);
262 MEDI_UNUSED(recvcounts);
263 MEDI_UNUSED(displs);
264
265 MPI_Iscatterv(recvbufAdjoints, recvbufCounts, recvbufDispls, recvtype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufSize, sendtype->getADTool().getAdjointMpiType(), root, comm, &request->request);
266 }
267#endif
268
269#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
270 template<typename SENDTYPE, typename RECVTYPE>
271 void AMPI_Allgather_adj(typename SENDTYPE::AdjointType* &sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE* sendtype, typename RECVTYPE::AdjointType* &recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE* recvtype, AMPI_Comm comm) {
272 MEDI_UNUSED(sendcount);
273 MEDI_UNUSED(recvcount);
274
275 MPI_Alltoall(recvbufAdjoints, recvbufSize, recvtype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufSize, sendtype->getADTool().getAdjointMpiType(), comm);
276 }
277#endif
278
279#if MEDI_MPI_VERSION_3_0 <= MEDI_MPI_TARGET
280 template<typename SENDTYPE, typename RECVTYPE>
281 void AMPI_Iallgather_adj(typename SENDTYPE::AdjointType* &sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE* sendtype, typename RECVTYPE::AdjointType* &recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE* recvtype, AMPI_Comm comm, AMPI_Request* request ) {
282 MEDI_UNUSED(sendcount);
283 MEDI_UNUSED(recvcount);
284
285 MPI_Ialltoall(recvbufAdjoints, recvbufSize, recvtype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufSize, sendtype->getADTool().getAdjointMpiType(), comm, &request->request);
286 }
287#endif
288
289#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
290 template<typename SENDTYPE, typename RECVTYPE>
291 void AMPI_Allgatherv_adj(typename SENDTYPE::AdjointType* &sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE* sendtype, typename RECVTYPE::AdjointType* &recvbufAdjoints, int* recvbufCounts, MEDI_OPTIONAL_CONST int* recvbufDispls, MEDI_OPTIONAL_CONST int* recvcounts, MEDI_OPTIONAL_CONST int* displs, RECVTYPE* recvtype, AMPI_Comm comm) {
292 MEDI_UNUSED(sendcount);
293 MEDI_UNUSED(recvcounts);
294 MEDI_UNUSED(displs);
295
296 LinearDisplacements linDis(getCommSize(comm), sendbufSize);
297
298 MPI_Alltoallv(recvbufAdjoints, recvbufCounts, recvbufDispls, recvtype->getADTool().getAdjointMpiType(), sendbufAdjoints, linDis.counts, linDis.displs, sendtype->getADTool().getAdjointMpiType(), comm);
299 }
300#endif
301
302#if MEDI_MPI_VERSION_3_0 <= MEDI_MPI_TARGET
303 template<typename SENDTYPE, typename RECVTYPE>
304 void AMPI_Iallgatherv_adj(typename SENDTYPE::AdjointType* &sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE* sendtype, typename RECVTYPE::AdjointType* &recvbufAdjoints, int* recvbufCounts, MEDI_OPTIONAL_CONST int* recvbufDispls, MEDI_OPTIONAL_CONST int* recvcounts, MEDI_OPTIONAL_CONST int* displs, RECVTYPE* recvtype, AMPI_Comm comm, AMPI_Request* request) {
305 MEDI_UNUSED(sendcount);
306 MEDI_UNUSED(recvcounts);
307 MEDI_UNUSED(displs);
308
309 LinearDisplacements* linDis = new LinearDisplacements(getCommSize(comm), sendbufSize);
310 request->setReverseData(reinterpret_cast<void*>(linDis), LinearDisplacements::deleteFunc);
311
312 MPI_Ialltoallv(recvbufAdjoints, recvbufCounts, recvbufDispls, recvtype->getADTool().getAdjointMpiType(), sendbufAdjoints, linDis->counts, linDis->displs, sendtype->getADTool().getAdjointMpiType(), comm, &request->request);
313 }
314#endif
315
316#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
317 template<typename SENDTYPE, typename RECVTYPE>
318 void AMPI_Alltoall_adj(typename SENDTYPE::AdjointType* &sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE* sendtype, typename RECVTYPE::AdjointType* &recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE* recvtype, AMPI_Comm comm) {
319 MEDI_UNUSED(sendcount);
320 MEDI_UNUSED(recvcount);
321
322 MPI_Alltoall(recvbufAdjoints, recvbufSize, recvtype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufSize, sendtype->getADTool().getAdjointMpiType(), comm);
323 }
324#endif
325
326#if MEDI_MPI_VERSION_3_0 <= MEDI_MPI_TARGET
327 template<typename SENDTYPE, typename RECVTYPE>
328 void AMPI_Ialltoall_adj(typename SENDTYPE::AdjointType* &sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE* sendtype, typename RECVTYPE::AdjointType* &recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE* recvtype, AMPI_Comm comm, AMPI_Request* request) {
329 MEDI_UNUSED(sendcount);
330 MEDI_UNUSED(recvcount);
331
332 MPI_Ialltoall(recvbufAdjoints, recvbufSize, recvtype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufSize, sendtype->getADTool().getAdjointMpiType(), comm, &request->request);
333 }
334#endif
335
336#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
337 template<typename SENDTYPE, typename RECVTYPE>
338 void AMPI_Alltoallv_adj(typename SENDTYPE::AdjointType* &sendbufAdjoints, int* sendbufCounts, MEDI_OPTIONAL_CONST int* sendbufDispls, MEDI_OPTIONAL_CONST int* sendcounts, MEDI_OPTIONAL_CONST int* sdispls, SENDTYPE* sendtype, typename RECVTYPE::AdjointType* &recvbufAdjoints, int* recvbufCounts, MEDI_OPTIONAL_CONST int* recvbufDispls, MEDI_OPTIONAL_CONST int* recvcounts, MEDI_OPTIONAL_CONST int* rdispls, RECVTYPE* recvtype, AMPI_Comm comm) {
339 MEDI_UNUSED(sendcounts);
340 MEDI_UNUSED(sdispls);
341 MEDI_UNUSED(recvcounts);
342 MEDI_UNUSED(rdispls);
343
344 MPI_Alltoallv(recvbufAdjoints, recvbufCounts, recvbufDispls, recvtype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufCounts, sendbufDispls, sendtype->getADTool().getAdjointMpiType(), comm);
345 }
346#endif
347
348#if MEDI_MPI_VERSION_3_0 <= MEDI_MPI_TARGET
349 template<typename SENDTYPE, typename RECVTYPE>
350 void AMPI_Ialltoallv_adj(typename SENDTYPE::AdjointType* &sendbufAdjoints, int* sendbufCounts, MEDI_OPTIONAL_CONST int* sendbufDispls, MEDI_OPTIONAL_CONST int* sendcounts, MEDI_OPTIONAL_CONST int* sdispls, SENDTYPE* sendtype, typename RECVTYPE::AdjointType* &recvbufAdjoints, int* recvbufCounts, MEDI_OPTIONAL_CONST int* recvbufDispls, MEDI_OPTIONAL_CONST int* recvcounts, MEDI_OPTIONAL_CONST int* rdispls, RECVTYPE* recvtype, AMPI_Comm comm, AMPI_Request* request) {
351 MEDI_UNUSED(sendcounts);
352 MEDI_UNUSED(sdispls);
353 MEDI_UNUSED(recvcounts);
354 MEDI_UNUSED(rdispls);
355
356 MPI_Ialltoallv(recvbufAdjoints, recvbufCounts, recvbufDispls, recvtype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufCounts, sendbufDispls, sendtype->getADTool().getAdjointMpiType(), comm, &request->request);
357 }
358#endif
359
360#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
361 template<typename DATATYPE>
362 void AMPI_Reduce_global_adj(typename DATATYPE::AdjointType* &sendbufAdjoints, int sendbufSize, typename DATATYPE::AdjointType* &recvbufAdjoints, int recvbufSize, int count, DATATYPE* datatype, AMPI_Op op, int root, AMPI_Comm comm) {
363 MEDI_UNUSED(count);
364 if(root == getCommRank(comm)) {
365 MPI_Bcast(recvbufAdjoints, recvbufSize, datatype->getADTool().getAdjointMpiType(), root, comm);
366 std::swap(sendbufAdjoints, recvbufAdjoints);
367 } else {
368 MPI_Bcast(sendbufAdjoints, sendbufSize, datatype->getADTool().getAdjointMpiType(), root, comm);
369 }
370 }
371#endif
372
373#if MEDI_MPI_VERSION_3_0 <= MEDI_MPI_TARGET
374 template<typename DATATYPE>
375 void AMPI_Ireduce_global_adj(typename DATATYPE::AdjointType* &sendbufAdjoints, int sendbufSize, typename DATATYPE::AdjointType* &recvbufAdjoints, int recvbufSize, int count, DATATYPE* datatype, AMPI_Op op, int root, AMPI_Comm comm, AMPI_Request* request) {
376 MEDI_UNUSED(count);
377 if(root == getCommRank(comm)) {
378 MPI_Ibcast(recvbufAdjoints, recvbufSize, datatype->getADTool().getAdjointMpiType(), root, comm, &request->request);
379 std::swap(sendbufAdjoints, recvbufAdjoints);
380 } else {
381 MPI_Ibcast(sendbufAdjoints, sendbufSize, datatype->getADTool().getAdjointMpiType(), root, comm, &request->request);
382 }
383 }
384#endif
385
386#if MEDI_MPI_VERSION_1_0 <= MEDI_MPI_TARGET
387 template<typename DATATYPE>
388 void AMPI_Allreduce_global_adj(typename DATATYPE::AdjointType* &sendbufAdjoints, int sendbufSize, typename DATATYPE::AdjointType* &recvbufAdjoints, int recvbufSize, int count, DATATYPE* datatype, AMPI_Op op, AMPI_Comm comm) {
389 MEDI_UNUSED(op);
390 MEDI_UNUSED(count);
391
392 MPI_Allgather(recvbufAdjoints, recvbufSize, datatype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufSize, datatype->getADTool().getAdjointMpiType(), comm);
393 }
394#endif
395
396#if MEDI_MPI_VERSION_3_0 <= MEDI_MPI_TARGET
397 template<typename DATATYPE>
398 void AMPI_Iallreduce_global_adj(typename DATATYPE::AdjointType* &sendbufAdjoints, int sendbufSize, typename DATATYPE::AdjointType* &recvbufAdjoints, int recvbufSize, int count, DATATYPE* datatype, AMPI_Op op, AMPI_Comm comm, AMPI_Request* request) {
399 MEDI_UNUSED(op);
400 MEDI_UNUSED(count);
401
402 MPI_Iallgather(recvbufAdjoints, recvbufSize, datatype->getADTool().getAdjointMpiType(), sendbufAdjoints, sendbufSize, datatype->getADTool().getAdjointMpiType(), comm, &request->request);
403 }
404#endif
405}
#define AMPI_Status
Definition ampiDefinitions.h:772
#define AMPI_Comm
Definition ampiDefinitions.h:778
#define MEDI_EXCEPTION(...)
Generates an exception.
Definition exceptions.hpp:44
#define MEDI_OPTIONAL_CONST
Definition macros.h:70
#define MEDI_UNUSED(name)
Definition macros.h:108
Global namespace for MeDiPack - Message Differentiation Package.
Definition adjointInterface.hpp:37
void AMPI_Recv_adj(typename DATATYPE::AdjointType *bufAdjoints, int bufSize, int count, DATATYPE *datatype, int src, int tag, MPI_Comm comm, MPI_Status *status, RecvAdjCall reverse_call)
Definition reverseFunctions.hpp:108
void AMPI_Irsend_adj(typename DATATYPE::AdjointType *bufAdjoints, int bufSize, int count, DATATYPE *datatype, int dest, int tag, MPI_Comm comm, AMPI_Request *request)
Definition reverseFunctions.hpp:100
void AMPI_Mrecv_adj(typename DATATYPE::AdjointType *bufAdjoints, int bufSize, int count, DATATYPE *datatype, AMPI_Message *message, MPI_Status *status, RecvAdjCall reverse_call)
Definition reverseFunctions.hpp:127
void AMPI_Irecv_adj(typename DATATYPE::AdjointType *bufAdjoints, int bufSize, int count, DATATYPE *datatype, int src, int tag, MPI_Comm comm, AMPI_Request *request, IrecvAdjCall reverse_call)
Definition reverseFunctions.hpp:134
void AMPI_Ireduce_global_adj(typename DATATYPE::AdjointType *&sendbufAdjoints, int sendbufSize, typename DATATYPE::AdjointType *&recvbufAdjoints, int recvbufSize, int count, DATATYPE *datatype, AMPI_Op op, int root, MPI_Comm comm, AMPI_Request *request)
Definition reverseFunctions.hpp:375
void AMPI_Iallgather_adj(typename SENDTYPE::AdjointType *&sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE *sendtype, typename RECVTYPE::AdjointType *&recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE *recvtype, MPI_Comm comm, AMPI_Request *request)
Definition reverseFunctions.hpp:281
void AMPI_Bcast_wrap_adj(typename DATATYPE::AdjointType *&sendbufAdjoints, int sendbufSize, typename DATATYPE::AdjointType *&recvbufAdjoints, int recvbufSize, int count, DATATYPE *datatype, int root, MPI_Comm comm)
Definition reverseFunctions.hpp:171
void AMPI_Igatherv_adj(typename SENDTYPE::AdjointType *&sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE *sendtype, typename RECVTYPE::AdjointType *&recvbufAdjoints, int *recvbufCounts, const int *recvbufDispls, const int *recvcounts, const int *displs, RECVTYPE *recvtype, int root, MPI_Comm comm, AMPI_Request *request)
Definition reverseFunctions.hpp:260
void AMPI_Iscatter_adj(typename SENDTYPE::AdjointType *&sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE *sendtype, typename RECVTYPE::AdjointType *&recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE *recvtype, int root, MPI_Comm comm, AMPI_Request *request)
Definition reverseFunctions.hpp:199
IrecvAdjCall
Definition enums.hpp:40
void AMPI_Gatherv_adj(typename SENDTYPE::AdjointType *&sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE *sendtype, typename RECVTYPE::AdjointType *&recvbufAdjoints, int *recvbufCounts, const int *recvbufDispls, const int *recvcounts, const int *displs, RECVTYPE *recvtype, int root, MPI_Comm comm)
Definition reverseFunctions.hpp:249
RecvAdjCall
Definition enums.hpp:44
void AMPI_Gather_adj(typename SENDTYPE::AdjointType *&sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE *sendtype, typename RECVTYPE::AdjointType *&recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE *recvtype, int root, MPI_Comm comm)
Definition reverseFunctions.hpp:229
void AMPI_Send_adj(typename DATATYPE::AdjointType *bufAdjoints, int bufSize, int count, DATATYPE *datatype, int dest, int tag, MPI_Comm comm)
Definition reverseFunctions.hpp:44
void AMPI_Ibcast_wrap_adj(typename DATATYPE::AdjointType *&sendbufAdjoints, int sendbufSize, typename DATATYPE::AdjointType *&recvbufAdjoints, int recvbufSize, int count, DATATYPE *datatype, int root, MPI_Comm comm, AMPI_Request *request)
Definition reverseFunctions.hpp:180
void AMPI_Iscatterv_adj(typename SENDTYPE::AdjointType *&sendbufAdjoints, int *sendbufCounts, const int *sendbufDispl, const int *sendcount, const int *displs, SENDTYPE *sendtype, typename RECVTYPE::AdjointType *&recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE *recvtype, int root, MPI_Comm comm, AMPI_Request *request)
Definition reverseFunctions.hpp:219
void AMPI_Scatterv_adj(typename SENDTYPE::AdjointType *&sendbufAdjoints, int *sendbufCounts, const int *sendbufDispl, const int *sendcount, const int *displs, SENDTYPE *sendtype, typename RECVTYPE::AdjointType *&recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE *recvtype, int root, MPI_Comm comm)
Definition reverseFunctions.hpp:209
void AMPI_Allgatherv_adj(typename SENDTYPE::AdjointType *&sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE *sendtype, typename RECVTYPE::AdjointType *&recvbufAdjoints, int *recvbufCounts, const int *recvbufDispls, const int *recvcounts, const int *displs, RECVTYPE *recvtype, MPI_Comm comm)
Definition reverseFunctions.hpp:291
void AMPI_Rsend_adj(typename DATATYPE::AdjointType *bufAdjoints, int bufSize, int count, DATATYPE *datatype, int dest, int tag, MPI_Comm comm)
Definition reverseFunctions.hpp:92
void AMPI_Allgather_adj(typename SENDTYPE::AdjointType *&sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE *sendtype, typename RECVTYPE::AdjointType *&recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE *recvtype, MPI_Comm comm)
Definition reverseFunctions.hpp:271
void AMPI_Ialltoall_adj(typename SENDTYPE::AdjointType *&sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE *sendtype, typename RECVTYPE::AdjointType *&recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE *recvtype, MPI_Comm comm, AMPI_Request *request)
Definition reverseFunctions.hpp:328
void AMPI_Allreduce_global_adj(typename DATATYPE::AdjointType *&sendbufAdjoints, int sendbufSize, typename DATATYPE::AdjointType *&recvbufAdjoints, int recvbufSize, int count, DATATYPE *datatype, AMPI_Op op, MPI_Comm comm)
Definition reverseFunctions.hpp:388
void AMPI_Isend_adj(typename DATATYPE::AdjointType *bufAdjoints, int bufSize, int count, DATATYPE *datatype, int dest, int tag, MPI_Comm comm, AMPI_Request *request)
Definition reverseFunctions.hpp:52
int getCommSize(MPI_Comm comm)
Helper function that gets the number of ranks from the communicator.
Definition mpiTools.h:64
void AMPI_Sendrecv_adj(typename SENDTYPE::AdjointType *sendbuf, int sendbufSize, int sendcount, SENDTYPE *sendtype, int dest, int sendtag, typename RECVTYPE::AdjointType *recvbuf, int recvbufSize, int recvcount, RECVTYPE *recvtype, int source, int recvtag, MPI_Comm comm, MPI_Status *status)
Definition reverseFunctions.hpp:160
int getCommRank(MPI_Comm comm)
Helper function that gets the own rank number from the communicator.
Definition mpiTools.h:52
void AMPI_Issend_adj(typename DATATYPE::AdjointType *bufAdjoints, int bufSize, int count, DATATYPE *datatype, int dest, int tag, MPI_Comm comm, AMPI_Request *request)
Definition reverseFunctions.hpp:84
void AMPI_Bsend_adj(typename DATATYPE::AdjointType *bufAdjoints, int bufSize, int count, DATATYPE *datatype, int dest, int tag, MPI_Comm comm)
Definition reverseFunctions.hpp:60
void AMPI_Iallgatherv_adj(typename SENDTYPE::AdjointType *&sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE *sendtype, typename RECVTYPE::AdjointType *&recvbufAdjoints, int *recvbufCounts, const int *recvbufDispls, const int *recvcounts, const int *displs, RECVTYPE *recvtype, MPI_Comm comm, AMPI_Request *request)
Definition reverseFunctions.hpp:304
void AMPI_Alltoall_adj(typename SENDTYPE::AdjointType *&sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE *sendtype, typename RECVTYPE::AdjointType *&recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE *recvtype, MPI_Comm comm)
Definition reverseFunctions.hpp:318
void AMPI_Igather_adj(typename SENDTYPE::AdjointType *&sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE *sendtype, typename RECVTYPE::AdjointType *&recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE *recvtype, int root, MPI_Comm comm, AMPI_Request *request)
Definition reverseFunctions.hpp:239
void AMPI_Ibsend_adj(typename DATATYPE::AdjointType *bufAdjoints, int bufSize, int count, DATATYPE *datatype, int dest, int tag, MPI_Comm comm, AMPI_Request *request)
Definition reverseFunctions.hpp:68
void AMPI_Iallreduce_global_adj(typename DATATYPE::AdjointType *&sendbufAdjoints, int sendbufSize, typename DATATYPE::AdjointType *&recvbufAdjoints, int recvbufSize, int count, DATATYPE *datatype, AMPI_Op op, MPI_Comm comm, AMPI_Request *request)
Definition reverseFunctions.hpp:398
void AMPI_Ssend_adj(typename DATATYPE::AdjointType *bufAdjoints, int bufSize, int count, DATATYPE *datatype, int dest, int tag, MPI_Comm comm)
Definition reverseFunctions.hpp:76
void AMPI_Ialltoallv_adj(typename SENDTYPE::AdjointType *&sendbufAdjoints, int *sendbufCounts, const int *sendbufDispls, const int *sendcounts, const int *sdispls, SENDTYPE *sendtype, typename RECVTYPE::AdjointType *&recvbufAdjoints, int *recvbufCounts, const int *recvbufDispls, const int *recvcounts, const int *rdispls, RECVTYPE *recvtype, MPI_Comm comm, AMPI_Request *request)
Definition reverseFunctions.hpp:350
void AMPI_Alltoallv_adj(typename SENDTYPE::AdjointType *&sendbufAdjoints, int *sendbufCounts, const int *sendbufDispls, const int *sendcounts, const int *sdispls, SENDTYPE *sendtype, typename RECVTYPE::AdjointType *&recvbufAdjoints, int *recvbufCounts, const int *recvbufDispls, const int *recvcounts, const int *rdispls, RECVTYPE *recvtype, MPI_Comm comm)
Definition reverseFunctions.hpp:338
void AMPI_Reduce_global_adj(typename DATATYPE::AdjointType *&sendbufAdjoints, int sendbufSize, typename DATATYPE::AdjointType *&recvbufAdjoints, int recvbufSize, int count, DATATYPE *datatype, AMPI_Op op, int root, MPI_Comm comm)
Definition reverseFunctions.hpp:362
void AMPI_Imrecv_adj(typename DATATYPE::AdjointType *bufAdjoints, int bufSize, int count, DATATYPE *datatype, AMPI_Message *message, AMPI_Request *request, IrecvAdjCall reverse_call)
Definition reverseFunctions.hpp:152
void AMPI_Scatter_adj(typename SENDTYPE::AdjointType *&sendbufAdjoints, int sendbufSize, int sendcount, SENDTYPE *sendtype, typename RECVTYPE::AdjointType *&recvbufAdjoints, int recvbufSize, int recvcount, RECVTYPE *recvtype, int root, MPI_Comm comm)
Definition reverseFunctions.hpp:189
Stores additional information for a MPI_Message.
Definition message.hpp:44
int src
Definition message.hpp:48
MPI_Comm comm
Definition message.hpp:49
int tag
Definition message.hpp:47
Structure for the special handling of the MPI_Op structure.
Definition op.hpp:50
Definition async.hpp:45
MPI_Request request
Definition async.hpp:46
void setReverseData(void *data, DeleteReverseData func)
Definition async.hpp:70
Helper structure for the easy creation of a linear displacement with a the same length.
Definition displacementTools.hpp:43
int * displs
The array with the displacements.
Definition displacementTools.hpp:48
int * counts
The array with the counts.
Definition displacementTools.hpp:53
static void deleteFunc(void *d)
Helper function for deleting a LinearDisplacements structure.
Definition displacementTools.hpp:90