Logo AND Algorithmique Numérique Distribuée

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