Logo AND Algorithmique Numérique Distribuée

Public GIT Repository
Update copyright lines with new year.
[simgrid.git] / src / smpi / bindings / smpi_pmpi_coll.cpp
1 /* Copyright (c) 2007-2020. The SimGrid Team. All rights reserved.          */
2
3 /* This program is free software; you can redistribute it and/or modify it
4  * under the terms of the license (GNU LGPL) which comes with this package. */
5
6 #include "private.hpp"
7 #include "smpi_coll.hpp"
8 #include "smpi_comm.hpp"
9 #include "smpi_request.hpp"
10 #include "smpi_datatype_derived.hpp"
11 #include "smpi_op.hpp"
12 #include "src/smpi/include/smpi_actor.hpp"
13
14 XBT_LOG_EXTERNAL_DEFAULT_CATEGORY(smpi_pmpi);
15
16   static const void* smpi_get_in_place_buf(const void* inplacebuf, const void* otherbuf,std::unique_ptr<unsigned char[]>& tmp_sendbuf, int count, MPI_Datatype datatype){
17   if (inplacebuf == MPI_IN_PLACE) {
18       tmp_sendbuf.reset(new unsigned char[count * datatype->get_extent()]);
19       simgrid::smpi::Datatype::copy(otherbuf, count, datatype, tmp_sendbuf.get(), count, datatype);
20     return tmp_sendbuf.get();
21   }else{
22     return inplacebuf;
23   }
24 }
25 /* PMPI User level calls */
26
27 int PMPI_Barrier(MPI_Comm comm)
28 {
29   return PMPI_Ibarrier(comm, MPI_REQUEST_IGNORED);
30 }
31
32 int PMPI_Ibarrier(MPI_Comm comm, MPI_Request *request)
33 {
34   CHECK_COMM(1)
35   CHECK_REQUEST(2)
36
37   smpi_bench_end();
38   int rank = simgrid::s4u::this_actor::get_pid();
39   TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Barrier" : "PMPI_Ibarrier",
40                      new simgrid::instr::NoOpTIData(request == MPI_REQUEST_IGNORED ? "barrier" : "ibarrier"));
41   if (request == MPI_REQUEST_IGNORED) {
42     simgrid::smpi::colls::barrier(comm);
43     // Barrier can be used to synchronize RMA calls. Finish all requests from comm before.
44     comm->finish_rma_calls();
45   } else
46     simgrid::smpi::colls::ibarrier(comm, request);
47
48   TRACE_smpi_comm_out(rank);
49   smpi_bench_begin();
50   return MPI_SUCCESS;
51 }
52
53 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
54 {
55   return PMPI_Ibcast(buf, count, datatype, root, comm, MPI_REQUEST_IGNORED);
56 }
57
58 int PMPI_Ibcast(void *buf, int count, MPI_Datatype datatype, 
59                    int root, MPI_Comm comm, MPI_Request* request)
60 {
61   CHECK_COMM(5)
62   CHECK_BUFFER(1, buf, count)
63   CHECK_COUNT(2, count)
64   CHECK_TYPE(3, datatype)
65   CHECK_ROOT(4)
66   CHECK_REQUEST(6)
67
68   smpi_bench_end();
69   int rank = simgrid::s4u::this_actor::get_pid();
70   TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Bcast" : "PMPI_Ibcast",
71                      new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "bcast" : "ibcast", root, -1.0,
72                                                     datatype->is_replayable() ? count : count * datatype->size(), -1,
73                                                     simgrid::smpi::Datatype::encode(datatype), ""));
74   if (comm->size() > 1) {
75     if (request == MPI_REQUEST_IGNORED)
76       simgrid::smpi::colls::bcast(buf, count, datatype, root, comm);
77     else
78       simgrid::smpi::colls::ibcast(buf, count, datatype, root, comm, request);
79   } else {
80     if (request != MPI_REQUEST_IGNORED)
81       *request = MPI_REQUEST_NULL;
82   }
83
84   TRACE_smpi_comm_out(rank);
85   smpi_bench_begin();
86   return MPI_SUCCESS;
87 }
88
89 int PMPI_Gather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
90                 int root, MPI_Comm comm){
91   return PMPI_Igather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
92 }
93
94 int PMPI_Igather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
95                  MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
96 {
97   CHECK_COMM(8)
98   if(sendbuf != MPI_IN_PLACE){
99     CHECK_BUFFER(1,sendbuf, sendcount)
100     CHECK_COUNT(2, sendcount)
101     CHECK_TYPE(3, sendtype)
102   }
103   if(comm->rank() == root){
104     CHECK_TYPE(6, recvtype)
105     CHECK_COUNT(5, recvcount)
106     CHECK_BUFFER(4, recvbuf, recvcount)
107   }
108   CHECK_ROOT(7)
109   CHECK_REQUEST(9)
110
111   smpi_bench_end();
112   const void* real_sendbuf   = sendbuf;
113   int real_sendcount         = sendcount;
114   MPI_Datatype real_sendtype = sendtype;
115   if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
116     real_sendcount = 0;
117     real_sendtype  = recvtype;
118   }
119   int rank = simgrid::s4u::this_actor::get_pid();
120
121   TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gather" : "PMPI_Igather",
122                      new simgrid::instr::CollTIData(
123                          request == MPI_REQUEST_IGNORED ? "gather" : "igather", root, -1.0,
124                          real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
125                          (comm->rank() != root || recvtype->is_replayable()) ? recvcount : recvcount * recvtype->size(),
126                          simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
127   if (request == MPI_REQUEST_IGNORED)
128     simgrid::smpi::colls::gather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm);
129   else
130     simgrid::smpi::colls::igather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm,
131                                   request);
132
133   TRACE_smpi_comm_out(rank);
134   smpi_bench_begin();
135   return MPI_SUCCESS;
136 }
137
138 int PMPI_Gatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, const int *recvcounts, const int *displs,
139                 MPI_Datatype recvtype, int root, MPI_Comm comm){
140   return PMPI_Igatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, root, comm, MPI_REQUEST_IGNORED);
141 }
142
143 int PMPI_Igatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
144                   MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
145 {
146   CHECK_COMM(9)
147   CHECK_BUFFER(1, sendbuf, sendcount)
148   if(sendbuf != MPI_IN_PLACE){
149     CHECK_TYPE(3, sendtype)
150     CHECK_COUNT(2, sendcount)
151   }
152   if(comm->rank() == root){
153     CHECK_TYPE(6, recvtype)
154     CHECK_NULL(5, MPI_ERR_COUNT, recvcounts)
155     CHECK_NULL(6, MPI_ERR_ARG, displs)
156   }
157   CHECK_ROOT(8)
158   CHECK_REQUEST(10)
159
160   if (comm->rank() == root){
161     for (int i = 0; i < comm->size(); i++) {
162       CHECK_COUNT(5, recvcounts[i])
163       CHECK_BUFFER(4,recvbuf,recvcounts[i])
164     }
165   }
166
167   smpi_bench_end();
168   const void* real_sendbuf   = sendbuf;
169   int real_sendcount         = sendcount;
170   MPI_Datatype real_sendtype = sendtype;
171   if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
172     real_sendcount = 0;
173     real_sendtype  = recvtype;
174   }
175
176   int rank         = simgrid::s4u::this_actor::get_pid();
177   int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
178
179   std::vector<int>* trace_recvcounts = new std::vector<int>();
180   if (comm->rank() == root) {
181     for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
182       trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
183   }
184
185   TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gatherv" : "PMPI_Igatherv",
186                      new simgrid::instr::VarCollTIData(
187                          request == MPI_REQUEST_IGNORED ? "gatherv" : "igatherv", root,
188                          real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
189                          nullptr, dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(real_sendtype),
190                          simgrid::smpi::Datatype::encode(recvtype)));
191   if (request == MPI_REQUEST_IGNORED)
192     simgrid::smpi::colls::gatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
193                                   root, comm);
194   else
195     simgrid::smpi::colls::igatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
196                                    root, comm, request);
197
198   TRACE_smpi_comm_out(rank);
199   smpi_bench_begin();
200   return MPI_SUCCESS;
201 }
202
203 int PMPI_Allgather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
204                    void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm){
205   return PMPI_Iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
206 }
207
208 int PMPI_Iallgather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
209                     MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
210 {
211   CHECK_COMM(7)
212   CHECK_BUFFER(1, sendbuf, sendcount)
213   CHECK_BUFFER(4, recvbuf, recvcount)
214   if(sendbuf != MPI_IN_PLACE){
215     CHECK_COUNT(2, sendcount)
216     CHECK_TYPE(3, sendtype)
217   }
218   CHECK_TYPE(6, recvtype)
219   CHECK_COUNT(5, recvcount)
220   CHECK_REQUEST(8)
221
222   smpi_bench_end();
223   if (sendbuf == MPI_IN_PLACE) {
224     sendbuf   = static_cast<char*>(recvbuf) + recvtype->get_extent() * recvcount * comm->rank();
225     sendcount = recvcount;
226     sendtype  = recvtype;
227   }
228   int rank = simgrid::s4u::this_actor::get_pid();
229
230   TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgather" : "PMPI_Iallggather",
231                      new simgrid::instr::CollTIData(
232                          request == MPI_REQUEST_IGNORED ? "allgather" : "iallgather", -1, -1.0,
233                          sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
234                          recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
235                          simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
236   if (request == MPI_REQUEST_IGNORED)
237     simgrid::smpi::colls::allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
238   else
239     simgrid::smpi::colls::iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, request);
240
241   TRACE_smpi_comm_out(rank);
242   smpi_bench_begin();
243   return MPI_SUCCESS;
244 }
245
246 int PMPI_Allgatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
247                    void *recvbuf, const int *recvcounts, const int *displs, MPI_Datatype recvtype, MPI_Comm comm){
248   return PMPI_Iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, MPI_REQUEST_IGNORED);
249 }
250
251 int PMPI_Iallgatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
252                      MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
253 {
254   CHECK_COMM(8)
255   CHECK_BUFFER(1, sendbuf, sendcount)
256   if(sendbuf != MPI_IN_PLACE)
257     CHECK_TYPE(3, sendtype)
258   CHECK_TYPE(6, recvtype)
259   CHECK_NULL(5, MPI_ERR_COUNT, recvcounts)
260   CHECK_NULL(6, MPI_ERR_ARG, displs)
261   if(sendbuf != MPI_IN_PLACE)
262     CHECK_COUNT(2, sendcount)
263   CHECK_REQUEST(9)
264   for (int i = 0; i < comm->size(); i++) {
265     CHECK_COUNT(5, recvcounts[i])
266     CHECK_BUFFER(4, recvbuf, recvcounts[i])
267   }
268
269   smpi_bench_end();
270   if (sendbuf == MPI_IN_PLACE) {
271     sendbuf   = static_cast<char*>(recvbuf) + recvtype->get_extent() * displs[comm->rank()];
272     sendcount = recvcounts[comm->rank()];
273     sendtype  = recvtype;
274   }
275   int rank         = simgrid::s4u::this_actor::get_pid();
276   int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
277
278   std::vector<int>* trace_recvcounts = new std::vector<int>();
279   for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
280     trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
281   }
282
283   TRACE_smpi_comm_in(
284       rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgatherv" : "PMPI_Iallgatherv",
285       new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "allgatherv" : "iallgatherv", -1,
286                                         sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(), nullptr,
287                                         dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
288                                         simgrid::smpi::Datatype::encode(recvtype)));
289   if (request == MPI_REQUEST_IGNORED)
290     simgrid::smpi::colls::allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
291   else
292     simgrid::smpi::colls::iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm,
293                                       request);
294
295   TRACE_smpi_comm_out(rank);
296   smpi_bench_begin();
297   return MPI_SUCCESS;
298 }
299
300 int PMPI_Scatter(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
301                 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
302   return PMPI_Iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
303 }
304
305 int PMPI_Iscatter(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
306                   MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
307 {
308   CHECK_COMM(8)
309   if(comm->rank() == root){
310     CHECK_BUFFER(1, sendbuf, sendcount)
311     CHECK_COUNT(2, sendcount)
312     CHECK_TYPE(3, sendtype)
313   }
314   if(recvbuf != MPI_IN_PLACE){
315     CHECK_BUFFER(4, recvbuf, recvcount)
316     CHECK_COUNT(5, recvcount)
317     CHECK_TYPE(6, recvtype)
318   }
319   CHECK_ROOT(8)
320   CHECK_REQUEST(9)
321
322   smpi_bench_end();
323   if (recvbuf == MPI_IN_PLACE) {
324     recvtype  = sendtype;
325     recvcount = sendcount;
326   }
327   int rank = simgrid::s4u::this_actor::get_pid();
328
329   TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatter" : "PMPI_Iscatter",
330                      new simgrid::instr::CollTIData(
331                          request == MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
332                          (comm->rank() != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
333                          recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
334                          simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
335   if (request == MPI_REQUEST_IGNORED)
336     simgrid::smpi::colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
337   else
338     simgrid::smpi::colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
339
340   TRACE_smpi_comm_out(rank);
341   smpi_bench_begin();
342   return MPI_SUCCESS;
343 }
344
345 int PMPI_Scatterv(const void *sendbuf, const int *sendcounts, const int *displs,
346                  MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
347   return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
348 }
349
350 int PMPI_Iscatterv(const void* sendbuf, const int* sendcounts, const int* displs, MPI_Datatype sendtype, void* recvbuf, int recvcount,
351                    MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
352 {
353   CHECK_COMM(9)
354   if(recvbuf != MPI_IN_PLACE){
355     CHECK_BUFFER(4, recvbuf, recvcount)
356     CHECK_COUNT(5, recvcount)
357     CHECK_TYPE(7, recvtype)
358   }
359   CHECK_ROOT(9)
360   CHECK_REQUEST(10)
361   if (comm->rank() == root) {
362     CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
363     CHECK_NULL(3, MPI_ERR_ARG, displs)
364     CHECK_TYPE(4, sendtype)
365     for (int i = 0; i < comm->size(); i++){
366       CHECK_BUFFER(1, sendbuf, sendcounts[i])
367       CHECK_COUNT(2, sendcounts[i])
368     }
369     if (recvbuf == MPI_IN_PLACE) {
370       recvtype  = sendtype;
371       recvcount = sendcounts[comm->rank()];
372     }
373   }
374
375   smpi_bench_end();
376
377   int rank         = simgrid::s4u::this_actor::get_pid();
378   int dt_size_send = sendtype->is_replayable() ? 1 : sendtype->size();
379
380   std::vector<int>* trace_sendcounts = new std::vector<int>();
381   if (comm->rank() == root) {
382     for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
383       trace_sendcounts->push_back(sendcounts[i] * dt_size_send);
384     }
385   }
386
387   TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatterv" : "PMPI_Iscatterv",
388                      new simgrid::instr::VarCollTIData(
389                          request == MPI_REQUEST_IGNORED ? "scatterv" : "iscatterv", root, dt_size_send,
390                          trace_sendcounts, recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
391                          nullptr, simgrid::smpi::Datatype::encode(sendtype),
392                          simgrid::smpi::Datatype::encode(recvtype)));
393   if (request == MPI_REQUEST_IGNORED)
394     simgrid::smpi::colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
395   else
396     simgrid::smpi::colls::iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm,
397                                     request);
398
399   TRACE_smpi_comm_out(rank);
400   smpi_bench_begin();
401   return MPI_SUCCESS;
402 }
403
404 int PMPI_Reduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
405 {
406   return PMPI_Ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, MPI_REQUEST_IGNORED);
407 }
408
409 int PMPI_Ireduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm, MPI_Request* request)
410 {
411   CHECK_COMM(7)
412   CHECK_BUFFER(1, sendbuf, count)
413   if(comm->rank() == root)
414     CHECK_BUFFER(5, recvbuf, count)
415   CHECK_TYPE(4, datatype)
416   CHECK_COUNT(3, count)
417   CHECK_OP(5)
418   CHECK_ROOT(7)
419   CHECK_REQUEST(8)
420
421   smpi_bench_end();
422   int rank = simgrid::s4u::this_actor::get_pid();
423
424   TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce" : "PMPI_Ireduce",
425                      new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "reduce" : "ireduce", root, 0,
426                                                     datatype->is_replayable() ? count : count * datatype->size(), -1,
427                                                     simgrid::smpi::Datatype::encode(datatype), ""));
428   if (request == MPI_REQUEST_IGNORED)
429     simgrid::smpi::colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
430   else
431     simgrid::smpi::colls::ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, request);
432
433   TRACE_smpi_comm_out(rank);
434   smpi_bench_begin();
435   return MPI_SUCCESS;
436 }
437
438 int PMPI_Reduce_local(const void* inbuf, void* inoutbuf, int count, MPI_Datatype datatype, MPI_Op op)
439 {
440   CHECK_BUFFER(1, inbuf, count)
441   CHECK_BUFFER(2, inoutbuf, count)
442   CHECK_TYPE(4, datatype)
443   CHECK_COUNT(3, count)
444   CHECK_OP(5)
445
446   smpi_bench_end();
447   op->apply(inbuf, inoutbuf, &count, datatype);
448   smpi_bench_begin();
449   return MPI_SUCCESS;
450 }
451
452 int PMPI_Allreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
453 {
454   return PMPI_Iallreduce(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
455 }
456
457 int PMPI_Iallreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
458 {
459   CHECK_COMM(6)
460   CHECK_BUFFER(1, sendbuf, count)
461   CHECK_BUFFER(2, recvbuf, count)
462   CHECK_TYPE(4, datatype)
463   CHECK_COUNT(3, count)
464   CHECK_REQUEST(7)
465   CHECK_OP(5)
466
467   smpi_bench_end();
468   std::unique_ptr<unsigned char[]> tmp_sendbuf;
469   const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
470
471   int rank = simgrid::s4u::this_actor::get_pid();
472
473   TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allreduce" : "PMPI_Iallreduce",
474                      new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "allreduce" : "iallreduce", -1, 0,
475                                                     datatype->is_replayable() ? count : count * datatype->size(), -1,
476                                                     simgrid::smpi::Datatype::encode(datatype), ""));
477
478   if (request == MPI_REQUEST_IGNORED)
479     simgrid::smpi::colls::allreduce(real_sendbuf, recvbuf, count, datatype, op, comm);
480   else
481     simgrid::smpi::colls::iallreduce(real_sendbuf, recvbuf, count, datatype, op, comm, request);
482
483   TRACE_smpi_comm_out(rank);
484   smpi_bench_begin();
485   return MPI_SUCCESS;
486 }
487
488 int PMPI_Scan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
489 {
490   return PMPI_Iscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
491 }
492
493 int PMPI_Iscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request)
494 {
495   CHECK_COMM(6)
496   CHECK_BUFFER(1,sendbuf,count)
497   CHECK_BUFFER(2,recvbuf,count)
498   CHECK_TYPE(4, datatype)
499   CHECK_COUNT(3, count)
500   CHECK_REQUEST(7)
501   CHECK_OP(5)
502
503   smpi_bench_end();
504   int rank         = simgrid::s4u::this_actor::get_pid();
505   std::unique_ptr<unsigned char[]> tmp_sendbuf;
506   const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
507
508   TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scan" : "PMPI_Iscan",
509                      new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "scan" : "iscan", -1,
510                                                      datatype->is_replayable() ? count : count * datatype->size(),
511                                                      simgrid::smpi::Datatype::encode(datatype)));
512
513   int retval;
514   if (request == MPI_REQUEST_IGNORED)
515     retval = simgrid::smpi::colls::scan(real_sendbuf, recvbuf, count, datatype, op, comm);
516   else
517     retval = simgrid::smpi::colls::iscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
518
519   TRACE_smpi_comm_out(rank);
520   smpi_bench_begin();
521   return retval;
522 }
523
524 int PMPI_Exscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
525 {
526   return PMPI_Iexscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
527 }
528
529 int PMPI_Iexscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request){
530   CHECK_COMM(6)
531   CHECK_BUFFER(1, sendbuf, count)
532   CHECK_BUFFER(2, recvbuf, count)
533   CHECK_TYPE(4, datatype)
534   CHECK_COUNT(3, count)
535   CHECK_REQUEST(7)
536   CHECK_OP(5)
537
538   smpi_bench_end();
539   int rank         = simgrid::s4u::this_actor::get_pid();
540   std::unique_ptr<unsigned char[]> tmp_sendbuf;
541   const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
542
543   TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Exscan" : "PMPI_Iexscan",
544                      new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "exscan" : "iexscan", -1,
545                                                      datatype->is_replayable() ? count : count * datatype->size(),
546                                                      simgrid::smpi::Datatype::encode(datatype)));
547
548   int retval;
549   if (request == MPI_REQUEST_IGNORED)
550     retval = simgrid::smpi::colls::exscan(real_sendbuf, recvbuf, count, datatype, op, comm);
551   else
552     retval = simgrid::smpi::colls::iexscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
553
554   TRACE_smpi_comm_out(rank);
555   smpi_bench_begin();
556   return retval;
557 }
558
559 int PMPI_Reduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
560 {
561   return PMPI_Ireduce_scatter(sendbuf, recvbuf, recvcounts, datatype, op, comm, MPI_REQUEST_IGNORED);
562 }
563
564 int PMPI_Ireduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
565 {
566   CHECK_COMM(6)
567   CHECK_TYPE(4, datatype)
568   CHECK_NULL(3, MPI_ERR_COUNT, recvcounts)
569   CHECK_REQUEST(7)
570   CHECK_OP(5)
571   for (int i = 0; i < comm->size(); i++) {
572     CHECK_COUNT(3, recvcounts[i])
573     CHECK_BUFFER(1, sendbuf, recvcounts[i])
574     CHECK_BUFFER(2, recvbuf, recvcounts[i])
575   }
576
577   smpi_bench_end();
578   int rank                           = simgrid::s4u::this_actor::get_pid();
579   std::vector<int>* trace_recvcounts = new std::vector<int>();
580   int dt_send_size                   = datatype->is_replayable() ? 1 : datatype->size();
581   int totalcount                     = 0;
582
583   for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
584     trace_recvcounts->push_back(recvcounts[i] * dt_send_size);
585     totalcount += recvcounts[i];
586   }
587   std::unique_ptr<unsigned char[]> tmp_sendbuf;
588   const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, totalcount, datatype);
589
590   TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter" : "PMPI_Ireduce_scatter",
591                      new simgrid::instr::VarCollTIData(
592                          request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, dt_send_size, nullptr,
593                          -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
594
595   if (request == MPI_REQUEST_IGNORED)
596     simgrid::smpi::colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm);
597   else
598     simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm, request);
599
600   TRACE_smpi_comm_out(rank);
601   smpi_bench_begin();
602   return MPI_SUCCESS;
603 }
604
605 int PMPI_Reduce_scatter_block(const void *sendbuf, void *recvbuf, int recvcount,
606                               MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
607 {
608   return PMPI_Ireduce_scatter_block(sendbuf, recvbuf, recvcount, datatype, op, comm, MPI_REQUEST_IGNORED);
609 }
610
611 int PMPI_Ireduce_scatter_block(const void* sendbuf, void* recvbuf, int recvcount, MPI_Datatype datatype, MPI_Op op,
612                                MPI_Comm comm, MPI_Request* request)
613 {
614   CHECK_COMM(6)
615   CHECK_BUFFER(1, sendbuf, recvcount)
616   CHECK_BUFFER(2, recvbuf, recvcount)
617   CHECK_TYPE(4, datatype)
618   CHECK_COUNT(3, recvcount)
619   CHECK_REQUEST(7)
620   CHECK_OP(5)
621
622   smpi_bench_end();
623   int count = comm->size();
624
625   int rank                           = simgrid::s4u::this_actor::get_pid();
626   int dt_send_size                   = datatype->is_replayable() ? 1 : datatype->size();
627   std::vector<int>* trace_recvcounts = new std::vector<int>(recvcount * dt_send_size); // copy data to avoid bad free
628   std::unique_ptr<unsigned char[]> tmp_sendbuf;
629   const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, recvcount * count, datatype);
630
631   TRACE_smpi_comm_in(
632       rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter_block" : "PMPI_Ireduce_scatter_block",
633       new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, 0,
634                                         nullptr, -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
635
636   int* recvcounts = new int[count];
637   for (int i      = 0; i < count; i++)
638     recvcounts[i] = recvcount;
639   if (request == MPI_REQUEST_IGNORED)
640     simgrid::smpi::colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm);
641   else
642     simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm, request);
643   delete[] recvcounts;
644
645   TRACE_smpi_comm_out(rank);
646   smpi_bench_begin();
647   return MPI_SUCCESS;
648 }
649
650 int PMPI_Alltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
651                   MPI_Datatype recvtype, MPI_Comm comm){
652   return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
653 }
654
655 int PMPI_Ialltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
656                    MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
657 {
658   CHECK_COMM(7)
659   CHECK_BUFFER(1, sendbuf, sendcount)
660   CHECK_BUFFER(4, recvbuf, recvcount)
661   if(sendbuf != MPI_IN_PLACE)
662     CHECK_TYPE(3, sendtype)
663   CHECK_TYPE(6, recvtype)
664   CHECK_COUNT(5, recvcount)
665   if(sendbuf != MPI_IN_PLACE)
666     CHECK_COUNT(2, sendcount)
667   CHECK_COUNT(5, recvcount)
668   CHECK_REQUEST(8)
669
670   smpi_bench_end();
671   int rank                 = simgrid::s4u::this_actor::get_pid();
672   int real_sendcount         = sendcount;
673   MPI_Datatype real_sendtype = sendtype;
674   
675   std::unique_ptr<unsigned char[]> tmp_sendbuf;
676   const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, recvcount * comm->size(), recvtype);
677
678   if (sendbuf == MPI_IN_PLACE) {
679     real_sendcount = recvcount;
680     real_sendtype  = recvtype;
681   }
682
683   TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoall" : "PMPI_Ialltoall",
684                      new simgrid::instr::CollTIData(
685                          request == MPI_REQUEST_IGNORED ? "alltoall" : "ialltoall", -1, -1.0,
686                          real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
687                          recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
688                          simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
689   int retval;
690   if (request == MPI_REQUEST_IGNORED)
691     retval =
692         simgrid::smpi::colls::alltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, comm);
693   else
694     retval = simgrid::smpi::colls::ialltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype,
695                                              comm, request);
696
697   TRACE_smpi_comm_out(rank);
698   smpi_bench_begin();
699   return retval;
700 }
701
702 int PMPI_Alltoallv(const void* sendbuf, const int* sendcounts, const int* senddispls, MPI_Datatype sendtype, void* recvbuf,
703                    const int* recvcounts, const int* recvdispls, MPI_Datatype recvtype, MPI_Comm comm)
704 {
705   return PMPI_Ialltoallv(sendbuf, sendcounts, senddispls, sendtype, recvbuf, recvcounts, recvdispls, recvtype, comm, MPI_REQUEST_IGNORED);
706 }
707
708 int PMPI_Ialltoallv(const void* sendbuf, const int* sendcounts, const int* senddispls, MPI_Datatype sendtype, void* recvbuf,
709                     const int* recvcounts, const int* recvdispls, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
710 {
711   CHECK_COMM(9)
712   if(sendbuf != MPI_IN_PLACE){
713     CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
714     CHECK_NULL(3, MPI_ERR_ARG, senddispls)
715     CHECK_TYPE(4, sendtype)
716   }
717   CHECK_TYPE(8, recvtype)
718   CHECK_NULL(6, MPI_ERR_COUNT, recvcounts)
719   CHECK_NULL(7, MPI_ERR_ARG, recvdispls)
720   CHECK_REQUEST(10)
721
722   int rank = simgrid::s4u::this_actor::get_pid();
723   int size = comm->size();
724   for (int i = 0; i < size; i++) {
725     if(sendbuf != MPI_IN_PLACE){
726       CHECK_BUFFER(1, sendbuf, sendcounts[i])
727       CHECK_COUNT(2, sendcounts[i])
728     }
729     CHECK_BUFFER(5, recvbuf, recvcounts[i])
730     CHECK_COUNT(6, recvcounts[i])
731   }
732
733   smpi_bench_end();
734   int send_size                      = 0;
735   int recv_size                      = 0;
736   std::vector<int>* trace_sendcounts = new std::vector<int>();
737   std::vector<int>* trace_recvcounts = new std::vector<int>();
738   int dt_size_recv                   = recvtype->size();
739
740   const int* real_sendcounts = sendcounts;
741   const int* real_senddispls  = senddispls;
742   MPI_Datatype real_sendtype = sendtype;
743   int maxsize              = 0;
744   for (int i = 0; i < size; i++) { // copy data to avoid bad free
745     recv_size += recvcounts[i] * dt_size_recv;
746     trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
747     if (((recvdispls[i] + recvcounts[i]) * dt_size_recv) > maxsize)
748       maxsize = (recvdispls[i] + recvcounts[i]) * dt_size_recv;
749   }
750
751   std::unique_ptr<unsigned char[]> tmp_sendbuf;
752   std::unique_ptr<int[]> tmp_sendcounts;
753   std::unique_ptr<int[]> tmp_senddispls;
754   const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
755   if (sendbuf == MPI_IN_PLACE) {
756     tmp_sendcounts.reset(new int[size]);
757     std::copy(recvcounts, recvcounts + size, tmp_sendcounts.get());
758     real_sendcounts = tmp_sendcounts.get();
759     tmp_senddispls.reset(new int[size]);
760     std::copy(recvdispls, recvdispls + size, tmp_senddispls.get());
761     real_senddispls = tmp_senddispls.get();
762     real_sendtype  = recvtype;
763   }
764
765   int dt_size_send = real_sendtype->size();
766
767   for (int i = 0; i < size; i++) { // copy data to avoid bad free
768     send_size += real_sendcounts[i] * dt_size_send;
769     trace_sendcounts->push_back(real_sendcounts[i] * dt_size_send);
770   }
771
772   TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallv" : "PMPI_Ialltoallv",
773                      new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
774                                                        send_size, trace_sendcounts, recv_size, trace_recvcounts,
775                                                        simgrid::smpi::Datatype::encode(real_sendtype),
776                                                        simgrid::smpi::Datatype::encode(recvtype)));
777
778   int retval;
779   if (request == MPI_REQUEST_IGNORED)
780     retval = simgrid::smpi::colls::alltoallv(real_sendbuf, real_sendcounts, real_senddispls, real_sendtype, recvbuf,
781                                              recvcounts, recvdispls, recvtype, comm);
782   else
783     retval = simgrid::smpi::colls::ialltoallv(real_sendbuf, real_sendcounts, real_senddispls, real_sendtype, recvbuf,
784                                               recvcounts, recvdispls, recvtype, comm, request);
785
786   TRACE_smpi_comm_out(rank);
787   smpi_bench_begin();
788   return retval;
789 }
790
791 int PMPI_Alltoallw(const void* sendbuf, const int* sendcounts, const int* senddispls, const MPI_Datatype* sendtypes, void* recvbuf,
792                    const int* recvcounts, const int* recvdispls, const MPI_Datatype* recvtypes, MPI_Comm comm)
793 {
794   return PMPI_Ialltoallw(sendbuf, sendcounts, senddispls, sendtypes, recvbuf, recvcounts, recvdispls, recvtypes, comm, MPI_REQUEST_IGNORED);
795 }
796
797 int PMPI_Ialltoallw(const void* sendbuf, const int* sendcounts, const int* senddispls, const MPI_Datatype* sendtypes, void* recvbuf,
798                     const int* recvcounts, const int* recvdispls, const MPI_Datatype* recvtypes, MPI_Comm comm, MPI_Request* request)
799 {
800   CHECK_COMM(9)
801   if(sendbuf != MPI_IN_PLACE){
802     CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
803     CHECK_NULL(3, MPI_ERR_ARG, senddispls)
804     CHECK_NULL(4, MPI_ERR_TYPE, sendtypes)
805   }
806   CHECK_NULL(6, MPI_ERR_COUNT, recvcounts)
807   CHECK_NULL(7, MPI_ERR_ARG, recvdispls)
808   CHECK_NULL(8, MPI_ERR_TYPE, recvtypes)
809   CHECK_REQUEST(10)
810   int rank = simgrid::s4u::this_actor::get_pid();
811   int size = comm->size();
812   for (int i = 0; i < size; i++) {
813     if(sendbuf != MPI_IN_PLACE){
814       CHECK_BUFFER(1, sendbuf, sendcounts[i])
815       CHECK_COUNT(2, sendcounts[i])
816       CHECK_TYPE(4, sendtypes[i])
817     }
818     CHECK_BUFFER(5, recvbuf, recvcounts[i])
819     CHECK_COUNT(6, recvcounts[i])
820     CHECK_TYPE(8, recvtypes[i])
821   }
822
823   smpi_bench_end();
824
825   int send_size                      = 0;
826   int recv_size                      = 0;
827   std::vector<int>* trace_sendcounts = new std::vector<int>();
828   std::vector<int>* trace_recvcounts = new std::vector<int>();
829
830   const int* real_sendcounts         = sendcounts;
831   const int* real_senddispls          = senddispls;
832   const MPI_Datatype* real_sendtypes = sendtypes;
833   unsigned long maxsize      = 0;
834   for (int i = 0; i < size; i++) { // copy data to avoid bad free
835     if (recvtypes[i] == MPI_DATATYPE_NULL) {
836       delete trace_recvcounts;
837       delete trace_sendcounts;
838       return MPI_ERR_TYPE;
839     }
840     recv_size += recvcounts[i] * recvtypes[i]->size();
841     trace_recvcounts->push_back(recvcounts[i] * recvtypes[i]->size());
842     if ((recvdispls[i] + (recvcounts[i] * recvtypes[i]->size())) > maxsize)
843       maxsize = recvdispls[i] + (recvcounts[i] * recvtypes[i]->size());
844   }
845
846   std::unique_ptr<unsigned char[]> tmp_sendbuf;
847   std::unique_ptr<int[]> tmp_sendcounts;
848   std::unique_ptr<int[]> tmp_senddispls;
849   std::unique_ptr<MPI_Datatype[]> tmp_sendtypes;
850   const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
851   if (sendbuf == MPI_IN_PLACE) {
852     tmp_sendcounts.reset(new int[size]);
853     std::copy(recvcounts, recvcounts + size, tmp_sendcounts.get());
854     real_sendcounts = tmp_sendcounts.get();
855     tmp_senddispls.reset(new int[size]);
856     std::copy(recvdispls, recvdispls + size, tmp_senddispls.get());
857     real_senddispls = tmp_senddispls.get();
858     tmp_sendtypes.reset(new MPI_Datatype[size]);
859     std::copy(recvtypes, recvtypes + size, tmp_sendtypes.get());
860     real_sendtypes = tmp_sendtypes.get();
861   }
862
863   for (int i = 0; i < size; i++) { // copy data to avoid bad free
864     send_size += real_sendcounts[i] * real_sendtypes[i]->size();
865     trace_sendcounts->push_back(real_sendcounts[i] * real_sendtypes[i]->size());
866   }
867
868   TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallw" : "PMPI_Ialltoallw",
869                      new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
870                                                        send_size, trace_sendcounts, recv_size, trace_recvcounts,
871                                                        simgrid::smpi::Datatype::encode(real_sendtypes[0]),
872                                                        simgrid::smpi::Datatype::encode(recvtypes[0])));
873
874   int retval;
875   if (request == MPI_REQUEST_IGNORED)
876     retval = simgrid::smpi::colls::alltoallw(real_sendbuf, real_sendcounts, real_senddispls, real_sendtypes, recvbuf,
877                                              recvcounts, recvdispls, recvtypes, comm);
878   else
879     retval = simgrid::smpi::colls::ialltoallw(real_sendbuf, real_sendcounts, real_senddispls, real_sendtypes, recvbuf,
880                                               recvcounts, recvdispls, recvtypes, comm, request);
881
882   TRACE_smpi_comm_out(rank);
883   smpi_bench_begin();
884   return retval;
885 }