1 /* Copyright (c) 2007-2020. 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"
16 XBT_LOG_EXTERNAL_DEFAULT_CATEGORY(smpi_pmpi);
18 static const void* smpi_get_in_place_buf(const void* inplacebuf, const void* otherbuf,
19 std::vector<unsigned char>& tmp_sendbuf, int count, MPI_Datatype datatype)
21 if (inplacebuf == MPI_IN_PLACE) {
22 tmp_sendbuf.resize(count * datatype->get_extent());
23 simgrid::smpi::Datatype::copy(otherbuf, count, datatype, tmp_sendbuf.data(), count, datatype);
24 return tmp_sendbuf.data();
29 /* PMPI User level calls */
31 int PMPI_Barrier(MPI_Comm comm)
33 return PMPI_Ibarrier(comm, MPI_REQUEST_IGNORED);
36 int PMPI_Ibarrier(MPI_Comm comm, MPI_Request *request)
42 int rank = simgrid::s4u::this_actor::get_pid();
43 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Barrier" : "PMPI_Ibarrier",
44 new simgrid::instr::NoOpTIData(request == MPI_REQUEST_IGNORED ? "barrier" : "ibarrier"));
45 if (request == MPI_REQUEST_IGNORED) {
46 simgrid::smpi::colls::barrier(comm);
47 // Barrier can be used to synchronize RMA calls. Finish all requests from comm before.
48 comm->finish_rma_calls();
50 simgrid::smpi::colls::ibarrier(comm, request);
52 TRACE_smpi_comm_out(rank);
57 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
59 return PMPI_Ibcast(buf, count, datatype, root, comm, MPI_REQUEST_IGNORED);
62 int PMPI_Ibcast(void *buf, int count, MPI_Datatype datatype,
63 int root, MPI_Comm comm, MPI_Request* request)
66 CHECK_BUFFER(1, buf, count)
68 CHECK_TYPE(3, datatype)
73 int rank = simgrid::s4u::this_actor::get_pid();
74 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Bcast" : "PMPI_Ibcast",
75 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "bcast" : "ibcast", root, -1.0,
76 datatype->is_replayable() ? count : count * datatype->size(), -1,
77 simgrid::smpi::Datatype::encode(datatype), ""));
78 if (comm->size() > 1) {
79 if (request == MPI_REQUEST_IGNORED)
80 simgrid::smpi::colls::bcast(buf, count, datatype, root, comm);
82 simgrid::smpi::colls::ibcast(buf, count, datatype, root, comm, request);
84 if (request != MPI_REQUEST_IGNORED)
85 *request = MPI_REQUEST_NULL;
88 TRACE_smpi_comm_out(rank);
93 int PMPI_Gather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
94 int root, MPI_Comm comm){
95 return PMPI_Igather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
98 int PMPI_Igather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
99 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
102 if(sendbuf != MPI_IN_PLACE){
103 CHECK_BUFFER(1,sendbuf, sendcount)
104 CHECK_COUNT(2, sendcount)
105 CHECK_TYPE(3, sendtype)
107 if(comm->rank() == root){
108 CHECK_TYPE(6, recvtype)
109 CHECK_COUNT(5, recvcount)
110 CHECK_BUFFER(4, recvbuf, recvcount)
116 const void* real_sendbuf = sendbuf;
117 int real_sendcount = sendcount;
118 MPI_Datatype real_sendtype = sendtype;
119 if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
121 real_sendtype = recvtype;
123 int rank = simgrid::s4u::this_actor::get_pid();
125 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gather" : "PMPI_Igather",
126 new simgrid::instr::CollTIData(
127 request == MPI_REQUEST_IGNORED ? "gather" : "igather", root, -1.0,
128 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
129 (comm->rank() != root || recvtype->is_replayable()) ? recvcount : recvcount * recvtype->size(),
130 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
131 if (request == MPI_REQUEST_IGNORED)
132 simgrid::smpi::colls::gather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm);
134 simgrid::smpi::colls::igather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm,
137 TRACE_smpi_comm_out(rank);
142 int PMPI_Gatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, const int *recvcounts, const int *displs,
143 MPI_Datatype recvtype, int root, MPI_Comm comm){
144 return PMPI_Igatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, root, comm, MPI_REQUEST_IGNORED);
147 int PMPI_Igatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
148 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
151 CHECK_BUFFER(1, sendbuf, sendcount)
152 if(sendbuf != MPI_IN_PLACE){
153 CHECK_TYPE(3, sendtype)
154 CHECK_COUNT(2, sendcount)
156 if(comm->rank() == root){
157 CHECK_TYPE(6, recvtype)
158 CHECK_NULL(5, MPI_ERR_COUNT, recvcounts)
159 CHECK_NULL(6, MPI_ERR_ARG, displs)
164 if (comm->rank() == root){
165 for (int i = 0; i < comm->size(); i++) {
166 CHECK_COUNT(5, recvcounts[i])
167 CHECK_BUFFER(4,recvbuf,recvcounts[i])
172 const void* real_sendbuf = sendbuf;
173 int real_sendcount = sendcount;
174 MPI_Datatype real_sendtype = sendtype;
175 if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
177 real_sendtype = recvtype;
180 int rank = simgrid::s4u::this_actor::get_pid();
181 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
183 auto* trace_recvcounts = new std::vector<int>();
184 if (comm->rank() == root) {
185 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
186 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
189 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gatherv" : "PMPI_Igatherv",
190 new simgrid::instr::VarCollTIData(
191 request == MPI_REQUEST_IGNORED ? "gatherv" : "igatherv", root,
192 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
193 nullptr, dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(real_sendtype),
194 simgrid::smpi::Datatype::encode(recvtype)));
195 if (request == MPI_REQUEST_IGNORED)
196 simgrid::smpi::colls::gatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
199 simgrid::smpi::colls::igatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
200 root, comm, request);
202 TRACE_smpi_comm_out(rank);
207 int PMPI_Allgather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
208 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm){
209 return PMPI_Iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
212 int PMPI_Iallgather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
213 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
216 CHECK_BUFFER(1, sendbuf, sendcount)
217 CHECK_BUFFER(4, recvbuf, recvcount)
218 if(sendbuf != MPI_IN_PLACE){
219 CHECK_COUNT(2, sendcount)
220 CHECK_TYPE(3, sendtype)
222 CHECK_TYPE(6, recvtype)
223 CHECK_COUNT(5, recvcount)
227 if (sendbuf == MPI_IN_PLACE) {
228 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * recvcount * comm->rank();
229 sendcount = recvcount;
232 int rank = simgrid::s4u::this_actor::get_pid();
234 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgather" : "PMPI_Iallggather",
235 new simgrid::instr::CollTIData(
236 request == MPI_REQUEST_IGNORED ? "allgather" : "iallgather", -1, -1.0,
237 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
238 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
239 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
240 if (request == MPI_REQUEST_IGNORED)
241 simgrid::smpi::colls::allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
243 simgrid::smpi::colls::iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, request);
245 TRACE_smpi_comm_out(rank);
250 int PMPI_Allgatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
251 void *recvbuf, const int *recvcounts, const int *displs, MPI_Datatype recvtype, MPI_Comm comm){
252 return PMPI_Iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, MPI_REQUEST_IGNORED);
255 int PMPI_Iallgatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
256 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
259 CHECK_BUFFER(1, sendbuf, sendcount)
260 if(sendbuf != MPI_IN_PLACE)
261 CHECK_TYPE(3, sendtype)
262 CHECK_TYPE(6, recvtype)
263 CHECK_NULL(5, MPI_ERR_COUNT, recvcounts)
264 CHECK_NULL(6, MPI_ERR_ARG, displs)
265 if(sendbuf != MPI_IN_PLACE)
266 CHECK_COUNT(2, sendcount)
268 for (int i = 0; i < comm->size(); i++) {
269 CHECK_COUNT(5, recvcounts[i])
270 CHECK_BUFFER(4, recvbuf, recvcounts[i])
274 if (sendbuf == MPI_IN_PLACE) {
275 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * displs[comm->rank()];
276 sendcount = recvcounts[comm->rank()];
279 int rank = simgrid::s4u::this_actor::get_pid();
280 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
282 auto* 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);
288 rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgatherv" : "PMPI_Iallgatherv",
289 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "allgatherv" : "iallgatherv", -1,
290 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(), nullptr,
291 dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
292 simgrid::smpi::Datatype::encode(recvtype)));
293 if (request == MPI_REQUEST_IGNORED)
294 simgrid::smpi::colls::allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
296 simgrid::smpi::colls::iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm,
299 TRACE_smpi_comm_out(rank);
304 int PMPI_Scatter(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
305 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
306 return PMPI_Iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
309 int PMPI_Iscatter(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
310 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
313 if(comm->rank() == root){
314 CHECK_BUFFER(1, sendbuf, sendcount)
315 CHECK_COUNT(2, sendcount)
316 CHECK_TYPE(3, sendtype)
318 if(recvbuf != MPI_IN_PLACE){
319 CHECK_BUFFER(4, recvbuf, recvcount)
320 CHECK_COUNT(5, recvcount)
321 CHECK_TYPE(6, recvtype)
327 if (recvbuf == MPI_IN_PLACE) {
329 recvcount = sendcount;
331 int rank = simgrid::s4u::this_actor::get_pid();
333 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatter" : "PMPI_Iscatter",
334 new simgrid::instr::CollTIData(
335 request == MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
336 (comm->rank() != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
337 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
338 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
339 if (request == MPI_REQUEST_IGNORED)
340 simgrid::smpi::colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
342 simgrid::smpi::colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
344 TRACE_smpi_comm_out(rank);
349 int PMPI_Scatterv(const void *sendbuf, const int *sendcounts, const int *displs,
350 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
351 return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
354 int PMPI_Iscatterv(const void* sendbuf, const int* sendcounts, const int* displs, MPI_Datatype sendtype, void* recvbuf, int recvcount,
355 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
358 if(recvbuf != MPI_IN_PLACE){
359 CHECK_BUFFER(4, recvbuf, recvcount)
360 CHECK_COUNT(5, recvcount)
361 CHECK_TYPE(7, recvtype)
365 if (comm->rank() == root) {
366 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
367 CHECK_NULL(3, MPI_ERR_ARG, displs)
368 CHECK_TYPE(4, sendtype)
369 for (int i = 0; i < comm->size(); i++){
370 CHECK_BUFFER(1, sendbuf, sendcounts[i])
371 CHECK_COUNT(2, sendcounts[i])
373 if (recvbuf == MPI_IN_PLACE) {
375 recvcount = sendcounts[comm->rank()];
381 int rank = simgrid::s4u::this_actor::get_pid();
382 int dt_size_send = sendtype->is_replayable() ? 1 : sendtype->size();
384 auto* trace_sendcounts = new std::vector<int>();
385 if (comm->rank() == root) {
386 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
387 trace_sendcounts->push_back(sendcounts[i] * dt_size_send);
391 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatterv" : "PMPI_Iscatterv",
392 new simgrid::instr::VarCollTIData(
393 request == MPI_REQUEST_IGNORED ? "scatterv" : "iscatterv", root, dt_size_send,
394 trace_sendcounts, recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
395 nullptr, simgrid::smpi::Datatype::encode(sendtype),
396 simgrid::smpi::Datatype::encode(recvtype)));
397 if (request == MPI_REQUEST_IGNORED)
398 simgrid::smpi::colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
400 simgrid::smpi::colls::iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm,
403 TRACE_smpi_comm_out(rank);
408 int PMPI_Reduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
410 return PMPI_Ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, MPI_REQUEST_IGNORED);
413 int PMPI_Ireduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm, MPI_Request* request)
416 CHECK_BUFFER(1, sendbuf, count)
417 if(comm->rank() == root)
418 CHECK_BUFFER(5, recvbuf, count)
419 CHECK_TYPE(4, datatype)
420 CHECK_COUNT(3, count)
426 int rank = simgrid::s4u::this_actor::get_pid();
428 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce" : "PMPI_Ireduce",
429 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "reduce" : "ireduce", root, 0,
430 datatype->is_replayable() ? count : count * datatype->size(), -1,
431 simgrid::smpi::Datatype::encode(datatype), ""));
432 if (request == MPI_REQUEST_IGNORED)
433 simgrid::smpi::colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
435 simgrid::smpi::colls::ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, request);
437 TRACE_smpi_comm_out(rank);
442 int PMPI_Reduce_local(const void* inbuf, void* inoutbuf, int count, MPI_Datatype datatype, MPI_Op op)
444 CHECK_BUFFER(1, inbuf, count)
445 CHECK_BUFFER(2, inoutbuf, count)
446 CHECK_TYPE(4, datatype)
447 CHECK_COUNT(3, count)
451 op->apply(inbuf, inoutbuf, &count, datatype);
456 int PMPI_Allreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
458 return PMPI_Iallreduce(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
461 int PMPI_Iallreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
464 CHECK_BUFFER(1, sendbuf, count)
465 CHECK_BUFFER(2, recvbuf, count)
466 CHECK_TYPE(4, datatype)
467 CHECK_COUNT(3, count)
472 std::vector<unsigned char> tmp_sendbuf;
473 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
475 int rank = simgrid::s4u::this_actor::get_pid();
477 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allreduce" : "PMPI_Iallreduce",
478 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "allreduce" : "iallreduce", -1, 0,
479 datatype->is_replayable() ? count : count * datatype->size(), -1,
480 simgrid::smpi::Datatype::encode(datatype), ""));
482 if (request == MPI_REQUEST_IGNORED)
483 simgrid::smpi::colls::allreduce(real_sendbuf, recvbuf, count, datatype, op, comm);
485 simgrid::smpi::colls::iallreduce(real_sendbuf, recvbuf, count, datatype, op, comm, request);
487 TRACE_smpi_comm_out(rank);
492 int PMPI_Scan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
494 return PMPI_Iscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
497 int PMPI_Iscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request)
500 CHECK_BUFFER(1,sendbuf,count)
501 CHECK_BUFFER(2,recvbuf,count)
502 CHECK_TYPE(4, datatype)
503 CHECK_COUNT(3, count)
508 int rank = simgrid::s4u::this_actor::get_pid();
509 std::vector<unsigned char> tmp_sendbuf;
510 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
512 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scan" : "PMPI_Iscan",
513 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "scan" : "iscan", -1,
514 datatype->is_replayable() ? count : count * datatype->size(),
515 simgrid::smpi::Datatype::encode(datatype)));
518 if (request == MPI_REQUEST_IGNORED)
519 retval = simgrid::smpi::colls::scan(real_sendbuf, recvbuf, count, datatype, op, comm);
521 retval = simgrid::smpi::colls::iscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
523 TRACE_smpi_comm_out(rank);
528 int PMPI_Exscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
530 return PMPI_Iexscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
533 int PMPI_Iexscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request){
535 CHECK_BUFFER(1, sendbuf, count)
536 CHECK_BUFFER(2, recvbuf, count)
537 CHECK_TYPE(4, datatype)
538 CHECK_COUNT(3, count)
543 int rank = simgrid::s4u::this_actor::get_pid();
544 std::vector<unsigned char> tmp_sendbuf;
545 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
547 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Exscan" : "PMPI_Iexscan",
548 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "exscan" : "iexscan", -1,
549 datatype->is_replayable() ? count : count * datatype->size(),
550 simgrid::smpi::Datatype::encode(datatype)));
553 if (request == MPI_REQUEST_IGNORED)
554 retval = simgrid::smpi::colls::exscan(real_sendbuf, recvbuf, count, datatype, op, comm);
556 retval = simgrid::smpi::colls::iexscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
558 TRACE_smpi_comm_out(rank);
563 int PMPI_Reduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
565 return PMPI_Ireduce_scatter(sendbuf, recvbuf, recvcounts, datatype, op, comm, MPI_REQUEST_IGNORED);
568 int PMPI_Ireduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
571 CHECK_TYPE(4, datatype)
572 CHECK_NULL(3, MPI_ERR_COUNT, recvcounts)
575 for (int i = 0; i < comm->size(); i++) {
576 CHECK_COUNT(3, recvcounts[i])
577 CHECK_BUFFER(1, sendbuf, recvcounts[i])
578 CHECK_BUFFER(2, recvbuf, recvcounts[i])
582 int rank = simgrid::s4u::this_actor::get_pid();
583 auto* trace_recvcounts = new std::vector<int>();
584 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
587 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
588 trace_recvcounts->push_back(recvcounts[i] * dt_send_size);
589 totalcount += recvcounts[i];
591 std::vector<unsigned char> tmp_sendbuf;
592 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, totalcount, datatype);
594 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter" : "PMPI_Ireduce_scatter",
595 new simgrid::instr::VarCollTIData(
596 request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, dt_send_size, nullptr,
597 -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
599 if (request == MPI_REQUEST_IGNORED)
600 simgrid::smpi::colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm);
602 simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm, request);
604 TRACE_smpi_comm_out(rank);
609 int PMPI_Reduce_scatter_block(const void *sendbuf, void *recvbuf, int recvcount,
610 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
612 return PMPI_Ireduce_scatter_block(sendbuf, recvbuf, recvcount, datatype, op, comm, MPI_REQUEST_IGNORED);
615 int PMPI_Ireduce_scatter_block(const void* sendbuf, void* recvbuf, int recvcount, MPI_Datatype datatype, MPI_Op op,
616 MPI_Comm comm, MPI_Request* request)
619 CHECK_BUFFER(1, sendbuf, recvcount)
620 CHECK_BUFFER(2, recvbuf, recvcount)
621 CHECK_TYPE(4, datatype)
622 CHECK_COUNT(3, recvcount)
627 int count = comm->size();
629 int rank = simgrid::s4u::this_actor::get_pid();
630 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
631 auto* trace_recvcounts = new std::vector<int>(recvcount * dt_send_size); // copy data to avoid bad free
632 std::vector<unsigned char> tmp_sendbuf;
633 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, recvcount * count, datatype);
636 rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter_block" : "PMPI_Ireduce_scatter_block",
637 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, 0,
638 nullptr, -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
640 std::vector<int> recvcounts(count);
641 for (int i = 0; i < count; i++)
642 recvcounts[i] = recvcount;
643 if (request == MPI_REQUEST_IGNORED)
644 simgrid::smpi::colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts.data(), datatype, op, comm);
646 simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts.data(), datatype, op, comm, request);
648 TRACE_smpi_comm_out(rank);
653 int PMPI_Alltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
654 MPI_Datatype recvtype, MPI_Comm comm){
655 return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
658 int PMPI_Ialltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
659 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
662 CHECK_BUFFER(1, sendbuf, sendcount)
663 CHECK_BUFFER(4, recvbuf, recvcount)
664 if(sendbuf != MPI_IN_PLACE)
665 CHECK_TYPE(3, sendtype)
666 CHECK_TYPE(6, recvtype)
667 CHECK_COUNT(5, recvcount)
668 if(sendbuf != MPI_IN_PLACE)
669 CHECK_COUNT(2, sendcount)
670 CHECK_COUNT(5, recvcount)
674 int rank = simgrid::s4u::this_actor::get_pid();
675 int real_sendcount = sendcount;
676 MPI_Datatype real_sendtype = sendtype;
678 std::vector<unsigned char> tmp_sendbuf;
679 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, recvcount * comm->size(), recvtype);
681 if (sendbuf == MPI_IN_PLACE) {
682 real_sendcount = recvcount;
683 real_sendtype = recvtype;
686 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoall" : "PMPI_Ialltoall",
687 new simgrid::instr::CollTIData(
688 request == MPI_REQUEST_IGNORED ? "alltoall" : "ialltoall", -1, -1.0,
689 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
690 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
691 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
693 if (request == MPI_REQUEST_IGNORED)
695 simgrid::smpi::colls::alltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, comm);
697 retval = simgrid::smpi::colls::ialltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype,
700 TRACE_smpi_comm_out(rank);
705 int PMPI_Alltoallv(const void* sendbuf, const int* sendcounts, const int* senddispls, MPI_Datatype sendtype, void* recvbuf,
706 const int* recvcounts, const int* recvdispls, MPI_Datatype recvtype, MPI_Comm comm)
708 return PMPI_Ialltoallv(sendbuf, sendcounts, senddispls, sendtype, recvbuf, recvcounts, recvdispls, recvtype, comm, MPI_REQUEST_IGNORED);
711 int PMPI_Ialltoallv(const void* sendbuf, const int* sendcounts, const int* senddispls, MPI_Datatype sendtype, void* recvbuf,
712 const int* recvcounts, const int* recvdispls, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
715 if(sendbuf != MPI_IN_PLACE){
716 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
717 CHECK_NULL(3, MPI_ERR_ARG, senddispls)
718 CHECK_TYPE(4, sendtype)
720 CHECK_TYPE(8, recvtype)
721 CHECK_NULL(6, MPI_ERR_COUNT, recvcounts)
722 CHECK_NULL(7, MPI_ERR_ARG, recvdispls)
725 int rank = simgrid::s4u::this_actor::get_pid();
726 int size = comm->size();
727 for (int i = 0; i < size; i++) {
728 if(sendbuf != MPI_IN_PLACE){
729 CHECK_BUFFER(1, sendbuf, sendcounts[i])
730 CHECK_COUNT(2, sendcounts[i])
732 CHECK_BUFFER(5, recvbuf, recvcounts[i])
733 CHECK_COUNT(6, recvcounts[i])
739 auto* trace_sendcounts = new std::vector<int>();
740 auto* trace_recvcounts = new std::vector<int>();
741 int dt_size_recv = recvtype->size();
743 const int* real_sendcounts = sendcounts;
744 const int* real_senddispls = senddispls;
745 MPI_Datatype real_sendtype = sendtype;
747 for (int i = 0; i < size; i++) { // copy data to avoid bad free
748 recv_size += recvcounts[i] * dt_size_recv;
749 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
750 if (((recvdispls[i] + recvcounts[i]) * dt_size_recv) > maxsize)
751 maxsize = (recvdispls[i] + recvcounts[i]) * dt_size_recv;
754 std::vector<unsigned char> tmp_sendbuf;
755 std::vector<int> tmp_sendcounts;
756 std::vector<int> tmp_senddispls;
757 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
758 if (sendbuf == MPI_IN_PLACE) {
759 tmp_sendcounts.assign(recvcounts, recvcounts + size);
760 real_sendcounts = tmp_sendcounts.data();
761 tmp_senddispls.assign(recvdispls, recvdispls + size);
762 real_senddispls = tmp_senddispls.data();
763 real_sendtype = recvtype;
766 int dt_size_send = real_sendtype->size();
768 for (int i = 0; i < size; i++) { // copy data to avoid bad free
769 send_size += real_sendcounts[i] * dt_size_send;
770 trace_sendcounts->push_back(real_sendcounts[i] * dt_size_send);
773 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallv" : "PMPI_Ialltoallv",
774 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
775 send_size, trace_sendcounts, recv_size, trace_recvcounts,
776 simgrid::smpi::Datatype::encode(real_sendtype),
777 simgrid::smpi::Datatype::encode(recvtype)));
780 if (request == MPI_REQUEST_IGNORED)
781 retval = simgrid::smpi::colls::alltoallv(real_sendbuf, real_sendcounts, real_senddispls, real_sendtype, recvbuf,
782 recvcounts, recvdispls, recvtype, comm);
784 retval = simgrid::smpi::colls::ialltoallv(real_sendbuf, real_sendcounts, real_senddispls, real_sendtype, recvbuf,
785 recvcounts, recvdispls, recvtype, comm, request);
787 TRACE_smpi_comm_out(rank);
792 int PMPI_Alltoallw(const void* sendbuf, const int* sendcounts, const int* senddispls, const MPI_Datatype* sendtypes, void* recvbuf,
793 const int* recvcounts, const int* recvdispls, const MPI_Datatype* recvtypes, MPI_Comm comm)
795 return PMPI_Ialltoallw(sendbuf, sendcounts, senddispls, sendtypes, recvbuf, recvcounts, recvdispls, recvtypes, comm, MPI_REQUEST_IGNORED);
798 int PMPI_Ialltoallw(const void* sendbuf, const int* sendcounts, const int* senddispls, const MPI_Datatype* sendtypes, void* recvbuf,
799 const int* recvcounts, const int* recvdispls, const MPI_Datatype* recvtypes, MPI_Comm comm, MPI_Request* request)
802 if(sendbuf != MPI_IN_PLACE){
803 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
804 CHECK_NULL(3, MPI_ERR_ARG, senddispls)
805 CHECK_NULL(4, MPI_ERR_TYPE, sendtypes)
807 CHECK_NULL(6, MPI_ERR_COUNT, recvcounts)
808 CHECK_NULL(7, MPI_ERR_ARG, recvdispls)
809 CHECK_NULL(8, MPI_ERR_TYPE, recvtypes)
811 int rank = simgrid::s4u::this_actor::get_pid();
812 int size = comm->size();
813 for (int i = 0; i < size; i++) {
814 if(sendbuf != MPI_IN_PLACE){
815 CHECK_BUFFER(1, sendbuf, sendcounts[i])
816 CHECK_COUNT(2, sendcounts[i])
817 CHECK_TYPE(4, sendtypes[i])
819 CHECK_BUFFER(5, recvbuf, recvcounts[i])
820 CHECK_COUNT(6, recvcounts[i])
821 CHECK_TYPE(8, recvtypes[i])
828 auto* trace_sendcounts = new std::vector<int>();
829 auto* trace_recvcounts = new std::vector<int>();
831 const int* real_sendcounts = sendcounts;
832 const int* real_senddispls = senddispls;
833 const MPI_Datatype* real_sendtypes = sendtypes;
834 unsigned long maxsize = 0;
835 for (int i = 0; i < size; i++) { // copy data to avoid bad free
836 if (recvtypes[i] == MPI_DATATYPE_NULL) {
837 delete trace_recvcounts;
838 delete trace_sendcounts;
841 recv_size += recvcounts[i] * recvtypes[i]->size();
842 trace_recvcounts->push_back(recvcounts[i] * recvtypes[i]->size());
843 if ((recvdispls[i] + (recvcounts[i] * recvtypes[i]->size())) > maxsize)
844 maxsize = recvdispls[i] + (recvcounts[i] * recvtypes[i]->size());
847 std::vector<unsigned char> tmp_sendbuf;
848 std::vector<int> tmp_sendcounts;
849 std::vector<int> tmp_senddispls;
850 std::vector<MPI_Datatype> tmp_sendtypes;
851 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
852 if (sendbuf == MPI_IN_PLACE) {
853 tmp_sendcounts.assign(recvcounts, recvcounts + size);
854 real_sendcounts = tmp_sendcounts.data();
855 tmp_senddispls.assign(recvdispls, recvdispls + size);
856 real_senddispls = tmp_senddispls.data();
857 tmp_sendtypes.assign(recvtypes, recvtypes + size);
858 real_sendtypes = tmp_sendtypes.data();
861 for (int i = 0; i < size; i++) { // copy data to avoid bad free
862 send_size += real_sendcounts[i] * real_sendtypes[i]->size();
863 trace_sendcounts->push_back(real_sendcounts[i] * real_sendtypes[i]->size());
866 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallw" : "PMPI_Ialltoallw",
867 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
868 send_size, trace_sendcounts, recv_size, trace_recvcounts,
869 simgrid::smpi::Datatype::encode(real_sendtypes[0]),
870 simgrid::smpi::Datatype::encode(recvtypes[0])));
873 if (request == MPI_REQUEST_IGNORED)
874 retval = simgrid::smpi::colls::alltoallw(real_sendbuf, real_sendcounts, real_senddispls, real_sendtypes, recvbuf,
875 recvcounts, recvdispls, recvtypes, comm);
877 retval = simgrid::smpi::colls::ialltoallw(real_sendbuf, real_sendcounts, real_senddispls, real_sendtypes, recvbuf,
878 recvcounts, recvdispls, recvtypes, comm, request);
880 TRACE_smpi_comm_out(rank);