1 /* Copyright (c) 2007-2017. The SimGrid Team. All rights reserved. */
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. */
7 #include "smpi_comm.hpp"
8 #include "smpi_coll.hpp"
9 #include "smpi_datatype_derived.hpp"
10 #include "smpi_op.hpp"
11 #include "smpi_process.hpp"
13 XBT_LOG_EXTERNAL_DEFAULT_CATEGORY(smpi_pmpi);
16 /* PMPI User level calls */
17 extern "C" { // Obviously, the C MPI interface should use the C linkage
19 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
25 if (comm == MPI_COMM_NULL) {
26 retval = MPI_ERR_COMM;
27 } else if (not datatype->is_valid()) {
30 int rank = comm != MPI_COMM_NULL ? smpi_process()->index() : -1;
31 int root_traced = comm->group()->index(root);
33 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
34 extra->type = TRACING_BCAST;
35 extra->root = root_traced;
37 extra->datatype1 = encode_datatype(datatype, &known);
40 dt_size_send = datatype->size();
41 extra->send_size = count * dt_size_send;
42 TRACE_smpi_collective_in(rank, __FUNCTION__, extra);
44 simgrid::smpi::Colls::bcast(buf, count, datatype, root, comm);
47 TRACE_smpi_collective_out(rank, __FUNCTION__);
53 int PMPI_Barrier(MPI_Comm comm)
59 if (comm == MPI_COMM_NULL) {
60 retval = MPI_ERR_COMM;
62 int rank = comm != MPI_COMM_NULL ? smpi_process()->index() : -1;
63 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
64 extra->type = TRACING_BARRIER;
65 TRACE_smpi_collective_in(rank, __FUNCTION__, extra);
67 simgrid::smpi::Colls::barrier(comm);
69 //Barrier can be used to synchronize RMA calls. Finish all requests from comm before.
70 comm->finish_rma_calls();
74 TRACE_smpi_collective_out(rank, __FUNCTION__);
81 int PMPI_Gather(void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
82 int root, MPI_Comm comm)
88 if (comm == MPI_COMM_NULL) {
89 retval = MPI_ERR_COMM;
90 } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
91 ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL))){
92 retval = MPI_ERR_TYPE;
93 } else if ((( sendbuf != MPI_IN_PLACE) && (sendcount <0)) || ((comm->rank() == root) && (recvcount <0))){
94 retval = MPI_ERR_COUNT;
97 char* sendtmpbuf = static_cast<char*>(sendbuf);
98 int sendtmpcount = sendcount;
99 MPI_Datatype sendtmptype = sendtype;
100 if( (comm->rank() == root) && (sendbuf == MPI_IN_PLACE )) {
102 sendtmptype=recvtype;
104 int rank = comm != MPI_COMM_NULL ? smpi_process()->index() : -1;
105 int root_traced = comm->group()->index(root);
106 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
107 extra->type = TRACING_GATHER;
108 extra->root = root_traced;
110 extra->datatype1 = encode_datatype(sendtmptype, &known);
111 int dt_size_send = 1;
113 dt_size_send = sendtmptype->size();
114 extra->send_size = sendtmpcount * dt_size_send;
115 extra->datatype2 = encode_datatype(recvtype, &known);
116 int dt_size_recv = 1;
117 if ((comm->rank() == root) && known == 0)
118 dt_size_recv = recvtype->size();
119 extra->recv_size = recvcount * dt_size_recv;
121 TRACE_smpi_collective_in(rank, __FUNCTION__, extra);
123 simgrid::smpi::Colls::gather(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, root, comm);
125 retval = MPI_SUCCESS;
126 TRACE_smpi_collective_out(rank, __FUNCTION__);
133 int PMPI_Gatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, int *recvcounts, int *displs,
134 MPI_Datatype recvtype, int root, MPI_Comm comm)
140 if (comm == MPI_COMM_NULL) {
141 retval = MPI_ERR_COMM;
142 } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
143 ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL))){
144 retval = MPI_ERR_TYPE;
145 } else if (( sendbuf != MPI_IN_PLACE) && (sendcount <0)){
146 retval = MPI_ERR_COUNT;
147 } else if (recvcounts == nullptr || displs == nullptr) {
148 retval = MPI_ERR_ARG;
150 char* sendtmpbuf = static_cast<char*>(sendbuf);
151 int sendtmpcount = sendcount;
152 MPI_Datatype sendtmptype = sendtype;
153 if( (comm->rank() == root) && (sendbuf == MPI_IN_PLACE )) {
155 sendtmptype=recvtype;
158 int rank = comm != MPI_COMM_NULL ? smpi_process()->index() : -1;
159 int root_traced = comm->group()->index(root);
160 int size = comm->size();
161 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
162 extra->type = TRACING_GATHERV;
163 extra->num_processes = size;
164 extra->root = root_traced;
166 extra->datatype1 = encode_datatype(sendtmptype, &known);
167 int dt_size_send = 1;
169 dt_size_send = sendtype->size();
170 extra->send_size = sendtmpcount * dt_size_send;
171 extra->datatype2 = encode_datatype(recvtype, &known);
172 int dt_size_recv = 1;
174 dt_size_recv = recvtype->size();
175 if (comm->rank() == root) {
176 extra->recvcounts = xbt_new(int, size);
177 for (int i = 0; i < size; i++) // copy data to avoid bad free
178 extra->recvcounts[i] = recvcounts[i] * dt_size_recv;
180 TRACE_smpi_collective_in(rank, __FUNCTION__, extra);
182 retval = simgrid::smpi::Colls::gatherv(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcounts, displs, recvtype, root, comm);
183 TRACE_smpi_collective_out(rank, __FUNCTION__);
190 int PMPI_Allgather(void *sendbuf, int sendcount, MPI_Datatype sendtype,
191 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm)
197 if (comm == MPI_COMM_NULL) {
198 retval = MPI_ERR_COMM;
199 } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
200 (recvtype == MPI_DATATYPE_NULL)){
201 retval = MPI_ERR_TYPE;
202 } else if ((( sendbuf != MPI_IN_PLACE) && (sendcount <0)) ||
204 retval = MPI_ERR_COUNT;
206 if(sendbuf == MPI_IN_PLACE) {
207 sendbuf=static_cast<char*>(recvbuf)+recvtype->get_extent()*recvcount*comm->rank();
211 int rank = comm != MPI_COMM_NULL ? smpi_process()->index() : -1;
212 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
213 extra->type = TRACING_ALLGATHER;
215 extra->datatype1 = encode_datatype(sendtype, &known);
216 int dt_size_send = 1;
218 dt_size_send = sendtype->size();
219 extra->send_size = sendcount * dt_size_send;
220 extra->datatype2 = encode_datatype(recvtype, &known);
221 int dt_size_recv = 1;
223 dt_size_recv = recvtype->size();
224 extra->recv_size = recvcount * dt_size_recv;
226 TRACE_smpi_collective_in(rank, __FUNCTION__, extra);
228 simgrid::smpi::Colls::allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
229 retval = MPI_SUCCESS;
230 TRACE_smpi_collective_out(rank, __FUNCTION__);
236 int PMPI_Allgatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype,
237 void *recvbuf, int *recvcounts, int *displs, MPI_Datatype recvtype, MPI_Comm comm)
243 if (comm == MPI_COMM_NULL) {
244 retval = MPI_ERR_COMM;
245 } else if (((sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) || (recvtype == MPI_DATATYPE_NULL)) {
246 retval = MPI_ERR_TYPE;
247 } else if (( sendbuf != MPI_IN_PLACE) && (sendcount <0)){
248 retval = MPI_ERR_COUNT;
249 } else if (recvcounts == nullptr || displs == nullptr) {
250 retval = MPI_ERR_ARG;
253 if(sendbuf == MPI_IN_PLACE) {
254 sendbuf=static_cast<char*>(recvbuf)+recvtype->get_extent()*displs[comm->rank()];
255 sendcount=recvcounts[comm->rank()];
258 int rank = comm != MPI_COMM_NULL ? smpi_process()->index() : -1;
260 int size = comm->size();
261 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
262 extra->type = TRACING_ALLGATHERV;
263 extra->num_processes = size;
265 extra->datatype1 = encode_datatype(sendtype, &known);
266 int dt_size_send = 1;
268 dt_size_send = sendtype->size();
269 extra->send_size = sendcount * dt_size_send;
270 extra->datatype2 = encode_datatype(recvtype, &known);
271 int dt_size_recv = 1;
273 dt_size_recv = recvtype->size();
274 extra->recvcounts = xbt_new(int, size);
275 for (i = 0; i < size; i++) // copy data to avoid bad free
276 extra->recvcounts[i] = recvcounts[i] * dt_size_recv;
278 TRACE_smpi_collective_in(rank, __FUNCTION__, extra);
280 simgrid::smpi::Colls::allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
281 retval = MPI_SUCCESS;
282 TRACE_smpi_collective_out(rank, __FUNCTION__);
289 int PMPI_Scatter(void *sendbuf, int sendcount, MPI_Datatype sendtype,
290 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm)
296 if (comm == MPI_COMM_NULL) {
297 retval = MPI_ERR_COMM;
298 } else if (((comm->rank() == root) && (not sendtype->is_valid())) ||
299 ((recvbuf != MPI_IN_PLACE) && (not recvtype->is_valid()))) {
300 retval = MPI_ERR_TYPE;
301 } else if ((sendbuf == recvbuf) ||
302 ((comm->rank()==root) && sendcount>0 && (sendbuf == nullptr))){
303 retval = MPI_ERR_BUFFER;
306 if (recvbuf == MPI_IN_PLACE) {
308 recvcount = sendcount;
310 int rank = comm != MPI_COMM_NULL ? smpi_process()->index() : -1;
311 int root_traced = comm->group()->index(root);
312 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
313 extra->type = TRACING_SCATTER;
314 extra->root = root_traced;
316 extra->datatype1 = encode_datatype(sendtype, &known);
317 int dt_size_send = 1;
318 if ((comm->rank() == root) && known == 0)
319 dt_size_send = sendtype->size();
320 extra->send_size = sendcount * dt_size_send;
321 extra->datatype2 = encode_datatype(recvtype, &known);
322 int dt_size_recv = 1;
324 dt_size_recv = recvtype->size();
325 extra->recv_size = recvcount * dt_size_recv;
326 TRACE_smpi_collective_in(rank, __FUNCTION__, extra);
328 simgrid::smpi::Colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
329 retval = MPI_SUCCESS;
330 TRACE_smpi_collective_out(rank, __FUNCTION__);
337 int PMPI_Scatterv(void *sendbuf, int *sendcounts, int *displs,
338 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm)
344 if (comm == MPI_COMM_NULL) {
345 retval = MPI_ERR_COMM;
346 } else if (sendcounts == nullptr || displs == nullptr) {
347 retval = MPI_ERR_ARG;
348 } else if (((comm->rank() == root) && (sendtype == MPI_DATATYPE_NULL)) ||
349 ((recvbuf != MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL))) {
350 retval = MPI_ERR_TYPE;
352 if (recvbuf == MPI_IN_PLACE) {
354 recvcount = sendcounts[comm->rank()];
356 int rank = comm != MPI_COMM_NULL ? smpi_process()->index() : -1;
357 int root_traced = comm->group()->index(root);
358 int size = comm->size();
359 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
360 extra->type = TRACING_SCATTERV;
361 extra->num_processes = size;
362 extra->root = root_traced;
364 extra->datatype1 = encode_datatype(sendtype, &known);
365 int dt_size_send = 1;
367 dt_size_send = sendtype->size();
368 if (comm->rank() == root) {
369 extra->sendcounts = xbt_new(int, size);
370 for (int i = 0; i < size; i++) // copy data to avoid bad free
371 extra->sendcounts[i] = sendcounts[i] * dt_size_send;
373 extra->datatype2 = encode_datatype(recvtype, &known);
374 int dt_size_recv = 1;
376 dt_size_recv = recvtype->size();
377 extra->recv_size = recvcount * dt_size_recv;
378 TRACE_smpi_collective_in(rank, __FUNCTION__, extra);
380 retval = simgrid::smpi::Colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
382 TRACE_smpi_collective_out(rank, __FUNCTION__);
389 int PMPI_Reduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
395 if (comm == MPI_COMM_NULL) {
396 retval = MPI_ERR_COMM;
397 } else if (not datatype->is_valid() || op == MPI_OP_NULL) {
398 retval = MPI_ERR_ARG;
400 int rank = comm != MPI_COMM_NULL ? smpi_process()->index() : -1;
401 int root_traced = comm->group()->index(root);
402 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
403 extra->type = TRACING_REDUCE;
405 extra->datatype1 = encode_datatype(datatype, &known);
406 int dt_size_send = 1;
408 dt_size_send = datatype->size();
409 extra->send_size = count * dt_size_send;
410 extra->root = root_traced;
412 TRACE_smpi_collective_in(rank, __FUNCTION__, extra);
414 simgrid::smpi::Colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
416 retval = MPI_SUCCESS;
417 TRACE_smpi_collective_out(rank, __FUNCTION__);
424 int PMPI_Reduce_local(void *inbuf, void *inoutbuf, int count, MPI_Datatype datatype, MPI_Op op){
428 if (not datatype->is_valid() || op == MPI_OP_NULL) {
429 retval = MPI_ERR_ARG;
431 op->apply(inbuf, inoutbuf, &count, datatype);
432 retval = MPI_SUCCESS;
438 int PMPI_Allreduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
444 if (comm == MPI_COMM_NULL) {
445 retval = MPI_ERR_COMM;
446 } else if (not datatype->is_valid()) {
447 retval = MPI_ERR_TYPE;
448 } else if (op == MPI_OP_NULL) {
452 char* sendtmpbuf = static_cast<char*>(sendbuf);
453 if( sendbuf == MPI_IN_PLACE ) {
454 sendtmpbuf = static_cast<char*>(xbt_malloc(count*datatype->get_extent()));
455 simgrid::smpi::Datatype::copy(recvbuf, count, datatype,sendtmpbuf, count, datatype);
457 int rank = comm != MPI_COMM_NULL ? smpi_process()->index() : -1;
458 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
459 extra->type = TRACING_ALLREDUCE;
461 extra->datatype1 = encode_datatype(datatype, &known);
462 int dt_size_send = 1;
464 dt_size_send = datatype->size();
465 extra->send_size = count * dt_size_send;
467 TRACE_smpi_collective_in(rank, __FUNCTION__, extra);
469 simgrid::smpi::Colls::allreduce(sendtmpbuf, recvbuf, count, datatype, op, comm);
471 if( sendbuf == MPI_IN_PLACE )
472 xbt_free(sendtmpbuf);
474 retval = MPI_SUCCESS;
475 TRACE_smpi_collective_out(rank, __FUNCTION__);
482 int PMPI_Scan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
488 if (comm == MPI_COMM_NULL) {
489 retval = MPI_ERR_COMM;
490 } else if (not datatype->is_valid()) {
491 retval = MPI_ERR_TYPE;
492 } else if (op == MPI_OP_NULL) {
495 int rank = comm != MPI_COMM_NULL ? smpi_process()->index() : -1;
496 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
497 extra->type = TRACING_SCAN;
499 extra->datatype1 = encode_datatype(datatype, &known);
500 int dt_size_send = 1;
502 dt_size_send = datatype->size();
503 extra->send_size = count * dt_size_send;
505 TRACE_smpi_collective_in(rank, __FUNCTION__, extra);
507 retval = simgrid::smpi::Colls::scan(sendbuf, recvbuf, count, datatype, op, comm);
509 TRACE_smpi_collective_out(rank, __FUNCTION__);
516 int PMPI_Exscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm){
521 if (comm == MPI_COMM_NULL) {
522 retval = MPI_ERR_COMM;
523 } else if (not datatype->is_valid()) {
524 retval = MPI_ERR_TYPE;
525 } else if (op == MPI_OP_NULL) {
528 int rank = comm != MPI_COMM_NULL ? smpi_process()->index() : -1;
529 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
530 extra->type = TRACING_EXSCAN;
532 extra->datatype1 = encode_datatype(datatype, &known);
533 int dt_size_send = 1;
535 dt_size_send = datatype->size();
536 extra->send_size = count * dt_size_send;
537 void* sendtmpbuf = sendbuf;
538 if (sendbuf == MPI_IN_PLACE) {
539 sendtmpbuf = static_cast<void*>(xbt_malloc(count * datatype->size()));
540 memcpy(sendtmpbuf, recvbuf, count * datatype->size());
542 TRACE_smpi_collective_in(rank, __FUNCTION__, extra);
544 retval = simgrid::smpi::Colls::exscan(sendtmpbuf, recvbuf, count, datatype, op, comm);
546 TRACE_smpi_collective_out(rank, __FUNCTION__);
547 if (sendbuf == MPI_IN_PLACE)
548 xbt_free(sendtmpbuf);
555 int PMPI_Reduce_scatter(void *sendbuf, void *recvbuf, int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
560 if (comm == MPI_COMM_NULL) {
561 retval = MPI_ERR_COMM;
562 } else if (not datatype->is_valid()) {
563 retval = MPI_ERR_TYPE;
564 } else if (op == MPI_OP_NULL) {
566 } else if (recvcounts == nullptr) {
567 retval = MPI_ERR_ARG;
569 int rank = comm != MPI_COMM_NULL ? smpi_process()->index() : -1;
571 int size = comm->size();
572 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
573 extra->type = TRACING_REDUCE_SCATTER;
574 extra->num_processes = size;
576 extra->datatype1 = encode_datatype(datatype, &known);
577 int dt_size_send = 1;
579 dt_size_send = datatype->size();
580 extra->send_size = 0;
581 extra->recvcounts = xbt_new(int, size);
583 for (i = 0; i < size; i++) { // copy data to avoid bad free
584 extra->recvcounts[i] = recvcounts[i] * dt_size_send;
585 totalcount += recvcounts[i];
587 void* sendtmpbuf = sendbuf;
588 if (sendbuf == MPI_IN_PLACE) {
589 sendtmpbuf = static_cast<void*>(xbt_malloc(totalcount * datatype->size()));
590 memcpy(sendtmpbuf, recvbuf, totalcount * datatype->size());
593 TRACE_smpi_collective_in(rank, __FUNCTION__, extra);
595 simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
596 retval = MPI_SUCCESS;
597 TRACE_smpi_collective_out(rank, __FUNCTION__);
599 if (sendbuf == MPI_IN_PLACE)
600 xbt_free(sendtmpbuf);
607 int PMPI_Reduce_scatter_block(void *sendbuf, void *recvbuf, int recvcount,
608 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
613 if (comm == MPI_COMM_NULL) {
614 retval = MPI_ERR_COMM;
615 } else if (not datatype->is_valid()) {
616 retval = MPI_ERR_TYPE;
617 } else if (op == MPI_OP_NULL) {
619 } else if (recvcount < 0) {
620 retval = MPI_ERR_ARG;
622 int count = comm->size();
624 int rank = comm != MPI_COMM_NULL ? smpi_process()->index() : -1;
625 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
626 extra->type = TRACING_REDUCE_SCATTER;
627 extra->num_processes = count;
629 extra->datatype1 = encode_datatype(datatype, &known);
630 int dt_size_send = 1;
632 dt_size_send = datatype->size();
633 extra->send_size = 0;
634 extra->recvcounts = xbt_new(int, count);
635 for (int i = 0; i < count; i++) // copy data to avoid bad free
636 extra->recvcounts[i] = recvcount * dt_size_send;
637 void* sendtmpbuf = sendbuf;
638 if (sendbuf == MPI_IN_PLACE) {
639 sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * count * datatype->size()));
640 memcpy(sendtmpbuf, recvbuf, recvcount * count * datatype->size());
643 TRACE_smpi_collective_in(rank, __FUNCTION__, extra);
645 int* recvcounts = static_cast<int*>(xbt_malloc(count * sizeof(int)));
646 for (int i = 0; i < count; i++)
647 recvcounts[i] = recvcount;
648 simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
649 xbt_free(recvcounts);
650 retval = MPI_SUCCESS;
652 TRACE_smpi_collective_out(rank, __FUNCTION__);
654 if (sendbuf == MPI_IN_PLACE)
655 xbt_free(sendtmpbuf);
662 int PMPI_Alltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
663 MPI_Datatype recvtype, MPI_Comm comm)
668 if (comm == MPI_COMM_NULL) {
669 retval = MPI_ERR_COMM;
670 } else if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL) {
671 retval = MPI_ERR_TYPE;
673 int rank = comm != MPI_COMM_NULL ? smpi_process()->index() : -1;
674 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
675 extra->type = TRACING_ALLTOALL;
677 void* sendtmpbuf = static_cast<char*>(sendbuf);
678 int sendtmpcount = sendcount;
679 MPI_Datatype sendtmptype = sendtype;
680 if (sendbuf == MPI_IN_PLACE) {
681 sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * comm->size() * recvtype->size()));
682 memcpy(sendtmpbuf, recvbuf, recvcount * comm->size() * recvtype->size());
683 sendtmpcount = recvcount;
684 sendtmptype = recvtype;
688 extra->datatype1 = encode_datatype(sendtmptype, &known);
690 extra->send_size = sendtmpcount * sendtmptype->size();
692 extra->send_size = sendtmpcount;
693 extra->datatype2 = encode_datatype(recvtype, &known);
695 extra->recv_size = recvcount * recvtype->size();
697 extra->recv_size = recvcount;
699 TRACE_smpi_collective_in(rank, __FUNCTION__, extra);
701 retval = simgrid::smpi::Colls::alltoall(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, comm);
703 TRACE_smpi_collective_out(rank, __FUNCTION__);
705 if (sendbuf == MPI_IN_PLACE)
706 xbt_free(sendtmpbuf);
713 int PMPI_Alltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf,
714 int* recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm)
720 if (comm == MPI_COMM_NULL) {
721 retval = MPI_ERR_COMM;
722 } else if (sendtype == MPI_DATATYPE_NULL || recvtype == MPI_DATATYPE_NULL) {
723 retval = MPI_ERR_TYPE;
724 } else if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
725 recvdisps == nullptr) {
726 retval = MPI_ERR_ARG;
728 int rank = comm != MPI_COMM_NULL ? smpi_process()->index() : -1;
730 int size = comm->size();
731 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
732 extra->type = TRACING_ALLTOALLV;
733 extra->send_size = 0;
734 extra->recv_size = 0;
735 extra->recvcounts = xbt_new(int, size);
736 extra->sendcounts = xbt_new(int, size);
738 extra->datatype2 = encode_datatype(recvtype, &known);
739 int dt_size_recv = recvtype->size();
741 void* sendtmpbuf = static_cast<char*>(sendbuf);
742 int* sendtmpcounts = sendcounts;
743 int* sendtmpdisps = senddisps;
744 MPI_Datatype sendtmptype = sendtype;
746 for (i = 0; i < size; i++) { // copy data to avoid bad free
747 extra->recv_size += recvcounts[i] * dt_size_recv;
748 extra->recvcounts[i] = recvcounts[i] * dt_size_recv;
749 if (((recvdisps[i] + recvcounts[i]) * dt_size_recv) > maxsize)
750 maxsize = (recvdisps[i] + recvcounts[i]) * dt_size_recv;
753 if (sendbuf == MPI_IN_PLACE) {
754 sendtmpbuf = static_cast<void*>(xbt_malloc(maxsize));
755 memcpy(sendtmpbuf, recvbuf, maxsize);
756 sendtmpcounts = static_cast<int*>(xbt_malloc(size * sizeof(int)));
757 memcpy(sendtmpcounts, recvcounts, size * sizeof(int));
758 sendtmpdisps = static_cast<int*>(xbt_malloc(size * sizeof(int)));
759 memcpy(sendtmpdisps, recvdisps, size * sizeof(int));
760 sendtmptype = recvtype;
763 extra->datatype1 = encode_datatype(sendtmptype, &known);
764 int dt_size_send = sendtmptype->size();
766 for (i = 0; i < size; i++) { // copy data to avoid bad free
767 extra->send_size += sendtmpcounts[i] * dt_size_send;
768 extra->sendcounts[i] = sendtmpcounts[i] * dt_size_send;
770 extra->num_processes = size;
771 TRACE_smpi_collective_in(rank, __FUNCTION__, extra);
772 retval = simgrid::smpi::Colls::alltoallv(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptype, recvbuf, recvcounts,
773 recvdisps, recvtype, comm);
774 TRACE_smpi_collective_out(rank, __FUNCTION__);
776 if (sendbuf == MPI_IN_PLACE) {
777 xbt_free(sendtmpbuf);
778 xbt_free(sendtmpcounts);
779 xbt_free(sendtmpdisps);