Logo AND Algorithmique Numérique Distribuée

Public GIT Repository
Remove extern "C" from cpp files (src/smpi/).
[simgrid.git] / src / smpi / bindings / smpi_pmpi_coll.cpp
1 /* Copyright (c) 2007-2018. 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_datatype_derived.hpp"
10 #include "smpi_op.hpp"
11 #include "smpi_process.hpp"
12
13 XBT_LOG_EXTERNAL_DEFAULT_CATEGORY(smpi_pmpi);
14
15 /* PMPI User level calls */
16
17 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
18 {
19   int retval = 0;
20
21   smpi_bench_end();
22
23   if (comm == MPI_COMM_NULL) {
24     retval = MPI_ERR_COMM;
25   } else if (not datatype->is_valid()) {
26     retval = MPI_ERR_ARG;
27   } else {
28     int rank = simgrid::s4u::Actor::self()->getPid();
29     TRACE_smpi_comm_in(rank, __FUNCTION__,
30                        new simgrid::instr::CollTIData("bcast", root, -1.0,
31                                                       datatype->is_replayable() ? count : count * datatype->size(), -1,
32                                                       encode_datatype(datatype), ""));
33     if (comm->size() > 1)
34       simgrid::smpi::Colls::bcast(buf, count, datatype, root, comm);
35     retval = MPI_SUCCESS;
36
37     TRACE_smpi_comm_out(rank);
38   }
39   smpi_bench_begin();
40   return retval;
41 }
42
43 int PMPI_Barrier(MPI_Comm comm)
44 {
45   int retval = 0;
46
47   smpi_bench_end();
48
49   if (comm == MPI_COMM_NULL) {
50     retval = MPI_ERR_COMM;
51   } else {
52     int rank = simgrid::s4u::Actor::self()->getPid();
53     TRACE_smpi_comm_in(rank, __FUNCTION__, new simgrid::instr::NoOpTIData("barrier"));
54
55     simgrid::smpi::Colls::barrier(comm);
56
57     //Barrier can be used to synchronize RMA calls. Finish all requests from comm before.
58     comm->finish_rma_calls();
59
60     retval = MPI_SUCCESS;
61
62     TRACE_smpi_comm_out(rank);
63   }
64
65   smpi_bench_begin();
66   return retval;
67 }
68
69 int PMPI_Gather(void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
70                 int root, MPI_Comm comm)
71 {
72   int retval = 0;
73
74   smpi_bench_end();
75
76   if (comm == MPI_COMM_NULL) {
77     retval = MPI_ERR_COMM;
78   } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
79             ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL))){
80     retval = MPI_ERR_TYPE;
81   } else if ((( sendbuf != MPI_IN_PLACE) && (sendcount <0)) || ((comm->rank() == root) && (recvcount <0))){
82     retval = MPI_ERR_COUNT;
83   } else {
84
85     char* sendtmpbuf = static_cast<char*>(sendbuf);
86     int sendtmpcount = sendcount;
87     MPI_Datatype sendtmptype = sendtype;
88     if( (comm->rank() == root) && (sendbuf == MPI_IN_PLACE )) {
89       sendtmpcount=0;
90       sendtmptype=recvtype;
91     }
92     int rank = simgrid::s4u::Actor::self()->getPid();
93
94     TRACE_smpi_comm_in(rank, __FUNCTION__,
95                        new simgrid::instr::CollTIData(
96                            "gather", root, -1.0,
97                            sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(),
98                            (comm->rank() != root || recvtype->is_replayable()) ? recvcount : recvcount * recvtype->size(),
99                            encode_datatype(sendtmptype), encode_datatype(recvtype)));
100
101     simgrid::smpi::Colls::gather(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, root, comm);
102
103     retval = MPI_SUCCESS;
104     TRACE_smpi_comm_out(rank);
105   }
106
107   smpi_bench_begin();
108   return retval;
109 }
110
111 int PMPI_Gatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, int *recvcounts, int *displs,
112                 MPI_Datatype recvtype, int root, MPI_Comm comm)
113 {
114   int retval = 0;
115
116   smpi_bench_end();
117
118   if (comm == MPI_COMM_NULL) {
119     retval = MPI_ERR_COMM;
120   } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
121             ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL))){
122     retval = MPI_ERR_TYPE;
123   } else if (( sendbuf != MPI_IN_PLACE) && (sendcount <0)){
124     retval = MPI_ERR_COUNT;
125   } else if (recvcounts == nullptr || displs == nullptr) {
126     retval = MPI_ERR_ARG;
127   } else {
128     char* sendtmpbuf = static_cast<char*>(sendbuf);
129     int sendtmpcount = sendcount;
130     MPI_Datatype sendtmptype = sendtype;
131     if( (comm->rank() == root) && (sendbuf == MPI_IN_PLACE )) {
132       sendtmpcount=0;
133       sendtmptype=recvtype;
134     }
135
136     int rank         = simgrid::s4u::Actor::self()->getPid();
137     int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
138
139     std::vector<int>* trace_recvcounts = new std::vector<int>;
140     if (comm->rank() == root) {
141       for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
142         trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
143     }
144
145     TRACE_smpi_comm_in(rank, __FUNCTION__,
146                        new simgrid::instr::VarCollTIData(
147                            "gatherV", root,
148                            sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(), nullptr,
149                            dt_size_recv, trace_recvcounts, encode_datatype(sendtmptype), encode_datatype(recvtype)));
150
151     retval = simgrid::smpi::Colls::gatherv(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcounts, displs, recvtype, root, comm);
152     TRACE_smpi_comm_out(rank);
153   }
154
155   smpi_bench_begin();
156   return retval;
157 }
158
159 int PMPI_Allgather(void *sendbuf, int sendcount, MPI_Datatype sendtype,
160                    void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm)
161 {
162   int retval = MPI_SUCCESS;
163
164   smpi_bench_end();
165
166   if (comm == MPI_COMM_NULL) {
167     retval = MPI_ERR_COMM;
168   } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
169             (recvtype == MPI_DATATYPE_NULL)){
170     retval = MPI_ERR_TYPE;
171   } else if ((( sendbuf != MPI_IN_PLACE) && (sendcount <0)) ||
172             (recvcount <0)){
173     retval = MPI_ERR_COUNT;
174   } else {
175     if(sendbuf == MPI_IN_PLACE) {
176       sendbuf=static_cast<char*>(recvbuf)+recvtype->get_extent()*recvcount*comm->rank();
177       sendcount=recvcount;
178       sendtype=recvtype;
179     }
180     int rank = simgrid::s4u::Actor::self()->getPid();
181
182     TRACE_smpi_comm_in(rank, __FUNCTION__,
183                        new simgrid::instr::CollTIData("allGather", -1, -1.0,
184                                                       sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
185                                                       recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
186                                                       encode_datatype(sendtype), encode_datatype(recvtype)));
187
188     simgrid::smpi::Colls::allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
189     TRACE_smpi_comm_out(rank);
190   }
191   smpi_bench_begin();
192   return retval;
193 }
194
195 int PMPI_Allgatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype,
196                    void *recvbuf, int *recvcounts, int *displs, MPI_Datatype recvtype, MPI_Comm comm)
197 {
198   int retval = 0;
199
200   smpi_bench_end();
201
202   if (comm == MPI_COMM_NULL) {
203     retval = MPI_ERR_COMM;
204   } else if (((sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) || (recvtype == MPI_DATATYPE_NULL)) {
205     retval = MPI_ERR_TYPE;
206   } else if (( sendbuf != MPI_IN_PLACE) && (sendcount <0)){
207     retval = MPI_ERR_COUNT;
208   } else if (recvcounts == nullptr || displs == nullptr) {
209     retval = MPI_ERR_ARG;
210   } else {
211
212     if(sendbuf == MPI_IN_PLACE) {
213       sendbuf=static_cast<char*>(recvbuf)+recvtype->get_extent()*displs[comm->rank()];
214       sendcount=recvcounts[comm->rank()];
215       sendtype=recvtype;
216     }
217     int rank               = simgrid::s4u::Actor::self()->getPid();
218     int dt_size_recv       = recvtype->is_replayable() ? 1 : recvtype->size();
219
220     std::vector<int>* trace_recvcounts = new std::vector<int>;
221     for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
222       trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
223
224     TRACE_smpi_comm_in(rank, __FUNCTION__,
225                        new simgrid::instr::VarCollTIData(
226                            "allGatherV", -1, sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(), nullptr,
227                            dt_size_recv, trace_recvcounts, encode_datatype(sendtype), encode_datatype(recvtype)));
228
229     simgrid::smpi::Colls::allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
230     retval = MPI_SUCCESS;
231     TRACE_smpi_comm_out(rank);
232   }
233
234   smpi_bench_begin();
235   return retval;
236 }
237
238 int PMPI_Scatter(void *sendbuf, int sendcount, MPI_Datatype sendtype,
239                 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm)
240 {
241   int retval = 0;
242
243   smpi_bench_end();
244
245   if (comm == MPI_COMM_NULL) {
246     retval = MPI_ERR_COMM;
247   } else if (((comm->rank() == root) && (not sendtype->is_valid())) ||
248              ((recvbuf != MPI_IN_PLACE) && (not recvtype->is_valid()))) {
249     retval = MPI_ERR_TYPE;
250   } else if ((sendbuf == recvbuf) ||
251       ((comm->rank()==root) && sendcount>0 && (sendbuf == nullptr))){
252     retval = MPI_ERR_BUFFER;
253   }else {
254
255     if (recvbuf == MPI_IN_PLACE) {
256       recvtype  = sendtype;
257       recvcount = sendcount;
258     }
259     int rank = simgrid::s4u::Actor::self()->getPid();
260
261     TRACE_smpi_comm_in(rank, __FUNCTION__,
262                        new simgrid::instr::CollTIData(
263                            "scatter", root, -1.0,
264                            (comm->rank() != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
265                            recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(), encode_datatype(sendtype),
266                            encode_datatype(recvtype)));
267
268     simgrid::smpi::Colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
269     retval = MPI_SUCCESS;
270     TRACE_smpi_comm_out(rank);
271   }
272
273   smpi_bench_begin();
274   return retval;
275 }
276
277 int PMPI_Scatterv(void *sendbuf, int *sendcounts, int *displs,
278                  MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm)
279 {
280   int retval = 0;
281
282   smpi_bench_end();
283
284   if (comm == MPI_COMM_NULL) {
285     retval = MPI_ERR_COMM;
286   } else if (sendcounts == nullptr || displs == nullptr) {
287     retval = MPI_ERR_ARG;
288   } else if (((comm->rank() == root) && (sendtype == MPI_DATATYPE_NULL)) ||
289              ((recvbuf != MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL))) {
290     retval = MPI_ERR_TYPE;
291   } else {
292     if (recvbuf == MPI_IN_PLACE) {
293       recvtype  = sendtype;
294       recvcount = sendcounts[comm->rank()];
295     }
296     int rank               = simgrid::s4u::Actor::self()->getPid();
297     int dt_size_send       = sendtype->is_replayable() ? 1 : sendtype->size();
298
299     std::vector<int>* trace_sendcounts = new std::vector<int>;
300     if (comm->rank() == root) {
301       for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
302         trace_sendcounts->push_back(sendcounts[i] * dt_size_send);
303     }
304
305     TRACE_smpi_comm_in(rank, __FUNCTION__, new simgrid::instr::VarCollTIData(
306                                                "scatterV", root, dt_size_send, trace_sendcounts,
307                                                recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(), nullptr,
308                                                encode_datatype(sendtype), encode_datatype(recvtype)));
309
310     retval = simgrid::smpi::Colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
311
312     TRACE_smpi_comm_out(rank);
313   }
314
315   smpi_bench_begin();
316   return retval;
317 }
318
319 int PMPI_Reduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
320 {
321   int retval = 0;
322
323   smpi_bench_end();
324
325   if (comm == MPI_COMM_NULL) {
326     retval = MPI_ERR_COMM;
327   } else if (not datatype->is_valid() || op == MPI_OP_NULL) {
328     retval = MPI_ERR_ARG;
329   } else {
330     int rank = simgrid::s4u::Actor::self()->getPid();
331
332     TRACE_smpi_comm_in(rank, __FUNCTION__,
333                        new simgrid::instr::CollTIData("reduce", root, 0,
334                                                       datatype->is_replayable() ? count : count * datatype->size(), -1,
335                                                       encode_datatype(datatype), ""));
336
337     simgrid::smpi::Colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
338
339     retval = MPI_SUCCESS;
340     TRACE_smpi_comm_out(rank);
341   }
342
343   smpi_bench_begin();
344   return retval;
345 }
346
347 int PMPI_Reduce_local(void *inbuf, void *inoutbuf, int count, MPI_Datatype datatype, MPI_Op op){
348   int retval = 0;
349
350   smpi_bench_end();
351   if (not datatype->is_valid() || op == MPI_OP_NULL) {
352     retval = MPI_ERR_ARG;
353   } else {
354     op->apply(inbuf, inoutbuf, &count, datatype);
355     retval = MPI_SUCCESS;
356   }
357   smpi_bench_begin();
358   return retval;
359 }
360
361 int PMPI_Allreduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
362 {
363   int retval = 0;
364
365   smpi_bench_end();
366
367   if (comm == MPI_COMM_NULL) {
368     retval = MPI_ERR_COMM;
369   } else if (not datatype->is_valid()) {
370     retval = MPI_ERR_TYPE;
371   } else if (op == MPI_OP_NULL) {
372     retval = MPI_ERR_OP;
373   } else {
374
375     char* sendtmpbuf = static_cast<char*>(sendbuf);
376     if( sendbuf == MPI_IN_PLACE ) {
377       sendtmpbuf = static_cast<char*>(xbt_malloc(count*datatype->get_extent()));
378       simgrid::smpi::Datatype::copy(recvbuf, count, datatype,sendtmpbuf, count, datatype);
379     }
380     int rank = simgrid::s4u::Actor::self()->getPid();
381
382     TRACE_smpi_comm_in(rank, __FUNCTION__,
383                        new simgrid::instr::CollTIData("allReduce", -1, 0,
384                                                       datatype->is_replayable() ? count : count * datatype->size(), -1,
385                                                       encode_datatype(datatype), ""));
386
387     simgrid::smpi::Colls::allreduce(sendtmpbuf, recvbuf, count, datatype, op, comm);
388
389     if( sendbuf == MPI_IN_PLACE )
390       xbt_free(sendtmpbuf);
391
392     retval = MPI_SUCCESS;
393     TRACE_smpi_comm_out(rank);
394   }
395
396   smpi_bench_begin();
397   return retval;
398 }
399
400 int PMPI_Scan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
401 {
402   int retval = 0;
403
404   smpi_bench_end();
405
406   if (comm == MPI_COMM_NULL) {
407     retval = MPI_ERR_COMM;
408   } else if (not datatype->is_valid()) {
409     retval = MPI_ERR_TYPE;
410   } else if (op == MPI_OP_NULL) {
411     retval = MPI_ERR_OP;
412   } else {
413     int rank = simgrid::s4u::Actor::self()->getPid();
414
415     TRACE_smpi_comm_in(rank, __FUNCTION__, new simgrid::instr::Pt2PtTIData(
416                                                "scan", -1, datatype->is_replayable() ? count : count * datatype->size(),
417                                                encode_datatype(datatype)));
418
419     retval = simgrid::smpi::Colls::scan(sendbuf, recvbuf, count, datatype, op, comm);
420
421     TRACE_smpi_comm_out(rank);
422   }
423
424   smpi_bench_begin();
425   return retval;
426 }
427
428 int PMPI_Exscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm){
429   int retval = 0;
430
431   smpi_bench_end();
432
433   if (comm == MPI_COMM_NULL) {
434     retval = MPI_ERR_COMM;
435   } else if (not datatype->is_valid()) {
436     retval = MPI_ERR_TYPE;
437   } else if (op == MPI_OP_NULL) {
438     retval = MPI_ERR_OP;
439   } else {
440     int rank         = simgrid::s4u::Actor::self()->getPid();
441     void* sendtmpbuf = sendbuf;
442     if (sendbuf == MPI_IN_PLACE) {
443       sendtmpbuf = static_cast<void*>(xbt_malloc(count * datatype->size()));
444       memcpy(sendtmpbuf, recvbuf, count * datatype->size());
445     }
446
447     TRACE_smpi_comm_in(rank, __FUNCTION__, new simgrid::instr::Pt2PtTIData(
448                                                "exscan", -1, datatype->is_replayable() ? count : count * datatype->size(),
449                                                encode_datatype(datatype)));
450
451     retval = simgrid::smpi::Colls::exscan(sendtmpbuf, recvbuf, count, datatype, op, comm);
452
453     TRACE_smpi_comm_out(rank);
454     if (sendbuf == MPI_IN_PLACE)
455       xbt_free(sendtmpbuf);
456   }
457
458   smpi_bench_begin();
459   return retval;
460 }
461
462 int PMPI_Reduce_scatter(void *sendbuf, void *recvbuf, int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
463 {
464   int retval = 0;
465   smpi_bench_end();
466
467   if (comm == MPI_COMM_NULL) {
468     retval = MPI_ERR_COMM;
469   } else if (not datatype->is_valid()) {
470     retval = MPI_ERR_TYPE;
471   } else if (op == MPI_OP_NULL) {
472     retval = MPI_ERR_OP;
473   } else if (recvcounts == nullptr) {
474     retval = MPI_ERR_ARG;
475   } else {
476     int rank                           = simgrid::s4u::Actor::self()->getPid();
477     std::vector<int>* trace_recvcounts = new std::vector<int>;
478     int dt_send_size                   = datatype->is_replayable() ? 1 : datatype->size();
479     int totalcount    = 0;
480
481     for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
482       trace_recvcounts->push_back(recvcounts[i] * dt_send_size);
483       totalcount += recvcounts[i];
484     }
485
486     void* sendtmpbuf = sendbuf;
487     if (sendbuf == MPI_IN_PLACE) {
488       sendtmpbuf = static_cast<void*>(xbt_malloc(totalcount * datatype->size()));
489       memcpy(sendtmpbuf, recvbuf, totalcount * datatype->size());
490     }
491
492     TRACE_smpi_comm_in(rank, __FUNCTION__,
493                        new simgrid::instr::VarCollTIData("reduceScatter", -1, dt_send_size, nullptr, -1,
494                                                          trace_recvcounts, encode_datatype(datatype), ""));
495
496     simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
497     retval = MPI_SUCCESS;
498     TRACE_smpi_comm_out(rank);
499
500     if (sendbuf == MPI_IN_PLACE)
501       xbt_free(sendtmpbuf);
502   }
503
504   smpi_bench_begin();
505   return retval;
506 }
507
508 int PMPI_Reduce_scatter_block(void *sendbuf, void *recvbuf, int recvcount,
509                               MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
510 {
511   int retval;
512   smpi_bench_end();
513
514   if (comm == MPI_COMM_NULL) {
515     retval = MPI_ERR_COMM;
516   } else if (not datatype->is_valid()) {
517     retval = MPI_ERR_TYPE;
518   } else if (op == MPI_OP_NULL) {
519     retval = MPI_ERR_OP;
520   } else if (recvcount < 0) {
521     retval = MPI_ERR_ARG;
522   } else {
523     int count = comm->size();
524
525     int rank                           = simgrid::s4u::Actor::self()->getPid();
526     int dt_send_size                   = datatype->is_replayable() ? 1 : datatype->size();
527     std::vector<int>* trace_recvcounts = new std::vector<int>(recvcount * dt_send_size); // copy data to avoid bad free
528
529     void* sendtmpbuf       = sendbuf;
530     if (sendbuf == MPI_IN_PLACE) {
531       sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * count * datatype->size()));
532       memcpy(sendtmpbuf, recvbuf, recvcount * count * datatype->size());
533     }
534
535     TRACE_smpi_comm_in(rank, __FUNCTION__,
536                        new simgrid::instr::VarCollTIData("reduceScatter", -1, 0, nullptr, -1, trace_recvcounts,
537                                                          encode_datatype(datatype), ""));
538
539     int* recvcounts = new int[count];
540     for (int i      = 0; i < count; i++)
541       recvcounts[i] = recvcount;
542     simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
543     delete[] recvcounts;
544     retval = MPI_SUCCESS;
545
546     TRACE_smpi_comm_out(rank);
547
548     if (sendbuf == MPI_IN_PLACE)
549       xbt_free(sendtmpbuf);
550   }
551
552   smpi_bench_begin();
553   return retval;
554 }
555
556 int PMPI_Alltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
557                   MPI_Datatype recvtype, MPI_Comm comm)
558 {
559   int retval = 0;
560   smpi_bench_end();
561
562   if (comm == MPI_COMM_NULL) {
563     retval = MPI_ERR_COMM;
564   } else if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL) {
565     retval = MPI_ERR_TYPE;
566   } else {
567     int rank                 = simgrid::s4u::Actor::self()->getPid();
568     void* sendtmpbuf         = static_cast<char*>(sendbuf);
569     int sendtmpcount         = sendcount;
570     MPI_Datatype sendtmptype = sendtype;
571     if (sendbuf == MPI_IN_PLACE) {
572       sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * comm->size() * recvtype->size()));
573       memcpy(sendtmpbuf, recvbuf, recvcount * comm->size() * recvtype->size());
574       sendtmpcount = recvcount;
575       sendtmptype  = recvtype;
576     }
577
578     TRACE_smpi_comm_in(
579         rank, __FUNCTION__,
580         new simgrid::instr::CollTIData("allToAll", -1, -1.0,
581                                        sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(),
582                                        recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
583                                        encode_datatype(sendtmptype), encode_datatype(recvtype)));
584
585     retval = simgrid::smpi::Colls::alltoall(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, comm);
586
587     TRACE_smpi_comm_out(rank);
588
589     if (sendbuf == MPI_IN_PLACE)
590       xbt_free(sendtmpbuf);
591   }
592
593   smpi_bench_begin();
594   return retval;
595 }
596
597 int PMPI_Alltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf,
598                    int* recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm)
599 {
600   int retval = 0;
601
602   smpi_bench_end();
603
604   if (comm == MPI_COMM_NULL) {
605     retval = MPI_ERR_COMM;
606   } else if (sendtype == MPI_DATATYPE_NULL || recvtype == MPI_DATATYPE_NULL) {
607     retval = MPI_ERR_TYPE;
608   } else if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
609              recvdisps == nullptr) {
610     retval = MPI_ERR_ARG;
611   } else {
612     int rank                           = simgrid::s4u::Actor::self()->getPid();
613     int size               = comm->size();
614     int send_size                      = 0;
615     int recv_size                      = 0;
616     std::vector<int>* trace_sendcounts = new std::vector<int>;
617     std::vector<int>* trace_recvcounts = new std::vector<int>;
618     int dt_size_recv       = recvtype->size();
619
620     void* sendtmpbuf         = static_cast<char*>(sendbuf);
621     int* sendtmpcounts       = sendcounts;
622     int* sendtmpdisps        = senddisps;
623     MPI_Datatype sendtmptype = sendtype;
624     int maxsize              = 0;
625     for (int i = 0; i < size; i++) { // copy data to avoid bad free
626       recv_size += recvcounts[i] * dt_size_recv;
627       trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
628       if (((recvdisps[i] + recvcounts[i]) * dt_size_recv) > maxsize)
629         maxsize = (recvdisps[i] + recvcounts[i]) * dt_size_recv;
630     }
631
632     if (sendbuf == MPI_IN_PLACE) {
633       sendtmpbuf = static_cast<void*>(xbt_malloc(maxsize));
634       memcpy(sendtmpbuf, recvbuf, maxsize);
635       sendtmpcounts = static_cast<int*>(xbt_malloc(size * sizeof(int)));
636       memcpy(sendtmpcounts, recvcounts, size * sizeof(int));
637       sendtmpdisps = static_cast<int*>(xbt_malloc(size * sizeof(int)));
638       memcpy(sendtmpdisps, recvdisps, size * sizeof(int));
639       sendtmptype = recvtype;
640     }
641
642     int dt_size_send = sendtmptype->size();
643
644     for (int i = 0; i < size; i++) { // copy data to avoid bad free
645       send_size += sendtmpcounts[i] * dt_size_send;
646       trace_sendcounts->push_back(sendtmpcounts[i] * dt_size_send);
647     }
648
649     TRACE_smpi_comm_in(rank, __FUNCTION__, new simgrid::instr::VarCollTIData(
650                                                "allToAllV", -1, send_size, trace_sendcounts, recv_size,
651                                                trace_recvcounts, encode_datatype(sendtype), encode_datatype(recvtype)));
652
653     retval = simgrid::smpi::Colls::alltoallv(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptype, recvbuf, recvcounts,
654                                     recvdisps, recvtype, comm);
655     TRACE_smpi_comm_out(rank);
656
657     if (sendbuf == MPI_IN_PLACE) {
658       xbt_free(sendtmpbuf);
659       xbt_free(sendtmpcounts);
660       xbt_free(sendtmpdisps);
661     }
662   }
663
664   smpi_bench_begin();
665   return retval;
666 }