1 /* Copyright (c) 2007-2021. 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)
68 CHECK_TYPE(3, datatype)
69 CHECK_BUFFER(1, buf, count, datatype)
74 int pid = simgrid::s4u::this_actor::get_pid();
75 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Bcast" : "PMPI_Ibcast",
76 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "bcast" : "ibcast", root, -1.0,
77 datatype->is_replayable() ? count : count * datatype->size(), -1,
78 simgrid::smpi::Datatype::encode(datatype), ""));
79 if (comm->size() > 1) {
80 if (request == MPI_REQUEST_IGNORED)
81 simgrid::smpi::colls::bcast(buf, count, datatype, root, comm);
83 simgrid::smpi::colls::ibcast(buf, count, datatype, root, comm, request);
85 if (request != MPI_REQUEST_IGNORED)
86 *request = MPI_REQUEST_NULL;
89 TRACE_smpi_comm_out(pid);
94 int PMPI_Gather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
95 int root, MPI_Comm comm){
96 return PMPI_Igather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
99 int PMPI_Igather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
100 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
105 int rank = comm->rank();
106 if(sendbuf != MPI_IN_PLACE){
107 CHECK_COUNT(2, sendcount)
108 CHECK_TYPE(3, sendtype)
109 CHECK_BUFFER(1,sendbuf, sendcount, sendtype)
112 CHECK_NOT_IN_PLACE_ROOT(4, recvbuf)
113 CHECK_TYPE(6, recvtype)
114 CHECK_COUNT(5, recvcount)
115 CHECK_BUFFER(4, recvbuf, recvcount, recvtype)
117 CHECK_NOT_IN_PLACE_ROOT(1, sendbuf)
122 const void* real_sendbuf = sendbuf;
123 int real_sendcount = sendcount;
124 MPI_Datatype real_sendtype = sendtype;
126 if (sendbuf == MPI_IN_PLACE) {
128 real_sendtype = recvtype;
129 } else if(recvtype->size() * recvcount != sendtype->size() * sendcount){
130 XBT_WARN("MPI_(I)Gather : received size at root differs from sent size : %zu vs %zu", recvtype->size() * recvcount , sendtype->size() * sendcount);
131 return MPI_ERR_TRUNCATE;
137 int pid = simgrid::s4u::this_actor::get_pid();
139 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Gather" : "PMPI_Igather",
140 new simgrid::instr::CollTIData(
141 request == MPI_REQUEST_IGNORED ? "gather" : "igather", root, -1.0,
142 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
143 (comm->rank() != root || recvtype->is_replayable()) ? recvcount : recvcount * recvtype->size(),
144 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
145 if (request == MPI_REQUEST_IGNORED)
146 simgrid::smpi::colls::gather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm);
148 simgrid::smpi::colls::igather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm,
151 TRACE_smpi_comm_out(pid);
156 int PMPI_Gatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, const int *recvcounts, const int *displs,
157 MPI_Datatype recvtype, int root, MPI_Comm comm){
158 return PMPI_Igatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, root, comm, MPI_REQUEST_IGNORED);
161 int PMPI_Igatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
162 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
167 int rank = comm->rank();
168 if(sendbuf != MPI_IN_PLACE){
169 CHECK_TYPE(3, sendtype)
170 CHECK_COUNT(2, sendcount)
172 CHECK_BUFFER(1, sendbuf, sendcount, sendtype)
174 CHECK_NOT_IN_PLACE_ROOT(4, recvbuf)
175 CHECK_TYPE(6, recvtype)
176 CHECK_NULL(5, MPI_ERR_COUNT, recvcounts)
177 CHECK_NULL(6, MPI_ERR_ARG, displs)
179 CHECK_NOT_IN_PLACE_ROOT(1, sendbuf)
185 for (int i = 0; i < comm->size(); i++) {
186 CHECK_COUNT(5, recvcounts[i])
187 CHECK_BUFFER(4,recvbuf,recvcounts[i], recvtype)
192 const void* real_sendbuf = sendbuf;
193 int real_sendcount = sendcount;
194 MPI_Datatype real_sendtype = sendtype;
195 if ((rank == root) && (sendbuf == MPI_IN_PLACE)) {
197 real_sendtype = recvtype;
200 int pid = simgrid::s4u::this_actor::get_pid();
201 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
203 auto trace_recvcounts = std::make_shared<std::vector<int>>();
205 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
206 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
209 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Gatherv" : "PMPI_Igatherv",
210 new simgrid::instr::VarCollTIData(
211 request == MPI_REQUEST_IGNORED ? "gatherv" : "igatherv", root,
212 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
213 nullptr, dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(real_sendtype),
214 simgrid::smpi::Datatype::encode(recvtype)));
215 if (request == MPI_REQUEST_IGNORED)
216 simgrid::smpi::colls::gatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
219 simgrid::smpi::colls::igatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
220 root, comm, request);
222 TRACE_smpi_comm_out(pid);
227 int PMPI_Allgather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
228 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm){
229 return PMPI_Iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
232 int PMPI_Iallgather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
233 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
238 int rank = comm->rank();
239 CHECK_NOT_IN_PLACE(4, recvbuf)
240 if(sendbuf != MPI_IN_PLACE){
241 CHECK_COUNT(2, sendcount)
242 CHECK_TYPE(3, sendtype)
244 CHECK_TYPE(6, recvtype)
245 CHECK_COUNT(5, recvcount)
246 CHECK_BUFFER(1, sendbuf, sendcount, sendtype)
247 CHECK_BUFFER(4, recvbuf, recvcount, recvtype)
250 if (sendbuf == MPI_IN_PLACE) {
251 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * recvcount * comm->rank();
252 sendcount = recvcount;
256 if(recvtype->size() * recvcount != sendtype->size() * sendcount){
257 XBT_WARN("MPI_(I)Allgather : received size from each process differs from sent size : %zu vs %zu", recvtype->size() * recvcount, sendtype->size() * sendcount);
258 return MPI_ERR_TRUNCATE;
263 int pid = simgrid::s4u::this_actor::get_pid();
265 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Allgather" : "PMPI_Iallggather",
266 new simgrid::instr::CollTIData(
267 request == MPI_REQUEST_IGNORED ? "allgather" : "iallgather", -1, -1.0,
268 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
269 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
270 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
271 if (request == MPI_REQUEST_IGNORED)
272 simgrid::smpi::colls::allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
274 simgrid::smpi::colls::iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, request);
276 TRACE_smpi_comm_out(pid);
281 int PMPI_Allgatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
282 void *recvbuf, const int *recvcounts, const int *displs, MPI_Datatype recvtype, MPI_Comm comm){
283 return PMPI_Iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, MPI_REQUEST_IGNORED);
286 int PMPI_Iallgatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
287 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
292 int rank = comm->rank();
293 if(sendbuf != MPI_IN_PLACE)
294 CHECK_TYPE(3, sendtype)
295 CHECK_TYPE(6, recvtype)
296 CHECK_NULL(5, MPI_ERR_COUNT, recvcounts)
297 CHECK_NULL(6, MPI_ERR_ARG, displs)
298 if(sendbuf != MPI_IN_PLACE){
299 CHECK_COUNT(2, sendcount)
300 CHECK_BUFFER(1, sendbuf, sendcount, sendtype)
303 CHECK_NOT_IN_PLACE(4, recvbuf)
304 for (int i = 0; i < comm->size(); i++) {
305 CHECK_COUNT(5, recvcounts[i])
306 CHECK_BUFFER(4, recvbuf, recvcounts[i], recvtype)
310 if (sendbuf == MPI_IN_PLACE) {
311 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * displs[comm->rank()];
312 sendcount = recvcounts[comm->rank()];
315 int pid = simgrid::s4u::this_actor::get_pid();
316 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
318 auto trace_recvcounts = std::make_shared<std::vector<int>>();
319 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
320 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
324 pid, request == MPI_REQUEST_IGNORED ? "PMPI_Allgatherv" : "PMPI_Iallgatherv",
325 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "allgatherv" : "iallgatherv", -1,
326 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(), nullptr,
327 dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
328 simgrid::smpi::Datatype::encode(recvtype)));
329 if (request == MPI_REQUEST_IGNORED)
330 simgrid::smpi::colls::allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
332 simgrid::smpi::colls::iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm,
335 TRACE_smpi_comm_out(pid);
340 int PMPI_Scatter(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
341 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
342 return PMPI_Iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
345 int PMPI_Iscatter(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
346 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
351 int rank = comm->rank();
353 CHECK_NOT_IN_PLACE_ROOT(1, sendbuf)
354 CHECK_COUNT(2, sendcount)
355 CHECK_TYPE(3, sendtype)
356 CHECK_BUFFER(1, sendbuf, sendcount, sendtype)
358 CHECK_NOT_IN_PLACE_ROOT(4, recvbuf)
360 if(recvbuf != MPI_IN_PLACE){
361 CHECK_COUNT(5, recvcount)
362 CHECK_TYPE(6, recvtype)
363 CHECK_BUFFER(4, recvbuf, recvcount, recvtype)
368 if (recvbuf == MPI_IN_PLACE) {
370 recvcount = sendcount;
373 if((rank == root) && (recvtype->size() * recvcount != sendtype->size() * sendcount)){
374 XBT_WARN("MPI_(I)Scatter : sent size to each process differs from receive size");
375 return MPI_ERR_TRUNCATE;
380 int pid = simgrid::s4u::this_actor::get_pid();
382 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Scatter" : "PMPI_Iscatter",
383 new simgrid::instr::CollTIData(
384 request == MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
385 (rank != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
386 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
387 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
388 if (request == MPI_REQUEST_IGNORED)
389 simgrid::smpi::colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
391 simgrid::smpi::colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
393 TRACE_smpi_comm_out(pid);
398 int PMPI_Scatterv(const void *sendbuf, const int *sendcounts, const int *displs,
399 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
400 return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
403 int PMPI_Iscatterv(const void* sendbuf, const int* sendcounts, const int* displs, MPI_Datatype sendtype, void* recvbuf, int recvcount,
404 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
409 int rank = comm->rank();
410 if(recvbuf != MPI_IN_PLACE){
411 CHECK_NOT_IN_PLACE_ROOT(1, sendbuf)
412 CHECK_COUNT(5, recvcount)
413 CHECK_TYPE(7, recvtype)
414 CHECK_BUFFER(4, recvbuf, recvcount, recvtype)
419 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
420 CHECK_NULL(3, MPI_ERR_ARG, displs)
421 CHECK_TYPE(4, sendtype)
422 for (int i = 0; i < comm->size(); i++){
423 CHECK_COUNT(2, sendcounts[i])
424 CHECK_BUFFER(1, sendbuf, sendcounts[i], sendtype)
426 if (recvbuf == MPI_IN_PLACE) {
428 recvcount = sendcounts[rank];
431 CHECK_NOT_IN_PLACE_ROOT(4, recvbuf)
436 int pid = simgrid::s4u::this_actor::get_pid();
437 int dt_size_send = sendtype->is_replayable() ? 1 : sendtype->size();
439 auto trace_sendcounts = std::make_shared<std::vector<int>>();
441 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
442 trace_sendcounts->push_back(sendcounts[i] * dt_size_send);
446 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Scatterv" : "PMPI_Iscatterv",
447 new simgrid::instr::VarCollTIData(
448 request == MPI_REQUEST_IGNORED ? "scatterv" : "iscatterv", root, dt_size_send,
449 trace_sendcounts, recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
450 nullptr, simgrid::smpi::Datatype::encode(sendtype),
451 simgrid::smpi::Datatype::encode(recvtype)));
452 if (request == MPI_REQUEST_IGNORED)
453 simgrid::smpi::colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
455 simgrid::smpi::colls::iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm,
458 TRACE_smpi_comm_out(pid);
463 int PMPI_Reduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
465 return PMPI_Ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, MPI_REQUEST_IGNORED);
468 int PMPI_Ireduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm, MPI_Request* request)
473 int rank = comm->rank();
474 CHECK_TYPE(4, datatype)
475 CHECK_COUNT(3, count)
476 CHECK_BUFFER(1, sendbuf, count, datatype)
478 CHECK_NOT_IN_PLACE(2, recvbuf)
479 CHECK_BUFFER(5, recvbuf, count, datatype)
481 CHECK_OP(5, op, datatype)
486 int pid = simgrid::s4u::this_actor::get_pid();
488 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce" : "PMPI_Ireduce",
489 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "reduce" : "ireduce", root, 0,
490 datatype->is_replayable() ? count : count * datatype->size(), -1,
491 simgrid::smpi::Datatype::encode(datatype), ""));
492 if (request == MPI_REQUEST_IGNORED)
493 simgrid::smpi::colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
495 simgrid::smpi::colls::ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, request);
497 TRACE_smpi_comm_out(pid);
502 int PMPI_Reduce_local(const void* inbuf, void* inoutbuf, int count, MPI_Datatype datatype, MPI_Op op)
506 CHECK_TYPE(4, datatype)
507 CHECK_COUNT(3, count)
508 CHECK_BUFFER(1, inbuf, count, datatype)
509 CHECK_BUFFER(2, inoutbuf, count, datatype)
510 CHECK_OP(5, op, datatype)
513 op->apply(inbuf, inoutbuf, &count, datatype);
518 int PMPI_Allreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
520 return PMPI_Iallreduce(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
523 int PMPI_Iallreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
528 int rank = comm->rank();
529 CHECK_NOT_IN_PLACE(2, recvbuf)
530 CHECK_TYPE(4, datatype)
531 CHECK_OP(5, op, datatype)
532 CHECK_COUNT(3, count)
533 CHECK_BUFFER(1, sendbuf, count, datatype)
534 CHECK_BUFFER(2, recvbuf, count, datatype)
538 std::vector<unsigned char> tmp_sendbuf;
539 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
541 int pid = simgrid::s4u::this_actor::get_pid();
543 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Allreduce" : "PMPI_Iallreduce",
544 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "allreduce" : "iallreduce", -1, 0,
545 datatype->is_replayable() ? count : count * datatype->size(), -1,
546 simgrid::smpi::Datatype::encode(datatype), ""));
548 if (request == MPI_REQUEST_IGNORED)
549 simgrid::smpi::colls::allreduce(real_sendbuf, recvbuf, count, datatype, op, comm);
551 simgrid::smpi::colls::iallreduce(real_sendbuf, recvbuf, count, datatype, op, comm, request);
553 TRACE_smpi_comm_out(pid);
558 int PMPI_Scan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
560 return PMPI_Iscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
563 int PMPI_Iscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request)
568 CHECK_TYPE(4, datatype)
569 CHECK_COUNT(3, count)
570 CHECK_BUFFER(1,sendbuf,count, datatype)
571 CHECK_BUFFER(2,recvbuf,count, datatype)
573 CHECK_OP(5, op, datatype)
576 int pid = simgrid::s4u::this_actor::get_pid();
577 std::vector<unsigned char> tmp_sendbuf;
578 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
580 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Scan" : "PMPI_Iscan",
581 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "scan" : "iscan", -1,
582 datatype->is_replayable() ? count : count * datatype->size(),
583 simgrid::smpi::Datatype::encode(datatype)));
586 if (request == MPI_REQUEST_IGNORED)
587 retval = simgrid::smpi::colls::scan(real_sendbuf, recvbuf, count, datatype, op, comm);
589 retval = simgrid::smpi::colls::iscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
591 TRACE_smpi_comm_out(pid);
596 int PMPI_Exscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
598 return PMPI_Iexscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
601 int PMPI_Iexscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request){
605 CHECK_TYPE(4, datatype)
606 CHECK_COUNT(3, count)
607 CHECK_BUFFER(1, sendbuf, count, datatype)
608 CHECK_BUFFER(2, recvbuf, count, datatype)
610 CHECK_OP(5, op, datatype)
613 int pid = simgrid::s4u::this_actor::get_pid();
614 std::vector<unsigned char> tmp_sendbuf;
615 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
617 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Exscan" : "PMPI_Iexscan",
618 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "exscan" : "iexscan", -1,
619 datatype->is_replayable() ? count : count * datatype->size(),
620 simgrid::smpi::Datatype::encode(datatype)));
623 if (request == MPI_REQUEST_IGNORED)
624 retval = simgrid::smpi::colls::exscan(real_sendbuf, recvbuf, count, datatype, op, comm);
626 retval = simgrid::smpi::colls::iexscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
628 TRACE_smpi_comm_out(pid);
633 int PMPI_Reduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
635 return PMPI_Ireduce_scatter(sendbuf, recvbuf, recvcounts, datatype, op, comm, MPI_REQUEST_IGNORED);
638 int PMPI_Ireduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
643 int rank = comm->rank();
644 CHECK_NOT_IN_PLACE(2, recvbuf)
645 CHECK_TYPE(4, datatype)
646 CHECK_NULL(3, MPI_ERR_COUNT, recvcounts)
648 CHECK_OP(5, op, datatype)
649 for (int i = 0; i < comm->size(); i++) {
650 CHECK_COUNT(3, recvcounts[i])
651 CHECK_BUFFER(1, sendbuf, recvcounts[i], datatype)
652 CHECK_BUFFER(2, recvbuf, recvcounts[i], datatype)
656 int pid = simgrid::s4u::this_actor::get_pid();
657 auto trace_recvcounts = std::make_shared<std::vector<int>>();
658 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
661 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
662 trace_recvcounts->push_back(recvcounts[i] * dt_send_size);
663 totalcount += recvcounts[i];
665 std::vector<unsigned char> tmp_sendbuf;
666 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, totalcount, datatype);
668 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter" : "PMPI_Ireduce_scatter",
669 new simgrid::instr::VarCollTIData(
670 request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, dt_send_size, nullptr,
671 -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
673 if (request == MPI_REQUEST_IGNORED)
674 simgrid::smpi::colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm);
676 simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm, request);
678 TRACE_smpi_comm_out(pid);
683 int PMPI_Reduce_scatter_block(const void *sendbuf, void *recvbuf, int recvcount,
684 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
686 return PMPI_Ireduce_scatter_block(sendbuf, recvbuf, recvcount, datatype, op, comm, MPI_REQUEST_IGNORED);
689 int PMPI_Ireduce_scatter_block(const void* sendbuf, void* recvbuf, int recvcount, MPI_Datatype datatype, MPI_Op op,
690 MPI_Comm comm, MPI_Request* request)
695 CHECK_TYPE(4, datatype)
696 CHECK_COUNT(3, recvcount)
697 CHECK_BUFFER(1, sendbuf, recvcount, datatype)
698 CHECK_BUFFER(2, recvbuf, recvcount, datatype)
700 CHECK_OP(5, op, datatype)
703 int count = comm->size();
705 int pid = simgrid::s4u::this_actor::get_pid();
706 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
707 auto trace_recvcounts = std::make_shared<std::vector<int>>(recvcount * dt_send_size); // copy data to avoid bad free
708 std::vector<unsigned char> tmp_sendbuf;
709 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, recvcount * count, datatype);
712 pid, 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,
714 nullptr, -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
716 std::vector<int> recvcounts(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(real_sendbuf, recvbuf, recvcounts.data(), datatype, op, comm);
722 simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts.data(), datatype, op, comm, request);
724 TRACE_smpi_comm_out(pid);
729 int PMPI_Alltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
730 MPI_Datatype recvtype, MPI_Comm comm){
731 return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
734 int PMPI_Ialltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
735 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
740 if(sendbuf != MPI_IN_PLACE){
741 CHECK_TYPE(3, sendtype)
742 CHECK_COUNT(2, sendcount)
743 CHECK_BUFFER(1, sendbuf, sendcount, sendtype)
745 CHECK_TYPE(6, recvtype)
746 CHECK_COUNT(5, recvcount)
747 CHECK_COUNT(5, recvcount)
748 CHECK_BUFFER(4, recvbuf, recvcount, recvtype)
751 int pid = simgrid::s4u::this_actor::get_pid();
752 int real_sendcount = sendcount;
753 MPI_Datatype real_sendtype = sendtype;
755 std::vector<unsigned char> tmp_sendbuf;
756 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, recvcount * comm->size(), recvtype);
758 if (sendbuf == MPI_IN_PLACE) {
759 real_sendcount = recvcount;
760 real_sendtype = recvtype;
763 if(recvtype->size() * recvcount != real_sendtype->size() * real_sendcount){
764 XBT_WARN("MPI_(I)Alltoall : receive size from each process differs from sent size : %zu vs %zu", recvtype->size() * recvcount, real_sendtype->size() * real_sendcount);
765 return MPI_ERR_TRUNCATE;
770 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoall" : "PMPI_Ialltoall",
771 new simgrid::instr::CollTIData(
772 request == MPI_REQUEST_IGNORED ? "alltoall" : "ialltoall", -1, -1.0,
773 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
774 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
775 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
777 if (request == MPI_REQUEST_IGNORED)
779 simgrid::smpi::colls::alltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, comm);
781 retval = simgrid::smpi::colls::ialltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype,
784 TRACE_smpi_comm_out(pid);
789 int PMPI_Alltoallv(const void* sendbuf, const int* sendcounts, const int* senddispls, MPI_Datatype sendtype, void* recvbuf,
790 const int* recvcounts, const int* recvdispls, MPI_Datatype recvtype, MPI_Comm comm)
792 return PMPI_Ialltoallv(sendbuf, sendcounts, senddispls, sendtype, recvbuf, recvcounts, recvdispls, recvtype, comm, MPI_REQUEST_IGNORED);
795 int PMPI_Ialltoallv(const void* sendbuf, const int* sendcounts, const int* senddispls, MPI_Datatype sendtype, void* recvbuf,
796 const int* recvcounts, const int* recvdispls, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
801 if(sendbuf != MPI_IN_PLACE){
802 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
803 CHECK_NULL(3, MPI_ERR_ARG, senddispls)
804 CHECK_TYPE(4, sendtype)
806 CHECK_TYPE(8, recvtype)
807 CHECK_NULL(6, MPI_ERR_COUNT, recvcounts)
808 CHECK_NULL(7, MPI_ERR_ARG, recvdispls)
811 int pid = 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], sendtype)
816 CHECK_COUNT(2, sendcounts[i])
818 CHECK_BUFFER(5, recvbuf, recvcounts[i], recvtype)
819 CHECK_COUNT(6, recvcounts[i])
825 auto trace_sendcounts = std::make_shared<std::vector<int>>();
826 auto trace_recvcounts = std::make_shared<std::vector<int>>();
827 int dt_size_recv = recvtype->size();
829 const int* real_sendcounts = sendcounts;
830 const int* real_senddispls = senddispls;
831 MPI_Datatype real_sendtype = sendtype;
833 for (int i = 0; i < size; i++) { // copy data to avoid bad free
834 recv_size += recvcounts[i] * dt_size_recv;
835 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
836 if (((recvdispls[i] + recvcounts[i]) * dt_size_recv) > maxsize)
837 maxsize = (recvdispls[i] + recvcounts[i]) * dt_size_recv;
840 std::vector<unsigned char> tmp_sendbuf;
841 std::vector<int> tmp_sendcounts;
842 std::vector<int> tmp_senddispls;
843 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
844 if (sendbuf == MPI_IN_PLACE) {
845 tmp_sendcounts.assign(recvcounts, recvcounts + size);
846 real_sendcounts = tmp_sendcounts.data();
847 tmp_senddispls.assign(recvdispls, recvdispls + size);
848 real_senddispls = tmp_senddispls.data();
849 real_sendtype = recvtype;
852 if(recvtype->size() * recvcounts[comm->rank()] != real_sendtype->size() * real_sendcounts[comm->rank()]){
853 XBT_WARN("MPI_(I)Alltoallv : receive size from me differs from sent size to me : %zu vs %zu", recvtype->size() * recvcounts[comm->rank()], real_sendtype->size() * real_sendcounts[comm->rank()]);
855 return MPI_ERR_TRUNCATE;
858 int dt_size_send = real_sendtype->size();
860 for (int i = 0; i < size; i++) { // copy data to avoid bad free
861 send_size += real_sendcounts[i] * dt_size_send;
862 trace_sendcounts->push_back(real_sendcounts[i] * dt_size_send);
865 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallv" : "PMPI_Ialltoallv",
866 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
867 send_size, trace_sendcounts, recv_size, trace_recvcounts,
868 simgrid::smpi::Datatype::encode(real_sendtype),
869 simgrid::smpi::Datatype::encode(recvtype)));
872 if (request == MPI_REQUEST_IGNORED)
873 retval = simgrid::smpi::colls::alltoallv(real_sendbuf, real_sendcounts, real_senddispls, real_sendtype, recvbuf,
874 recvcounts, recvdispls, recvtype, comm);
876 retval = simgrid::smpi::colls::ialltoallv(real_sendbuf, real_sendcounts, real_senddispls, real_sendtype, recvbuf,
877 recvcounts, recvdispls, recvtype, comm, request);
879 TRACE_smpi_comm_out(pid);
884 int PMPI_Alltoallw(const void* sendbuf, const int* sendcounts, const int* senddispls, const MPI_Datatype* sendtypes, void* recvbuf,
885 const int* recvcounts, const int* recvdispls, const MPI_Datatype* recvtypes, MPI_Comm comm)
887 return PMPI_Ialltoallw(sendbuf, sendcounts, senddispls, sendtypes, recvbuf, recvcounts, recvdispls, recvtypes, comm, MPI_REQUEST_IGNORED);
890 int PMPI_Ialltoallw(const void* sendbuf, const int* sendcounts, const int* senddispls, const MPI_Datatype* sendtypes, void* recvbuf,
891 const int* recvcounts, const int* recvdispls, const MPI_Datatype* recvtypes, MPI_Comm comm, MPI_Request* request)
896 if(sendbuf != MPI_IN_PLACE){
897 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
898 CHECK_NULL(3, MPI_ERR_ARG, senddispls)
899 CHECK_NULL(4, MPI_ERR_TYPE, sendtypes)
901 CHECK_NULL(6, MPI_ERR_COUNT, recvcounts)
902 CHECK_NULL(7, MPI_ERR_ARG, recvdispls)
903 CHECK_NULL(8, MPI_ERR_TYPE, recvtypes)
905 int pid = simgrid::s4u::this_actor::get_pid();
906 int size = comm->size();
907 for (int i = 0; i < size; i++) {
908 if(sendbuf != MPI_IN_PLACE){
909 CHECK_COUNT(2, sendcounts[i])
910 CHECK_TYPE(4, sendtypes[i])
911 CHECK_BUFFER(1, sendbuf, sendcounts[i], sendtypes[i])
913 CHECK_COUNT(6, recvcounts[i])
914 CHECK_TYPE(8, recvtypes[i])
915 CHECK_BUFFER(5, recvbuf, recvcounts[i], recvtypes[i])
922 auto trace_sendcounts = std::make_shared<std::vector<int>>();
923 auto trace_recvcounts = std::make_shared<std::vector<int>>();
925 const int* real_sendcounts = sendcounts;
926 const int* real_senddispls = senddispls;
927 const MPI_Datatype* real_sendtypes = sendtypes;
928 unsigned long maxsize = 0;
929 for (int i = 0; i < size; i++) { // copy data to avoid bad free
930 if (recvtypes[i] == MPI_DATATYPE_NULL)
932 recv_size += recvcounts[i] * recvtypes[i]->size();
933 trace_recvcounts->push_back(recvcounts[i] * recvtypes[i]->size());
934 if ((recvdispls[i] + (recvcounts[i] * recvtypes[i]->size())) > maxsize)
935 maxsize = recvdispls[i] + (recvcounts[i] * recvtypes[i]->size());
938 std::vector<unsigned char> tmp_sendbuf;
939 std::vector<int> tmp_sendcounts;
940 std::vector<int> tmp_senddispls;
941 std::vector<MPI_Datatype> tmp_sendtypes;
942 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
943 if (sendbuf == MPI_IN_PLACE) {
944 tmp_sendcounts.assign(recvcounts, recvcounts + size);
945 real_sendcounts = tmp_sendcounts.data();
946 tmp_senddispls.assign(recvdispls, recvdispls + size);
947 real_senddispls = tmp_senddispls.data();
948 tmp_sendtypes.assign(recvtypes, recvtypes + size);
949 real_sendtypes = tmp_sendtypes.data();
953 if(recvtypes[comm->rank()]->size() * recvcounts[comm->rank()] != real_sendtypes[comm->rank()]->size() * real_sendcounts[comm->rank()]){
954 XBT_WARN("MPI_(I)Alltoallw : receive size from me differs from sent size to me : %zu vs %zu", recvtypes[comm->rank()]->size() * recvcounts[comm->rank()], real_sendtypes[comm->rank()]->size() * real_sendcounts[comm->rank()]);
956 return MPI_ERR_TRUNCATE;
959 for (int i = 0; i < size; i++) { // copy data to avoid bad free
960 send_size += real_sendcounts[i] * real_sendtypes[i]->size();
961 trace_sendcounts->push_back(real_sendcounts[i] * real_sendtypes[i]->size());
964 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallw" : "PMPI_Ialltoallw",
965 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
966 send_size, trace_sendcounts, recv_size, trace_recvcounts,
967 simgrid::smpi::Datatype::encode(real_sendtypes[0]),
968 simgrid::smpi::Datatype::encode(recvtypes[0])));
971 if (request == MPI_REQUEST_IGNORED)
972 retval = simgrid::smpi::colls::alltoallw(real_sendbuf, real_sendcounts, real_senddispls, real_sendtypes, recvbuf,
973 recvcounts, recvdispls, recvtypes, comm);
975 retval = simgrid::smpi::colls::ialltoallw(real_sendbuf, real_sendcounts, real_senddispls, real_sendtypes, recvbuf,
976 recvcounts, recvdispls, recvtypes, comm, request);
978 TRACE_smpi_comm_out(pid);