1 /* Copyright (c) 2007-2019. 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_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"
14 XBT_LOG_EXTERNAL_DEFAULT_CATEGORY(smpi_pmpi);
16 /* PMPI User level calls */
18 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
20 return PMPI_Ibcast(buf, count, datatype, root, comm, MPI_REQUEST_IGNORED);
23 int PMPI_Barrier(MPI_Comm comm)
25 return PMPI_Ibarrier(comm, MPI_REQUEST_IGNORED);
28 int PMPI_Ibarrier(MPI_Comm comm, MPI_Request *request)
32 if (comm == MPI_COMM_NULL) {
33 retval = MPI_ERR_COMM;
34 } else if(request == nullptr){
37 int rank = simgrid::s4u::this_actor::get_pid();
38 TRACE_smpi_comm_in(rank, __func__, new simgrid::instr::NoOpTIData(request==MPI_REQUEST_IGNORED? "barrier" : "ibarrier"));
39 if(request==MPI_REQUEST_IGNORED){
40 simgrid::smpi::Colls::barrier(comm);
41 //Barrier can be used to synchronize RMA calls. Finish all requests from comm before.
42 comm->finish_rma_calls();
44 simgrid::smpi::Colls::ibarrier(comm, request);
45 TRACE_smpi_comm_out(rank);
51 int PMPI_Ibcast(void *buf, int count, MPI_Datatype datatype,
52 int root, MPI_Comm comm, MPI_Request* request)
56 if (comm == MPI_COMM_NULL) {
57 retval = MPI_ERR_COMM;
58 } else if (not datatype->is_valid()) {
60 } else if(request == nullptr){
63 int rank = simgrid::s4u::this_actor::get_pid();
64 TRACE_smpi_comm_in(rank, __func__,
65 new simgrid::instr::CollTIData(request==MPI_REQUEST_IGNORED?"bcast":"ibcast", root, -1.0,
66 datatype->is_replayable() ? count : count * datatype->size(), -1,
67 simgrid::smpi::Datatype::encode(datatype), ""));
68 if (comm->size() > 1){
69 if(request==MPI_REQUEST_IGNORED)
70 simgrid::smpi::Colls::bcast(buf, count, datatype, root, comm);
72 simgrid::smpi::Colls::ibcast(buf, count, datatype, root, comm, request);
76 TRACE_smpi_comm_out(rank);
82 int PMPI_Gather(void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
83 int root, MPI_Comm comm){
84 return PMPI_Igather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
87 int PMPI_Igather(void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
88 int root, MPI_Comm comm, MPI_Request *request)
94 if (comm == MPI_COMM_NULL) {
95 retval = MPI_ERR_COMM;
96 } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
97 ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL))){
98 retval = MPI_ERR_TYPE;
99 } else if ((( sendbuf != MPI_IN_PLACE) && (sendcount <0)) || ((comm->rank() == root) && (recvcount <0))){
100 retval = MPI_ERR_COUNT;
101 } else if (request == nullptr){
102 retval = MPI_ERR_ARG;
105 char* sendtmpbuf = static_cast<char*>(sendbuf);
106 int sendtmpcount = sendcount;
107 MPI_Datatype sendtmptype = sendtype;
108 if( (comm->rank() == root) && (sendbuf == MPI_IN_PLACE )) {
110 sendtmptype=recvtype;
112 int rank = simgrid::s4u::this_actor::get_pid();
116 new simgrid::instr::CollTIData(
117 request==MPI_REQUEST_IGNORED ? "gather":"igather", root, -1.0, sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(),
118 (comm->rank() != root || recvtype->is_replayable()) ? recvcount : recvcount * recvtype->size(),
119 simgrid::smpi::Datatype::encode(sendtmptype), simgrid::smpi::Datatype::encode(recvtype)));
120 if(request == MPI_REQUEST_IGNORED)
121 simgrid::smpi::Colls::gather(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, root, comm);
123 simgrid::smpi::Colls::igather(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, root, comm, request);
125 retval = MPI_SUCCESS;
126 TRACE_smpi_comm_out(rank);
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){
135 return PMPI_Igatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, root, comm, MPI_REQUEST_IGNORED);
138 int PMPI_Igatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, int *recvcounts, int *displs,
139 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request *request)
145 if (comm == MPI_COMM_NULL) {
146 retval = MPI_ERR_COMM;
147 } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
148 ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL))){
149 retval = MPI_ERR_TYPE;
150 } else if (( sendbuf != MPI_IN_PLACE) && (sendcount <0)){
151 retval = MPI_ERR_COUNT;
152 } else if ((comm->rank() == root) && (recvcounts == nullptr || displs == nullptr)) {
153 retval = MPI_ERR_ARG;
154 } else if (request == nullptr){
155 retval = MPI_ERR_ARG;
157 char* sendtmpbuf = static_cast<char*>(sendbuf);
158 int sendtmpcount = sendcount;
159 MPI_Datatype sendtmptype = sendtype;
160 if( (comm->rank() == root) && (sendbuf == MPI_IN_PLACE )) {
162 sendtmptype=recvtype;
165 int rank = simgrid::s4u::this_actor::get_pid();
166 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
168 std::vector<int>* trace_recvcounts = new std::vector<int>;
169 if (comm->rank() == root) {
170 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
171 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
174 TRACE_smpi_comm_in(rank, __func__,
175 new simgrid::instr::VarCollTIData(
176 request==MPI_REQUEST_IGNORED ? "gatherv":"igatherv", root,
177 sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(), nullptr,
178 dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtmptype),
179 simgrid::smpi::Datatype::encode(recvtype)));
180 if(request == MPI_REQUEST_IGNORED)
181 retval = simgrid::smpi::Colls::gatherv(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcounts, displs, recvtype, root, comm);
183 retval = simgrid::smpi::Colls::igatherv(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcounts, displs, recvtype, root, comm, request);
184 TRACE_smpi_comm_out(rank);
191 int PMPI_Allgather(void *sendbuf, int sendcount, MPI_Datatype sendtype,
192 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm){
193 return PMPI_Iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
196 int PMPI_Iallgather(void *sendbuf, int sendcount, MPI_Datatype sendtype,
197 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
199 int retval = MPI_SUCCESS;
203 if (comm == MPI_COMM_NULL) {
204 retval = MPI_ERR_COMM;
205 } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
206 (recvtype == MPI_DATATYPE_NULL)){
207 retval = MPI_ERR_TYPE;
208 } else if ((( sendbuf != MPI_IN_PLACE) && (sendcount <0)) ||
210 retval = MPI_ERR_COUNT;
211 } else if (request == nullptr){
212 retval = MPI_ERR_ARG;
214 if(sendbuf == MPI_IN_PLACE) {
215 sendbuf=static_cast<char*>(recvbuf)+recvtype->get_extent()*recvcount*comm->rank();
219 int rank = simgrid::s4u::this_actor::get_pid();
221 TRACE_smpi_comm_in(rank, __func__,
222 new simgrid::instr::CollTIData(
223 request==MPI_REQUEST_IGNORED ? "allgather" : "iallgather", -1, -1.0, sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
224 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
225 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
226 if(request == MPI_REQUEST_IGNORED)
227 simgrid::smpi::Colls::allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
229 simgrid::smpi::Colls::iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, request);
230 TRACE_smpi_comm_out(rank);
236 int PMPI_Allgatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype,
237 void *recvbuf, int *recvcounts, int *displs, MPI_Datatype recvtype, MPI_Comm comm){
238 return PMPI_Iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, MPI_REQUEST_IGNORED);
241 int PMPI_Iallgatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype,
242 void *recvbuf, int *recvcounts, int *displs, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
248 if (comm == MPI_COMM_NULL) {
249 retval = MPI_ERR_COMM;
250 } else if (((sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) || (recvtype == MPI_DATATYPE_NULL)) {
251 retval = MPI_ERR_TYPE;
252 } else if (( sendbuf != MPI_IN_PLACE) && (sendcount <0)){
253 retval = MPI_ERR_COUNT;
254 } else if (recvcounts == nullptr || displs == nullptr) {
255 retval = MPI_ERR_ARG;
256 } else if (request == nullptr){
257 retval = MPI_ERR_ARG;
260 if(sendbuf == MPI_IN_PLACE) {
261 sendbuf=static_cast<char*>(recvbuf)+recvtype->get_extent()*displs[comm->rank()];
262 sendcount=recvcounts[comm->rank()];
265 int rank = simgrid::s4u::this_actor::get_pid();
266 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
268 std::vector<int>* trace_recvcounts = new std::vector<int>;
269 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
270 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
272 TRACE_smpi_comm_in(rank, __func__,
273 new simgrid::instr::VarCollTIData(
274 request==MPI_REQUEST_IGNORED ? "allgatherv" : "iallgatherv", -1, sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
275 nullptr, dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
276 simgrid::smpi::Datatype::encode(recvtype)));
277 if(request == MPI_REQUEST_IGNORED)
278 simgrid::smpi::Colls::allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
280 simgrid::smpi::Colls::iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, request);
281 retval = MPI_SUCCESS;
282 TRACE_smpi_comm_out(rank);
289 int PMPI_Scatter(void *sendbuf, int sendcount, MPI_Datatype sendtype,
290 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
291 return PMPI_Iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
294 int PMPI_Iscatter(void *sendbuf, int sendcount, MPI_Datatype sendtype,
295 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
301 if (comm == MPI_COMM_NULL) {
302 retval = MPI_ERR_COMM;
303 } else if (((comm->rank() == root) && (not sendtype->is_valid())) ||
304 ((recvbuf != MPI_IN_PLACE) && (not recvtype->is_valid()))) {
305 retval = MPI_ERR_TYPE;
306 } else if ((sendbuf == recvbuf) ||
307 ((comm->rank()==root) && sendcount>0 && (sendbuf == nullptr))){
308 retval = MPI_ERR_BUFFER;
311 if (recvbuf == MPI_IN_PLACE) {
313 recvcount = sendcount;
315 int rank = simgrid::s4u::this_actor::get_pid();
319 new simgrid::instr::CollTIData(
320 request==MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
321 (comm->rank() != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
322 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
323 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
324 if(request == MPI_REQUEST_IGNORED)
325 simgrid::smpi::Colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
327 simgrid::smpi::Colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
328 retval = MPI_SUCCESS;
329 TRACE_smpi_comm_out(rank);
336 int PMPI_Scatterv(void *sendbuf, int *sendcounts, int *displs,
337 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
338 return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
341 int PMPI_Iscatterv(void *sendbuf, int *sendcounts, int *displs,
342 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request *request)
348 if (comm == MPI_COMM_NULL) {
349 retval = MPI_ERR_COMM;
350 } else if (sendcounts == nullptr || displs == nullptr) {
351 retval = MPI_ERR_ARG;
352 } else if (((comm->rank() == root) && (sendtype == MPI_DATATYPE_NULL)) ||
353 ((recvbuf != MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL))) {
354 retval = MPI_ERR_TYPE;
355 } else if (request == nullptr){
356 retval = MPI_ERR_ARG;
358 if (recvbuf == MPI_IN_PLACE) {
360 recvcount = sendcounts[comm->rank()];
362 int rank = simgrid::s4u::this_actor::get_pid();
363 int dt_size_send = sendtype->is_replayable() ? 1 : sendtype->size();
365 std::vector<int>* trace_sendcounts = new std::vector<int>;
366 if (comm->rank() == root) {
367 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
368 trace_sendcounts->push_back(sendcounts[i] * dt_size_send);
371 TRACE_smpi_comm_in(rank, __func__,
372 new simgrid::instr::VarCollTIData(
373 request==MPI_REQUEST_IGNORED ? "scatterv":"iscatterv", root, dt_size_send, trace_sendcounts,
374 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(), nullptr,
375 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
376 if(request == MPI_REQUEST_IGNORED)
377 retval = simgrid::smpi::Colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
379 retval = simgrid::smpi::Colls::iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
381 TRACE_smpi_comm_out(rank);
388 int PMPI_Reduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
394 if (comm == MPI_COMM_NULL) {
395 retval = MPI_ERR_COMM;
396 } else if (not datatype->is_valid() || op == MPI_OP_NULL) {
397 retval = MPI_ERR_ARG;
399 int rank = simgrid::s4u::this_actor::get_pid();
401 TRACE_smpi_comm_in(rank, __func__,
402 new simgrid::instr::CollTIData("reduce", root, 0,
403 datatype->is_replayable() ? count : count * datatype->size(), -1,
404 simgrid::smpi::Datatype::encode(datatype), ""));
406 simgrid::smpi::Colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
408 retval = MPI_SUCCESS;
409 TRACE_smpi_comm_out(rank);
416 int PMPI_Reduce_local(void *inbuf, void *inoutbuf, int count, MPI_Datatype datatype, MPI_Op op){
420 if (not datatype->is_valid() || op == MPI_OP_NULL) {
421 retval = MPI_ERR_ARG;
423 op->apply(inbuf, inoutbuf, &count, datatype);
424 retval = MPI_SUCCESS;
430 int PMPI_Allreduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
432 return PMPI_Iallreduce(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
435 int PMPI_Iallreduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
441 if (comm == MPI_COMM_NULL) {
442 retval = MPI_ERR_COMM;
443 } else if (not datatype->is_valid()) {
444 retval = MPI_ERR_TYPE;
445 } else if (op == MPI_OP_NULL) {
447 } else if (request != MPI_REQUEST_IGNORED) {
448 xbt_die("Iallreduce is not yet implemented. WIP");
449 retval = MPI_ERR_ARG;
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 = simgrid::s4u::this_actor::get_pid();
459 TRACE_smpi_comm_in(rank, __func__,
460 new simgrid::instr::CollTIData(request==MPI_REQUEST_IGNORED ? "allreduce":"iallreduce", -1, 0,
461 datatype->is_replayable() ? count : count * datatype->size(), -1,
462 simgrid::smpi::Datatype::encode(datatype), ""));
464 // if(request == MPI_REQUEST_IGNORED)
465 simgrid::smpi::Colls::allreduce(sendtmpbuf, recvbuf, count, datatype, op, comm);
467 // simgrid::smpi::Colls::iallreduce(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
469 if( sendbuf == MPI_IN_PLACE )
470 xbt_free(sendtmpbuf);
472 retval = MPI_SUCCESS;
473 TRACE_smpi_comm_out(rank);
480 int PMPI_Scan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
486 if (comm == MPI_COMM_NULL) {
487 retval = MPI_ERR_COMM;
488 } else if (not datatype->is_valid()) {
489 retval = MPI_ERR_TYPE;
490 } else if (op == MPI_OP_NULL) {
493 int rank = simgrid::s4u::this_actor::get_pid();
494 void* sendtmpbuf = sendbuf;
495 if (sendbuf == MPI_IN_PLACE) {
496 sendtmpbuf = static_cast<void*>(xbt_malloc(count * datatype->size()));
497 memcpy(sendtmpbuf, recvbuf, count * datatype->size());
499 TRACE_smpi_comm_in(rank, __func__, new simgrid::instr::Pt2PtTIData(
500 "scan", -1, datatype->is_replayable() ? count : count * datatype->size(),
501 simgrid::smpi::Datatype::encode(datatype)));
503 retval = simgrid::smpi::Colls::scan(sendtmpbuf, recvbuf, count, datatype, op, comm);
505 TRACE_smpi_comm_out(rank);
506 if (sendbuf == MPI_IN_PLACE)
507 xbt_free(sendtmpbuf);
514 int PMPI_Exscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm){
519 if (comm == MPI_COMM_NULL) {
520 retval = MPI_ERR_COMM;
521 } else if (not datatype->is_valid()) {
522 retval = MPI_ERR_TYPE;
523 } else if (op == MPI_OP_NULL) {
526 int rank = simgrid::s4u::this_actor::get_pid();
527 void* sendtmpbuf = sendbuf;
528 if (sendbuf == MPI_IN_PLACE) {
529 sendtmpbuf = static_cast<void*>(xbt_malloc(count * datatype->size()));
530 memcpy(sendtmpbuf, recvbuf, count * datatype->size());
533 TRACE_smpi_comm_in(rank, __func__, new simgrid::instr::Pt2PtTIData(
534 "exscan", -1, datatype->is_replayable() ? count : count * datatype->size(),
535 simgrid::smpi::Datatype::encode(datatype)));
537 retval = simgrid::smpi::Colls::exscan(sendtmpbuf, recvbuf, count, datatype, op, comm);
539 TRACE_smpi_comm_out(rank);
540 if (sendbuf == MPI_IN_PLACE)
541 xbt_free(sendtmpbuf);
548 int PMPI_Reduce_scatter(void *sendbuf, void *recvbuf, int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
553 if (comm == MPI_COMM_NULL) {
554 retval = MPI_ERR_COMM;
555 } else if (not datatype->is_valid()) {
556 retval = MPI_ERR_TYPE;
557 } else if (op == MPI_OP_NULL) {
559 } else if (recvcounts == nullptr) {
560 retval = MPI_ERR_ARG;
562 int rank = simgrid::s4u::this_actor::get_pid();
563 std::vector<int>* trace_recvcounts = new std::vector<int>;
564 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
567 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
568 trace_recvcounts->push_back(recvcounts[i] * dt_send_size);
569 totalcount += recvcounts[i];
572 void* sendtmpbuf = sendbuf;
573 if (sendbuf == MPI_IN_PLACE) {
574 sendtmpbuf = static_cast<void*>(xbt_malloc(totalcount * datatype->size()));
575 memcpy(sendtmpbuf, recvbuf, totalcount * datatype->size());
578 TRACE_smpi_comm_in(rank, __func__, new simgrid::instr::VarCollTIData(
579 "reducescatter", -1, dt_send_size, nullptr, -1, trace_recvcounts,
580 simgrid::smpi::Datatype::encode(datatype), ""));
582 simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
583 retval = MPI_SUCCESS;
584 TRACE_smpi_comm_out(rank);
586 if (sendbuf == MPI_IN_PLACE)
587 xbt_free(sendtmpbuf);
594 int PMPI_Reduce_scatter_block(void *sendbuf, void *recvbuf, int recvcount,
595 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
600 if (comm == MPI_COMM_NULL) {
601 retval = MPI_ERR_COMM;
602 } else if (not datatype->is_valid()) {
603 retval = MPI_ERR_TYPE;
604 } else if (op == MPI_OP_NULL) {
606 } else if (recvcount < 0) {
607 retval = MPI_ERR_ARG;
609 int count = comm->size();
611 int rank = simgrid::s4u::this_actor::get_pid();
612 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
613 std::vector<int>* trace_recvcounts = new std::vector<int>(recvcount * dt_send_size); // copy data to avoid bad free
615 void* sendtmpbuf = sendbuf;
616 if (sendbuf == MPI_IN_PLACE) {
617 sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * count * datatype->size()));
618 memcpy(sendtmpbuf, recvbuf, recvcount * count * datatype->size());
621 TRACE_smpi_comm_in(rank, __func__,
622 new simgrid::instr::VarCollTIData("reducescatter", -1, 0, nullptr, -1, trace_recvcounts,
623 simgrid::smpi::Datatype::encode(datatype), ""));
625 int* recvcounts = new int[count];
626 for (int i = 0; i < count; i++)
627 recvcounts[i] = recvcount;
628 simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
630 retval = MPI_SUCCESS;
632 TRACE_smpi_comm_out(rank);
634 if (sendbuf == MPI_IN_PLACE)
635 xbt_free(sendtmpbuf);
641 int PMPI_Alltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
642 MPI_Datatype recvtype, MPI_Comm comm){
643 return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
646 int PMPI_Ialltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
647 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request *request)
652 if (comm == MPI_COMM_NULL) {
653 retval = MPI_ERR_COMM;
654 } else if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL) {
655 retval = MPI_ERR_TYPE;
656 } else if (request == nullptr){
657 retval = MPI_ERR_ARG;
659 int rank = simgrid::s4u::this_actor::get_pid();
660 void* sendtmpbuf = static_cast<char*>(sendbuf);
661 int sendtmpcount = sendcount;
662 MPI_Datatype sendtmptype = sendtype;
663 if (sendbuf == MPI_IN_PLACE) {
664 sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * comm->size() * recvtype->size()));
665 memcpy(sendtmpbuf, recvbuf, recvcount * comm->size() * recvtype->size());
666 sendtmpcount = recvcount;
667 sendtmptype = recvtype;
670 TRACE_smpi_comm_in(rank, __func__,
671 new simgrid::instr::CollTIData(
672 request==MPI_REQUEST_IGNORED ? "alltoall" : "ialltoall", -1, -1.0,
673 sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(),
674 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
675 simgrid::smpi::Datatype::encode(sendtmptype), simgrid::smpi::Datatype::encode(recvtype)));
676 if(request == MPI_REQUEST_IGNORED)
677 retval = simgrid::smpi::Colls::alltoall(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, comm);
679 retval = simgrid::smpi::Colls::ialltoall(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, comm, request);
681 TRACE_smpi_comm_out(rank);
683 if (sendbuf == MPI_IN_PLACE)
684 xbt_free(sendtmpbuf);
691 int PMPI_Alltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf,
692 int* recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm)
694 return PMPI_Ialltoallv(sendbuf, sendcounts, senddisps, sendtype, recvbuf, recvcounts, recvdisps, recvtype, comm, MPI_REQUEST_IGNORED);
697 int PMPI_Ialltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf,
698 int* recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request *request)
704 if (comm == MPI_COMM_NULL) {
705 retval = MPI_ERR_COMM;
706 } else if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL) {
707 retval = MPI_ERR_TYPE;
708 } else if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
709 recvdisps == nullptr) {
710 retval = MPI_ERR_ARG;
711 } else if (request == nullptr){
712 retval = MPI_ERR_ARG;
714 int rank = simgrid::s4u::this_actor::get_pid();
715 int size = comm->size();
718 std::vector<int>* trace_sendcounts = new std::vector<int>;
719 std::vector<int>* trace_recvcounts = new std::vector<int>;
720 int dt_size_recv = recvtype->size();
722 void* sendtmpbuf = static_cast<char*>(sendbuf);
723 int* sendtmpcounts = sendcounts;
724 int* sendtmpdisps = senddisps;
725 MPI_Datatype sendtmptype = sendtype;
727 for (int i = 0; i < size; i++) { // copy data to avoid bad free
728 recv_size += recvcounts[i] * dt_size_recv;
729 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
730 if (((recvdisps[i] + recvcounts[i]) * dt_size_recv) > maxsize)
731 maxsize = (recvdisps[i] + recvcounts[i]) * dt_size_recv;
734 if (sendbuf == MPI_IN_PLACE) {
735 sendtmpbuf = static_cast<void*>(xbt_malloc(maxsize));
736 memcpy(sendtmpbuf, recvbuf, maxsize);
737 sendtmpcounts = static_cast<int*>(xbt_malloc(size * sizeof(int)));
738 memcpy(sendtmpcounts, recvcounts, size * sizeof(int));
739 sendtmpdisps = static_cast<int*>(xbt_malloc(size * sizeof(int)));
740 memcpy(sendtmpdisps, recvdisps, size * sizeof(int));
741 sendtmptype = recvtype;
744 int dt_size_send = sendtmptype->size();
746 for (int i = 0; i < size; i++) { // copy data to avoid bad free
747 send_size += sendtmpcounts[i] * dt_size_send;
748 trace_sendcounts->push_back(sendtmpcounts[i] * dt_size_send);
751 TRACE_smpi_comm_in(rank, __func__,
752 new simgrid::instr::VarCollTIData(request==MPI_REQUEST_IGNORED ? "alltoallv":"ialltoallv", -1, send_size, trace_sendcounts, recv_size,
753 trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
754 simgrid::smpi::Datatype::encode(recvtype)));
756 if(request == MPI_REQUEST_IGNORED)
757 retval = simgrid::smpi::Colls::alltoallv(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptype, recvbuf, recvcounts,
758 recvdisps, recvtype, comm);
760 retval = simgrid::smpi::Colls::ialltoallv(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptype, recvbuf, recvcounts,
761 recvdisps, recvtype, comm, request);
762 TRACE_smpi_comm_out(rank);
764 if (sendbuf == MPI_IN_PLACE) {
765 xbt_free(sendtmpbuf);
766 xbt_free(sendtmpcounts);
767 xbt_free(sendtmpdisps);