1 /* Copyright (c) 2007-2022. 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)
40 CHECK_COLLECTIVE(comm, request == MPI_REQUEST_IGNORED ? "PMPI_Barrier" : "PMPI_Ibarrier")
41 const SmpiBenchGuard suspend_bench;
42 aid_t pid = simgrid::s4u::this_actor::get_pid();
43 TRACE_smpi_comm_in(pid, 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(pid);
56 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
58 return PMPI_Ibcast(buf, count, datatype, root, comm, MPI_REQUEST_IGNORED);
61 int PMPI_Ibcast(void* buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm, MPI_Request* request)
66 CHECK_TYPE(3, datatype)
67 CHECK_BUFFER(1, buf, count, datatype)
70 CHECK_COLLECTIVE(comm, std::string(request == MPI_REQUEST_IGNORED ? "PMPI_Bcast" : "PMPI_Ibcast") + " with root " +
73 const SmpiBenchGuard suspend_bench;
74 aid_t 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,
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);
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)
103 int rank = comm->rank();
104 if(sendbuf != MPI_IN_PLACE){
105 CHECK_COUNT(2, sendcount)
106 CHECK_TYPE(3, sendtype)
107 CHECK_BUFFER(1,sendbuf, sendcount, sendtype)
111 CHECK_NOT_IN_PLACE_ROOT(4, recvbuf)
112 CHECK_TYPE(6, recvtype)
113 CHECK_COUNT(5, recvcount)
114 CHECK_BUFFER(4, recvbuf, recvcount, recvtype)
116 CHECK_NOT_IN_PLACE_ROOT(1, sendbuf)
120 CHECK_COLLECTIVE(comm, std::string(request == MPI_REQUEST_IGNORED ? "PMPI_Gather" : "PMPI_Igather") + +" with root " +
121 std::to_string(root))
123 const void* real_sendbuf = sendbuf;
124 int real_sendcount = sendcount;
125 MPI_Datatype real_sendtype = sendtype;
127 if (sendbuf == MPI_IN_PLACE) {
129 real_sendtype = recvtype;
130 } else if(recvtype->size() * recvcount != sendtype->size() * sendcount){
131 XBT_WARN("MPI_(I)Gather : received size at root differs from sent size : %zu vs %zu", recvtype->size() * recvcount , sendtype->size() * sendcount);
132 return MPI_ERR_TRUNCATE;
136 const SmpiBenchGuard suspend_bench;
138 aid_t pid = simgrid::s4u::this_actor::get_pid();
140 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Gather" : "PMPI_Igather",
141 new simgrid::instr::CollTIData(
142 request == MPI_REQUEST_IGNORED ? "gather" : "igather", root, -1.0,
143 real_sendcount, recvcount,
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);
155 int PMPI_Gatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, const int *recvcounts, const int *displs,
156 MPI_Datatype recvtype, int root, MPI_Comm comm){
157 return PMPI_Igatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, root, comm, MPI_REQUEST_IGNORED);
160 int PMPI_Igatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
161 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
165 int rank = comm->rank();
166 if(sendbuf != MPI_IN_PLACE){
167 CHECK_TYPE(3, sendtype)
168 CHECK_COUNT(2, sendcount)
170 CHECK_BUFFER(1, sendbuf, sendcount, sendtype)
173 CHECK_NOT_IN_PLACE_ROOT(4, recvbuf)
174 CHECK_TYPE(6, recvtype)
175 CHECK_NULL(5, MPI_ERR_COUNT, recvcounts)
176 CHECK_NULL(6, MPI_ERR_ARG, displs)
178 CHECK_NOT_IN_PLACE_ROOT(1, sendbuf)
182 CHECK_COLLECTIVE(comm, std::string(request == MPI_REQUEST_IGNORED ? "PMPI_Gatherv" : "PMPI_Igatherv") +
183 " with root " + std::to_string(root))
186 for (int i = 0; i < comm->size(); i++) {
187 CHECK_COUNT(5, recvcounts[i])
188 CHECK_BUFFER(4,recvbuf,recvcounts[i], recvtype)
192 const SmpiBenchGuard suspend_bench;
193 const void* real_sendbuf = sendbuf;
194 int real_sendcount = sendcount;
195 MPI_Datatype real_sendtype = sendtype;
196 if ((rank == root) && (sendbuf == MPI_IN_PLACE)) {
198 real_sendtype = recvtype;
201 aid_t pid = simgrid::s4u::this_actor::get_pid();
203 auto trace_recvcounts = std::make_shared<std::vector<int>>();
205 trace_recvcounts->insert(trace_recvcounts->end(), &recvcounts[0], &recvcounts[comm->size()]);
206 else //this is not significant outside of root, put 0 as we don't know if recvcounts is initialized
207 trace_recvcounts->insert(trace_recvcounts->end(), comm->size(), 0);
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,
213 nullptr, -1, 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);
226 int PMPI_Allgather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
227 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm){
228 return PMPI_Iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
231 int PMPI_Iallgather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
232 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
237 int rank = comm->rank();
238 CHECK_NOT_IN_PLACE(4, recvbuf)
239 if(sendbuf != MPI_IN_PLACE){
240 CHECK_COUNT(2, sendcount)
241 CHECK_TYPE(3, sendtype)
243 CHECK_TYPE(6, recvtype)
244 CHECK_COUNT(5, recvcount)
245 CHECK_BUFFER(1, sendbuf, sendcount, sendtype)
246 CHECK_BUFFER(4, recvbuf, recvcount, recvtype)
248 CHECK_COLLECTIVE(comm, request == MPI_REQUEST_IGNORED ? "PMPI_Allgather" : "PMPI_Iallggather")
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;
261 const SmpiBenchGuard suspend_bench;
263 aid_t 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 sendcount, recvcount,
269 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
270 if (request == MPI_REQUEST_IGNORED)
271 simgrid::smpi::colls::allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
273 simgrid::smpi::colls::iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, request);
275 TRACE_smpi_comm_out(pid);
279 int PMPI_Allgatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
280 void *recvbuf, const int *recvcounts, const int *displs, MPI_Datatype recvtype, MPI_Comm comm){
281 return PMPI_Iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, MPI_REQUEST_IGNORED);
284 int PMPI_Iallgatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
285 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
290 int rank = comm->rank();
291 if(sendbuf != MPI_IN_PLACE)
292 CHECK_TYPE(3, sendtype)
293 CHECK_TYPE(6, recvtype)
294 CHECK_NULL(5, MPI_ERR_COUNT, recvcounts)
295 CHECK_NULL(6, MPI_ERR_ARG, displs)
296 if(sendbuf != MPI_IN_PLACE){
297 CHECK_COUNT(2, sendcount)
298 CHECK_BUFFER(1, sendbuf, sendcount, sendtype)
301 CHECK_NOT_IN_PLACE(4, recvbuf)
302 for (int i = 0; i < comm->size(); i++) {
303 CHECK_COUNT(5, recvcounts[i])
304 CHECK_BUFFER(4, recvbuf, recvcounts[i], recvtype)
306 CHECK_COLLECTIVE(comm, MPI_REQUEST_IGNORED ? "PMPI_Allgatherv" : "PMPI_Iallgatherv")
308 const SmpiBenchGuard suspend_bench;
309 if (sendbuf == MPI_IN_PLACE) {
310 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * displs[comm->rank()];
311 sendcount = recvcounts[comm->rank()];
314 aid_t pid = simgrid::s4u::this_actor::get_pid();
316 auto trace_recvcounts = std::make_shared<std::vector<int>>();
317 trace_recvcounts->insert(trace_recvcounts->end(), &recvcounts[0], &recvcounts[comm->size()]);
320 pid, request == MPI_REQUEST_IGNORED ? "PMPI_Allgatherv" : "PMPI_Iallgatherv",
321 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "allgatherv" : "iallgatherv", -1,
323 -1, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
324 simgrid::smpi::Datatype::encode(recvtype)));
325 if (request == MPI_REQUEST_IGNORED)
326 simgrid::smpi::colls::allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
328 simgrid::smpi::colls::iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm,
331 TRACE_smpi_comm_out(pid);
335 int PMPI_Scatter(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
336 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
337 return PMPI_Iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
340 int PMPI_Iscatter(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
341 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
345 int rank = comm->rank();
348 CHECK_NOT_IN_PLACE_ROOT(1, sendbuf)
349 CHECK_COUNT(2, sendcount)
350 CHECK_TYPE(3, sendtype)
351 CHECK_BUFFER(1, sendbuf, sendcount, sendtype)
353 CHECK_NOT_IN_PLACE_ROOT(4, recvbuf)
355 if(recvbuf != MPI_IN_PLACE){
356 CHECK_COUNT(5, recvcount)
357 CHECK_TYPE(6, recvtype)
358 CHECK_BUFFER(4, recvbuf, recvcount, recvtype)
362 CHECK_COLLECTIVE(comm, std::string(request == MPI_REQUEST_IGNORED ? "PMPI_Scatter" : "PMPI_Iscatter") +
363 " with root " + std::to_string(root))
365 if (recvbuf == MPI_IN_PLACE) {
367 recvcount = sendcount;
370 if((rank == root) && (recvtype->size() * recvcount != sendtype->size() * sendcount)){
371 XBT_WARN("MPI_(I)Scatter : sent size to each process differs from receive size");
372 return MPI_ERR_TRUNCATE;
375 const SmpiBenchGuard suspend_bench;
377 aid_t pid = simgrid::s4u::this_actor::get_pid();
379 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Scatter" : "PMPI_Iscatter",
380 new simgrid::instr::CollTIData(
381 request == MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
382 sendcount, recvcount,
383 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
384 if (request == MPI_REQUEST_IGNORED)
385 simgrid::smpi::colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
387 simgrid::smpi::colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
389 TRACE_smpi_comm_out(pid);
393 int PMPI_Scatterv(const void *sendbuf, const int *sendcounts, const int *displs,
394 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
395 return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
398 int PMPI_Iscatterv(const void* sendbuf, const int* sendcounts, const int* displs, MPI_Datatype sendtype, void* recvbuf, int recvcount,
399 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
403 int rank = comm->rank();
404 if(recvbuf != MPI_IN_PLACE){
405 CHECK_COUNT(5, recvcount)
406 CHECK_TYPE(7, recvtype)
407 CHECK_BUFFER(4, recvbuf, recvcount, recvtype)
413 CHECK_NOT_IN_PLACE_ROOT(1, sendbuf)
414 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
415 CHECK_NULL(3, MPI_ERR_ARG, displs)
416 CHECK_TYPE(4, sendtype)
417 for (int i = 0; i < comm->size(); i++){
418 CHECK_COUNT(2, sendcounts[i])
419 CHECK_BUFFER(1, sendbuf, sendcounts[i], sendtype)
421 if (recvbuf == MPI_IN_PLACE) {
423 recvcount = sendcounts[rank];
426 CHECK_NOT_IN_PLACE_ROOT(4, recvbuf)
428 CHECK_COLLECTIVE(comm, std::string(request == MPI_REQUEST_IGNORED ? "PMPI_Scatterv" : "PMPI_Iscatterv") +
429 " with root " + std::to_string(root))
431 const SmpiBenchGuard suspend_bench;
433 aid_t pid = simgrid::s4u::this_actor::get_pid();
435 auto trace_sendcounts = std::make_shared<std::vector<int>>();
437 trace_sendcounts->insert(trace_sendcounts->end(), &sendcounts[0], &sendcounts[comm->size()]);
438 else //this is not significant outside of root, put 0 as we don't know if sendcounts is initialized
439 trace_sendcounts->insert(trace_sendcounts->end(), comm->size(), 0);
442 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Scatterv" : "PMPI_Iscatterv",
443 new simgrid::instr::VarCollTIData(
444 request == MPI_REQUEST_IGNORED ? "scatterv" : "iscatterv", root, -1,
445 trace_sendcounts, recvcount,
446 nullptr, simgrid::smpi::Datatype::encode(sendtype),
447 simgrid::smpi::Datatype::encode(recvtype)));
448 if (request == MPI_REQUEST_IGNORED)
449 simgrid::smpi::colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
451 simgrid::smpi::colls::iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm,
454 TRACE_smpi_comm_out(pid);
458 int PMPI_Reduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
460 return PMPI_Ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, MPI_REQUEST_IGNORED);
463 int PMPI_Ireduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm, MPI_Request* request)
467 int rank = comm->rank();
468 CHECK_TYPE(4, datatype)
469 CHECK_COUNT(3, count)
470 CHECK_BUFFER(1, sendbuf, count, datatype)
473 CHECK_NOT_IN_PLACE(2, recvbuf)
474 CHECK_BUFFER(5, recvbuf, count, datatype)
476 CHECK_OP(5, op, datatype)
479 CHECK_COLLECTIVE(comm, std::string(request == MPI_REQUEST_IGNORED ? "PMPI_Reduce" : "PMPI_Ireduce") + " with op " +
480 op->name() + " and root " + std::to_string(root))
482 const SmpiBenchGuard suspend_bench;
483 aid_t pid = simgrid::s4u::this_actor::get_pid();
485 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce" : "PMPI_Ireduce",
486 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "reduce" : "ireduce", root, 0,
488 simgrid::smpi::Datatype::encode(datatype), ""));
489 if (request == MPI_REQUEST_IGNORED)
490 simgrid::smpi::colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
492 simgrid::smpi::colls::ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, request);
494 TRACE_smpi_comm_out(pid);
498 int PMPI_Reduce_local(const void* inbuf, void* inoutbuf, int count, MPI_Datatype datatype, MPI_Op op)
502 CHECK_TYPE(4, datatype)
503 CHECK_COUNT(3, count)
504 CHECK_BUFFER(1, inbuf, count, datatype)
505 CHECK_BUFFER(2, inoutbuf, count, datatype)
506 CHECK_OP(5, op, datatype)
508 const SmpiBenchGuard suspend_bench;
509 op->apply(inbuf, inoutbuf, &count, datatype);
513 int PMPI_Allreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
515 return PMPI_Iallreduce(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
518 int PMPI_Iallreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
523 int rank = comm->rank();
524 CHECK_NOT_IN_PLACE(2, recvbuf)
525 CHECK_TYPE(4, datatype)
526 CHECK_OP(5, op, datatype)
527 CHECK_COUNT(3, count)
528 CHECK_BUFFER(1, sendbuf, count, datatype)
529 CHECK_BUFFER(2, recvbuf, count, datatype)
531 CHECK_COLLECTIVE(comm, std::string(request == MPI_REQUEST_IGNORED ? "PMPI_Alleduce" : "PMPI_Iallreduce") +
532 " with op " + op->name())
534 const SmpiBenchGuard suspend_bench;
535 std::vector<unsigned char> tmp_sendbuf;
536 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
538 aid_t pid = simgrid::s4u::this_actor::get_pid();
540 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Allreduce" : "PMPI_Iallreduce",
541 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "allreduce" : "iallreduce", -1, 0,
543 simgrid::smpi::Datatype::encode(datatype), ""));
545 if (request == MPI_REQUEST_IGNORED)
546 simgrid::smpi::colls::allreduce(real_sendbuf, recvbuf, count, datatype, op, comm);
548 simgrid::smpi::colls::iallreduce(real_sendbuf, recvbuf, count, datatype, op, comm, request);
550 TRACE_smpi_comm_out(pid);
554 int PMPI_Scan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
556 return PMPI_Iscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
559 int PMPI_Iscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request)
564 CHECK_TYPE(4, datatype)
565 CHECK_COUNT(3, count)
566 CHECK_BUFFER(1,sendbuf,count, datatype)
567 CHECK_BUFFER(2,recvbuf,count, datatype)
569 CHECK_OP(5, op, datatype)
570 CHECK_COLLECTIVE(comm,
571 std::string(request == MPI_REQUEST_IGNORED ? "PMPI_Scan" : "PMPI_Iscan") + " with op " + op->name())
573 const SmpiBenchGuard suspend_bench;
574 aid_t pid = simgrid::s4u::this_actor::get_pid();
575 std::vector<unsigned char> tmp_sendbuf;
576 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
578 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Scan" : "PMPI_Iscan",
579 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "scan" : "iscan", -1, 0.0,
580 count, 0, simgrid::smpi::Datatype::encode(datatype), ""));
583 if (request == MPI_REQUEST_IGNORED)
584 retval = simgrid::smpi::colls::scan(real_sendbuf, recvbuf, count, datatype, op, comm);
586 retval = simgrid::smpi::colls::iscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
588 TRACE_smpi_comm_out(pid);
592 int PMPI_Exscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
594 return PMPI_Iexscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
597 int PMPI_Iexscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request){
601 CHECK_TYPE(4, datatype)
602 CHECK_COUNT(3, count)
603 CHECK_BUFFER(1, sendbuf, count, datatype)
604 CHECK_BUFFER(2, recvbuf, count, datatype)
606 CHECK_OP(5, op, datatype)
607 CHECK_COLLECTIVE(comm, std::string(request == MPI_REQUEST_IGNORED ? "PMPI_Exscan" : "PMPI_Iexscan") + " with op " +
610 const SmpiBenchGuard suspend_bench;
611 aid_t pid = simgrid::s4u::this_actor::get_pid();
612 std::vector<unsigned char> tmp_sendbuf;
613 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
615 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Exscan" : "PMPI_Iexscan",
616 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "exscan" : "iexscan", -1, 0.0,
617 count, 0, simgrid::smpi::Datatype::encode(datatype), ""));
620 if (request == MPI_REQUEST_IGNORED)
621 retval = simgrid::smpi::colls::exscan(real_sendbuf, recvbuf, count, datatype, op, comm);
623 retval = simgrid::smpi::colls::iexscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
625 TRACE_smpi_comm_out(pid);
629 int PMPI_Reduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
631 return PMPI_Ireduce_scatter(sendbuf, recvbuf, recvcounts, datatype, op, comm, MPI_REQUEST_IGNORED);
634 int PMPI_Ireduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
639 int rank = comm->rank();
640 CHECK_NOT_IN_PLACE(2, recvbuf)
641 CHECK_TYPE(4, datatype)
642 CHECK_NULL(3, MPI_ERR_COUNT, recvcounts)
644 CHECK_OP(5, op, datatype)
645 for (int i = 0; i < comm->size(); i++) {
646 CHECK_COUNT(3, recvcounts[i])
647 CHECK_BUFFER(1, sendbuf, recvcounts[i], datatype)
648 CHECK_BUFFER(2, recvbuf, recvcounts[i], datatype)
650 CHECK_COLLECTIVE(comm, std::string(request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter" : "PMPI_Ireduce_scatter") +
651 " with op " + op->name())
653 const SmpiBenchGuard suspend_bench;
654 aid_t pid = simgrid::s4u::this_actor::get_pid();
655 auto trace_recvcounts = std::make_shared<std::vector<int>>();
656 trace_recvcounts->insert(trace_recvcounts->end(), &recvcounts[0], &recvcounts[comm->size()]);
660 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
661 totalcount += recvcounts[i];
663 std::vector<unsigned char> tmp_sendbuf;
664 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, totalcount, datatype);
666 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter" : "PMPI_Ireduce_scatter",
667 new simgrid::instr::VarCollTIData(
668 request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, -1, nullptr,
669 -1 , trace_recvcounts, std::to_string(0), simgrid::smpi::Datatype::encode(datatype)));
671 if (request == MPI_REQUEST_IGNORED)
672 simgrid::smpi::colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm);
674 simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm, request);
676 TRACE_smpi_comm_out(pid);
680 int PMPI_Reduce_scatter_block(const void *sendbuf, void *recvbuf, int recvcount,
681 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
683 return PMPI_Ireduce_scatter_block(sendbuf, recvbuf, recvcount, datatype, op, comm, MPI_REQUEST_IGNORED);
686 int PMPI_Ireduce_scatter_block(const void* sendbuf, void* recvbuf, int recvcount, MPI_Datatype datatype, MPI_Op op,
687 MPI_Comm comm, MPI_Request* request)
692 CHECK_TYPE(4, datatype)
693 CHECK_COUNT(3, recvcount)
694 CHECK_BUFFER(1, sendbuf, recvcount, datatype)
695 CHECK_BUFFER(2, recvbuf, recvcount, datatype)
697 CHECK_OP(5, op, datatype)
699 comm, std::string(request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter_block" : "PMPI_Ireduce_scatter_block") +
700 " with op " + op->name())
702 const SmpiBenchGuard suspend_bench;
703 int count = comm->size();
705 aid_t pid = simgrid::s4u::this_actor::get_pid();
706 auto trace_recvcounts = std::make_shared<std::vector<int>>(recvcount);
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, -1,
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);
728 int PMPI_Alltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
729 MPI_Datatype recvtype, MPI_Comm comm){
730 return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
733 int PMPI_Ialltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
734 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
739 if(sendbuf != MPI_IN_PLACE){
740 CHECK_TYPE(3, sendtype)
741 CHECK_COUNT(2, sendcount)
742 CHECK_BUFFER(1, sendbuf, sendcount, sendtype)
744 CHECK_TYPE(6, recvtype)
745 CHECK_COUNT(5, recvcount)
746 CHECK_COUNT(5, recvcount)
747 CHECK_BUFFER(4, recvbuf, recvcount, recvtype)
749 CHECK_COLLECTIVE(comm, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoall" : "PMPI_Ialltoall")
751 aid_t 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;
768 const SmpiBenchGuard suspend_bench;
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_sendcount, recvcount,
774 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
776 if (request == MPI_REQUEST_IGNORED)
778 simgrid::smpi::colls::alltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, comm);
780 retval = simgrid::smpi::colls::ialltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype,
783 TRACE_smpi_comm_out(pid);
787 int PMPI_Alltoallv(const void* sendbuf, const int* sendcounts, const int* senddispls, MPI_Datatype sendtype, void* recvbuf,
788 const int* recvcounts, const int* recvdispls, MPI_Datatype recvtype, MPI_Comm comm)
790 return PMPI_Ialltoallv(sendbuf, sendcounts, senddispls, sendtype, recvbuf, recvcounts, recvdispls, recvtype, comm, MPI_REQUEST_IGNORED);
793 int PMPI_Ialltoallv(const void* sendbuf, const int* sendcounts, const int* senddispls, MPI_Datatype sendtype, void* recvbuf,
794 const int* recvcounts, const int* recvdispls, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
799 if(sendbuf != MPI_IN_PLACE){
800 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
801 CHECK_NULL(3, MPI_ERR_ARG, senddispls)
802 CHECK_TYPE(4, sendtype)
804 CHECK_TYPE(8, recvtype)
805 CHECK_NULL(6, MPI_ERR_COUNT, recvcounts)
806 CHECK_NULL(7, MPI_ERR_ARG, recvdispls)
808 CHECK_COLLECTIVE(comm, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallv" : "PMPI_Ialltoallv")
810 aid_t pid = simgrid::s4u::this_actor::get_pid();
811 int size = comm->size();
812 for (int i = 0; i < size; i++) {
813 if(sendbuf != MPI_IN_PLACE){
814 CHECK_BUFFER(1, sendbuf, sendcounts[i], sendtype)
815 CHECK_COUNT(2, sendcounts[i])
817 CHECK_BUFFER(5, recvbuf, recvcounts[i], recvtype)
818 CHECK_COUNT(6, recvcounts[i])
821 const SmpiBenchGuard suspend_bench;
824 auto trace_sendcounts = std::make_shared<std::vector<int>>();
825 auto trace_recvcounts = std::make_shared<std::vector<int>>();
826 trace_recvcounts->insert(trace_recvcounts->end(), &recvcounts[0], &recvcounts[size]);
828 int dt_size_recv = recvtype->size();
830 const int* real_sendcounts = sendcounts;
831 const int* real_senddispls = senddispls;
832 MPI_Datatype real_sendtype = sendtype;
834 std::vector<unsigned char> tmp_sendbuf;
835 std::vector<int> tmp_sendcounts;
836 std::vector<int> tmp_senddispls;
837 const void* real_sendbuf;
839 if (sendbuf == MPI_IN_PLACE) {
840 tmp_sendcounts.assign(recvcounts, recvcounts + size);
841 real_sendcounts = tmp_sendcounts.data();
842 tmp_senddispls.assign(recvdispls, recvdispls + size);
843 real_senddispls = tmp_senddispls.data();
844 real_sendtype = recvtype;
847 for (int i = 0; i < size; i++) { // copy data to avoid bad free
848 send_size += real_sendcounts[i] ;
849 recv_size += recvcounts[i];
850 if (((recvdispls[i] + recvcounts[i]) * dt_size_recv) > maxsize)
851 maxsize = (recvdispls[i] + recvcounts[i]) * dt_size_recv;
853 real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
855 if(recvtype->size() * recvcounts[comm->rank()] != real_sendtype->size() * real_sendcounts[comm->rank()]){
856 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()]);
857 return MPI_ERR_TRUNCATE;
860 trace_sendcounts->insert(trace_sendcounts->end(), &real_sendcounts[0], &real_sendcounts[size]);
862 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallv" : "PMPI_Ialltoallv",
863 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
864 send_size, trace_sendcounts, recv_size, trace_recvcounts,
865 simgrid::smpi::Datatype::encode(real_sendtype),
866 simgrid::smpi::Datatype::encode(recvtype)));
869 if (request == MPI_REQUEST_IGNORED)
870 retval = simgrid::smpi::colls::alltoallv(real_sendbuf, real_sendcounts, real_senddispls, real_sendtype, recvbuf,
871 recvcounts, recvdispls, recvtype, comm);
873 retval = simgrid::smpi::colls::ialltoallv(real_sendbuf, real_sendcounts, real_senddispls, real_sendtype, recvbuf,
874 recvcounts, recvdispls, recvtype, comm, request);
876 TRACE_smpi_comm_out(pid);
880 int PMPI_Alltoallw(const void* sendbuf, const int* sendcounts, const int* senddispls, const MPI_Datatype* sendtypes, void* recvbuf,
881 const int* recvcounts, const int* recvdispls, const MPI_Datatype* recvtypes, MPI_Comm comm)
883 return PMPI_Ialltoallw(sendbuf, sendcounts, senddispls, sendtypes, recvbuf, recvcounts, recvdispls, recvtypes, comm, MPI_REQUEST_IGNORED);
886 int PMPI_Ialltoallw(const void* sendbuf, const int* sendcounts, const int* senddispls, const MPI_Datatype* sendtypes, void* recvbuf,
887 const int* recvcounts, const int* recvdispls, const MPI_Datatype* recvtypes, MPI_Comm comm, MPI_Request* request)
892 if(sendbuf != MPI_IN_PLACE){
893 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
894 CHECK_NULL(3, MPI_ERR_ARG, senddispls)
895 CHECK_NULL(4, MPI_ERR_TYPE, sendtypes)
897 CHECK_NULL(6, MPI_ERR_COUNT, recvcounts)
898 CHECK_NULL(7, MPI_ERR_ARG, recvdispls)
899 CHECK_NULL(8, MPI_ERR_TYPE, recvtypes)
901 aid_t pid = simgrid::s4u::this_actor::get_pid();
902 int size = comm->size();
903 for (int i = 0; i < size; i++) {
904 if(sendbuf != MPI_IN_PLACE){
905 CHECK_COUNT(2, sendcounts[i])
906 CHECK_TYPE(4, sendtypes[i])
907 CHECK_BUFFER(1, sendbuf, sendcounts[i], sendtypes[i])
909 CHECK_COUNT(6, recvcounts[i])
910 CHECK_TYPE(8, recvtypes[i])
911 CHECK_BUFFER(5, recvbuf, recvcounts[i], recvtypes[i])
913 CHECK_COLLECTIVE(comm, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallw" : "PMPI_Ialltoallw")
915 const SmpiBenchGuard suspend_bench;
919 auto trace_sendcounts = std::make_shared<std::vector<int>>();
920 auto trace_recvcounts = std::make_shared<std::vector<int>>();
921 trace_recvcounts->insert(trace_recvcounts->end(), &recvcounts[0], &recvcounts[size]);
923 const int* real_sendcounts = sendcounts;
924 const int* real_senddispls = senddispls;
925 const MPI_Datatype* real_sendtypes = sendtypes;
927 unsigned long maxsize = 0;
928 for (int i = 0; i < size; i++) { // copy data to avoid bad free
929 if (recvtypes[i] == MPI_DATATYPE_NULL)
931 recv_size += recvcounts[i] * recvtypes[i]->size();
932 if ((recvdispls[i] + (recvcounts[i] * recvtypes[i]->size())) > maxsize)
933 maxsize = recvdispls[i] + (recvcounts[i] * recvtypes[i]->size());
936 std::vector<unsigned char> tmp_sendbuf;
937 std::vector<int> tmp_sendcounts;
938 std::vector<int> tmp_senddispls;
939 std::vector<MPI_Datatype> tmp_sendtypes;
940 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
941 if (sendbuf == MPI_IN_PLACE) {
942 tmp_sendcounts.assign(recvcounts, recvcounts + size);
943 real_sendcounts = tmp_sendcounts.data();
944 tmp_senddispls.assign(recvdispls, recvdispls + size);
945 real_senddispls = tmp_senddispls.data();
946 tmp_sendtypes.assign(recvtypes, recvtypes + size);
947 real_sendtypes = tmp_sendtypes.data();
951 if(recvtypes[comm->rank()]->size() * recvcounts[comm->rank()] != real_sendtypes[comm->rank()]->size() * real_sendcounts[comm->rank()]){
952 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()]);
953 return MPI_ERR_TRUNCATE;
956 trace_sendcounts->insert(trace_sendcounts->end(), &real_sendcounts[0], &real_sendcounts[size]);
957 for (int i = 0; i < size; i++) { // copy data to avoid bad free
958 send_size += real_sendcounts[i] * real_sendtypes[i]->size();
961 TRACE_smpi_comm_in(pid, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallw" : "PMPI_Ialltoallw",
962 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
963 send_size, trace_sendcounts, recv_size, trace_recvcounts,
964 simgrid::smpi::Datatype::encode(real_sendtypes[0]),
965 simgrid::smpi::Datatype::encode(recvtypes[0])));
968 if (request == MPI_REQUEST_IGNORED)
969 retval = simgrid::smpi::colls::alltoallw(real_sendbuf, real_sendcounts, real_senddispls, real_sendtypes, recvbuf,
970 recvcounts, recvdispls, recvtypes, comm);
972 retval = simgrid::smpi::colls::ialltoallw(real_sendbuf, real_sendcounts, real_senddispls, real_sendtypes, recvbuf,
973 recvcounts, recvdispls, recvtypes, comm, request);
975 TRACE_smpi_comm_out(pid);