Logo AND Algorithmique Numérique Distribuée

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