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,std::unique_ptr<unsigned char[]>& tmp_sendbuf, int count, MPI_Datatype datatype){
19 if (inplacebuf == MPI_IN_PLACE) {
20 tmp_sendbuf = std::make_unique<unsigned char[]>(count * datatype->get_extent());
21 simgrid::smpi::Datatype::copy(otherbuf, count, datatype, tmp_sendbuf.get(), count, datatype);
22 return tmp_sendbuf.get();
27 /* PMPI User level calls */
29 int PMPI_Barrier(MPI_Comm comm)
31 return PMPI_Ibarrier(comm, MPI_REQUEST_IGNORED);
34 int PMPI_Ibarrier(MPI_Comm comm, MPI_Request *request)
40 int rank = simgrid::s4u::this_actor::get_pid();
41 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Barrier" : "PMPI_Ibarrier",
42 new simgrid::instr::NoOpTIData(request == MPI_REQUEST_IGNORED ? "barrier" : "ibarrier"));
43 if (request == MPI_REQUEST_IGNORED) {
44 simgrid::smpi::colls::barrier(comm);
45 // Barrier can be used to synchronize RMA calls. Finish all requests from comm before.
46 comm->finish_rma_calls();
48 simgrid::smpi::colls::ibarrier(comm, request);
50 TRACE_smpi_comm_out(rank);
55 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
57 return PMPI_Ibcast(buf, count, datatype, root, comm, MPI_REQUEST_IGNORED);
60 int PMPI_Ibcast(void *buf, int count, MPI_Datatype datatype,
61 int root, MPI_Comm comm, MPI_Request* request)
64 CHECK_BUFFER(1, buf, count)
66 CHECK_TYPE(3, datatype)
71 int rank = simgrid::s4u::this_actor::get_pid();
72 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Bcast" : "PMPI_Ibcast",
73 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "bcast" : "ibcast", root, -1.0,
74 datatype->is_replayable() ? count : count * datatype->size(), -1,
75 simgrid::smpi::Datatype::encode(datatype), ""));
76 if (comm->size() > 1) {
77 if (request == MPI_REQUEST_IGNORED)
78 simgrid::smpi::colls::bcast(buf, count, datatype, root, comm);
80 simgrid::smpi::colls::ibcast(buf, count, datatype, root, comm, request);
82 if (request != MPI_REQUEST_IGNORED)
83 *request = MPI_REQUEST_NULL;
86 TRACE_smpi_comm_out(rank);
91 int PMPI_Gather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
92 int root, MPI_Comm comm){
93 return PMPI_Igather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
96 int PMPI_Igather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
97 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
100 if(sendbuf != MPI_IN_PLACE){
101 CHECK_BUFFER(1,sendbuf, sendcount)
102 CHECK_COUNT(2, sendcount)
103 CHECK_TYPE(3, sendtype)
105 if(comm->rank() == root){
106 CHECK_TYPE(6, recvtype)
107 CHECK_COUNT(5, recvcount)
108 CHECK_BUFFER(4, recvbuf, recvcount)
114 const void* real_sendbuf = sendbuf;
115 int real_sendcount = sendcount;
116 MPI_Datatype real_sendtype = sendtype;
117 if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
119 real_sendtype = recvtype;
121 int rank = simgrid::s4u::this_actor::get_pid();
123 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gather" : "PMPI_Igather",
124 new simgrid::instr::CollTIData(
125 request == MPI_REQUEST_IGNORED ? "gather" : "igather", root, -1.0,
126 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
127 (comm->rank() != root || recvtype->is_replayable()) ? recvcount : recvcount * recvtype->size(),
128 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
129 if (request == MPI_REQUEST_IGNORED)
130 simgrid::smpi::colls::gather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm);
132 simgrid::smpi::colls::igather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm,
135 TRACE_smpi_comm_out(rank);
140 int PMPI_Gatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, const int *recvcounts, const int *displs,
141 MPI_Datatype recvtype, int root, MPI_Comm comm){
142 return PMPI_Igatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, root, comm, MPI_REQUEST_IGNORED);
145 int PMPI_Igatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
146 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
149 CHECK_BUFFER(1, sendbuf, sendcount)
150 if(sendbuf != MPI_IN_PLACE){
151 CHECK_TYPE(3, sendtype)
152 CHECK_COUNT(2, sendcount)
154 if(comm->rank() == root){
155 CHECK_TYPE(6, recvtype)
156 CHECK_NULL(5, MPI_ERR_COUNT, recvcounts)
157 CHECK_NULL(6, MPI_ERR_ARG, displs)
162 if (comm->rank() == root){
163 for (int i = 0; i < comm->size(); i++) {
164 CHECK_COUNT(5, recvcounts[i])
165 CHECK_BUFFER(4,recvbuf,recvcounts[i])
170 const void* real_sendbuf = sendbuf;
171 int real_sendcount = sendcount;
172 MPI_Datatype real_sendtype = sendtype;
173 if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
175 real_sendtype = recvtype;
178 int rank = simgrid::s4u::this_actor::get_pid();
179 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
181 auto* trace_recvcounts = new std::vector<int>();
182 if (comm->rank() == root) {
183 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
184 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
187 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gatherv" : "PMPI_Igatherv",
188 new simgrid::instr::VarCollTIData(
189 request == MPI_REQUEST_IGNORED ? "gatherv" : "igatherv", root,
190 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
191 nullptr, dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(real_sendtype),
192 simgrid::smpi::Datatype::encode(recvtype)));
193 if (request == MPI_REQUEST_IGNORED)
194 simgrid::smpi::colls::gatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
197 simgrid::smpi::colls::igatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
198 root, comm, request);
200 TRACE_smpi_comm_out(rank);
205 int PMPI_Allgather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
206 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm){
207 return PMPI_Iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
210 int PMPI_Iallgather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
211 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
214 CHECK_BUFFER(1, sendbuf, sendcount)
215 CHECK_BUFFER(4, recvbuf, recvcount)
216 if(sendbuf != MPI_IN_PLACE){
217 CHECK_COUNT(2, sendcount)
218 CHECK_TYPE(3, sendtype)
220 CHECK_TYPE(6, recvtype)
221 CHECK_COUNT(5, recvcount)
225 if (sendbuf == MPI_IN_PLACE) {
226 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * recvcount * comm->rank();
227 sendcount = recvcount;
230 int rank = simgrid::s4u::this_actor::get_pid();
232 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgather" : "PMPI_Iallggather",
233 new simgrid::instr::CollTIData(
234 request == MPI_REQUEST_IGNORED ? "allgather" : "iallgather", -1, -1.0,
235 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
236 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
237 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
238 if (request == MPI_REQUEST_IGNORED)
239 simgrid::smpi::colls::allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
241 simgrid::smpi::colls::iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, request);
243 TRACE_smpi_comm_out(rank);
248 int PMPI_Allgatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
249 void *recvbuf, const int *recvcounts, const int *displs, MPI_Datatype recvtype, MPI_Comm comm){
250 return PMPI_Iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, MPI_REQUEST_IGNORED);
253 int PMPI_Iallgatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
254 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
257 CHECK_BUFFER(1, sendbuf, sendcount)
258 if(sendbuf != MPI_IN_PLACE)
259 CHECK_TYPE(3, sendtype)
260 CHECK_TYPE(6, recvtype)
261 CHECK_NULL(5, MPI_ERR_COUNT, recvcounts)
262 CHECK_NULL(6, MPI_ERR_ARG, displs)
263 if(sendbuf != MPI_IN_PLACE)
264 CHECK_COUNT(2, sendcount)
266 for (int i = 0; i < comm->size(); i++) {
267 CHECK_COUNT(5, recvcounts[i])
268 CHECK_BUFFER(4, recvbuf, recvcounts[i])
272 if (sendbuf == MPI_IN_PLACE) {
273 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * displs[comm->rank()];
274 sendcount = recvcounts[comm->rank()];
277 int rank = simgrid::s4u::this_actor::get_pid();
278 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
280 auto* trace_recvcounts = new std::vector<int>();
281 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
282 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
286 rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgatherv" : "PMPI_Iallgatherv",
287 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "allgatherv" : "iallgatherv", -1,
288 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(), nullptr,
289 dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
290 simgrid::smpi::Datatype::encode(recvtype)));
291 if (request == MPI_REQUEST_IGNORED)
292 simgrid::smpi::colls::allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
294 simgrid::smpi::colls::iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm,
297 TRACE_smpi_comm_out(rank);
302 int PMPI_Scatter(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
303 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
304 return PMPI_Iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
307 int PMPI_Iscatter(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
308 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
311 if(comm->rank() == root){
312 CHECK_BUFFER(1, sendbuf, sendcount)
313 CHECK_COUNT(2, sendcount)
314 CHECK_TYPE(3, sendtype)
316 if(recvbuf != MPI_IN_PLACE){
317 CHECK_BUFFER(4, recvbuf, recvcount)
318 CHECK_COUNT(5, recvcount)
319 CHECK_TYPE(6, recvtype)
325 if (recvbuf == MPI_IN_PLACE) {
327 recvcount = sendcount;
329 int rank = simgrid::s4u::this_actor::get_pid();
331 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatter" : "PMPI_Iscatter",
332 new simgrid::instr::CollTIData(
333 request == MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
334 (comm->rank() != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
335 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
336 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
337 if (request == MPI_REQUEST_IGNORED)
338 simgrid::smpi::colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
340 simgrid::smpi::colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
342 TRACE_smpi_comm_out(rank);
347 int PMPI_Scatterv(const void *sendbuf, const int *sendcounts, const int *displs,
348 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
349 return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
352 int PMPI_Iscatterv(const void* sendbuf, const int* sendcounts, const int* displs, MPI_Datatype sendtype, void* recvbuf, int recvcount,
353 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
356 if(recvbuf != MPI_IN_PLACE){
357 CHECK_BUFFER(4, recvbuf, recvcount)
358 CHECK_COUNT(5, recvcount)
359 CHECK_TYPE(7, recvtype)
363 if (comm->rank() == root) {
364 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
365 CHECK_NULL(3, MPI_ERR_ARG, displs)
366 CHECK_TYPE(4, sendtype)
367 for (int i = 0; i < comm->size(); i++){
368 CHECK_BUFFER(1, sendbuf, sendcounts[i])
369 CHECK_COUNT(2, sendcounts[i])
371 if (recvbuf == MPI_IN_PLACE) {
373 recvcount = sendcounts[comm->rank()];
379 int rank = simgrid::s4u::this_actor::get_pid();
380 int dt_size_send = sendtype->is_replayable() ? 1 : sendtype->size();
382 auto* trace_sendcounts = new std::vector<int>();
383 if (comm->rank() == root) {
384 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
385 trace_sendcounts->push_back(sendcounts[i] * dt_size_send);
389 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatterv" : "PMPI_Iscatterv",
390 new simgrid::instr::VarCollTIData(
391 request == MPI_REQUEST_IGNORED ? "scatterv" : "iscatterv", root, dt_size_send,
392 trace_sendcounts, recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
393 nullptr, simgrid::smpi::Datatype::encode(sendtype),
394 simgrid::smpi::Datatype::encode(recvtype)));
395 if (request == MPI_REQUEST_IGNORED)
396 simgrid::smpi::colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
398 simgrid::smpi::colls::iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm,
401 TRACE_smpi_comm_out(rank);
406 int PMPI_Reduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
408 return PMPI_Ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, MPI_REQUEST_IGNORED);
411 int PMPI_Ireduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm, MPI_Request* request)
414 CHECK_BUFFER(1, sendbuf, count)
415 if(comm->rank() == root)
416 CHECK_BUFFER(5, recvbuf, count)
417 CHECK_TYPE(4, datatype)
418 CHECK_COUNT(3, count)
424 int rank = simgrid::s4u::this_actor::get_pid();
426 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce" : "PMPI_Ireduce",
427 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "reduce" : "ireduce", root, 0,
428 datatype->is_replayable() ? count : count * datatype->size(), -1,
429 simgrid::smpi::Datatype::encode(datatype), ""));
430 if (request == MPI_REQUEST_IGNORED)
431 simgrid::smpi::colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
433 simgrid::smpi::colls::ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, request);
435 TRACE_smpi_comm_out(rank);
440 int PMPI_Reduce_local(const void* inbuf, void* inoutbuf, int count, MPI_Datatype datatype, MPI_Op op)
442 CHECK_BUFFER(1, inbuf, count)
443 CHECK_BUFFER(2, inoutbuf, count)
444 CHECK_TYPE(4, datatype)
445 CHECK_COUNT(3, count)
449 op->apply(inbuf, inoutbuf, &count, datatype);
454 int PMPI_Allreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
456 return PMPI_Iallreduce(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
459 int PMPI_Iallreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
462 CHECK_BUFFER(1, sendbuf, count)
463 CHECK_BUFFER(2, recvbuf, count)
464 CHECK_TYPE(4, datatype)
465 CHECK_COUNT(3, count)
470 std::unique_ptr<unsigned char[]> tmp_sendbuf;
471 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
473 int rank = simgrid::s4u::this_actor::get_pid();
475 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allreduce" : "PMPI_Iallreduce",
476 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "allreduce" : "iallreduce", -1, 0,
477 datatype->is_replayable() ? count : count * datatype->size(), -1,
478 simgrid::smpi::Datatype::encode(datatype), ""));
480 if (request == MPI_REQUEST_IGNORED)
481 simgrid::smpi::colls::allreduce(real_sendbuf, recvbuf, count, datatype, op, comm);
483 simgrid::smpi::colls::iallreduce(real_sendbuf, recvbuf, count, datatype, op, comm, request);
485 TRACE_smpi_comm_out(rank);
490 int PMPI_Scan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
492 return PMPI_Iscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
495 int PMPI_Iscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request)
498 CHECK_BUFFER(1,sendbuf,count)
499 CHECK_BUFFER(2,recvbuf,count)
500 CHECK_TYPE(4, datatype)
501 CHECK_COUNT(3, count)
506 int rank = simgrid::s4u::this_actor::get_pid();
507 std::unique_ptr<unsigned char[]> tmp_sendbuf;
508 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
510 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scan" : "PMPI_Iscan",
511 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "scan" : "iscan", -1,
512 datatype->is_replayable() ? count : count * datatype->size(),
513 simgrid::smpi::Datatype::encode(datatype)));
516 if (request == MPI_REQUEST_IGNORED)
517 retval = simgrid::smpi::colls::scan(real_sendbuf, recvbuf, count, datatype, op, comm);
519 retval = simgrid::smpi::colls::iscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
521 TRACE_smpi_comm_out(rank);
526 int PMPI_Exscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
528 return PMPI_Iexscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
531 int PMPI_Iexscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request){
533 CHECK_BUFFER(1, sendbuf, count)
534 CHECK_BUFFER(2, recvbuf, count)
535 CHECK_TYPE(4, datatype)
536 CHECK_COUNT(3, count)
541 int rank = simgrid::s4u::this_actor::get_pid();
542 std::unique_ptr<unsigned char[]> tmp_sendbuf;
543 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
545 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Exscan" : "PMPI_Iexscan",
546 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "exscan" : "iexscan", -1,
547 datatype->is_replayable() ? count : count * datatype->size(),
548 simgrid::smpi::Datatype::encode(datatype)));
551 if (request == MPI_REQUEST_IGNORED)
552 retval = simgrid::smpi::colls::exscan(real_sendbuf, recvbuf, count, datatype, op, comm);
554 retval = simgrid::smpi::colls::iexscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
556 TRACE_smpi_comm_out(rank);
561 int PMPI_Reduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
563 return PMPI_Ireduce_scatter(sendbuf, recvbuf, recvcounts, datatype, op, comm, MPI_REQUEST_IGNORED);
566 int PMPI_Ireduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
569 CHECK_TYPE(4, datatype)
570 CHECK_NULL(3, MPI_ERR_COUNT, recvcounts)
573 for (int i = 0; i < comm->size(); i++) {
574 CHECK_COUNT(3, recvcounts[i])
575 CHECK_BUFFER(1, sendbuf, recvcounts[i])
576 CHECK_BUFFER(2, recvbuf, recvcounts[i])
580 int rank = simgrid::s4u::this_actor::get_pid();
581 auto* trace_recvcounts = new std::vector<int>();
582 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
585 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
586 trace_recvcounts->push_back(recvcounts[i] * dt_send_size);
587 totalcount += recvcounts[i];
589 std::unique_ptr<unsigned char[]> tmp_sendbuf;
590 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, totalcount, datatype);
592 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter" : "PMPI_Ireduce_scatter",
593 new simgrid::instr::VarCollTIData(
594 request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, dt_send_size, nullptr,
595 -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
597 if (request == MPI_REQUEST_IGNORED)
598 simgrid::smpi::colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm);
600 simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm, request);
602 TRACE_smpi_comm_out(rank);
607 int PMPI_Reduce_scatter_block(const void *sendbuf, void *recvbuf, int recvcount,
608 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
610 return PMPI_Ireduce_scatter_block(sendbuf, recvbuf, recvcount, datatype, op, comm, MPI_REQUEST_IGNORED);
613 int PMPI_Ireduce_scatter_block(const void* sendbuf, void* recvbuf, int recvcount, MPI_Datatype datatype, MPI_Op op,
614 MPI_Comm comm, MPI_Request* request)
617 CHECK_BUFFER(1, sendbuf, recvcount)
618 CHECK_BUFFER(2, recvbuf, recvcount)
619 CHECK_TYPE(4, datatype)
620 CHECK_COUNT(3, recvcount)
625 int count = comm->size();
627 int rank = simgrid::s4u::this_actor::get_pid();
628 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
629 auto* trace_recvcounts = new std::vector<int>(recvcount * dt_send_size); // copy data to avoid bad free
630 std::unique_ptr<unsigned char[]> tmp_sendbuf;
631 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, recvcount * count, datatype);
634 rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter_block" : "PMPI_Ireduce_scatter_block",
635 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, 0,
636 nullptr, -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
638 auto* recvcounts = new int[count];
639 for (int i = 0; i < count; i++)
640 recvcounts[i] = recvcount;
641 if (request == MPI_REQUEST_IGNORED)
642 simgrid::smpi::colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm);
644 simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm, request);
647 TRACE_smpi_comm_out(rank);
652 int PMPI_Alltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
653 MPI_Datatype recvtype, MPI_Comm comm){
654 return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
657 int PMPI_Ialltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
658 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
661 CHECK_BUFFER(1, sendbuf, sendcount)
662 CHECK_BUFFER(4, recvbuf, recvcount)
663 if(sendbuf != MPI_IN_PLACE)
664 CHECK_TYPE(3, sendtype)
665 CHECK_TYPE(6, recvtype)
666 CHECK_COUNT(5, recvcount)
667 if(sendbuf != MPI_IN_PLACE)
668 CHECK_COUNT(2, sendcount)
669 CHECK_COUNT(5, recvcount)
673 int rank = simgrid::s4u::this_actor::get_pid();
674 int real_sendcount = sendcount;
675 MPI_Datatype real_sendtype = sendtype;
677 std::unique_ptr<unsigned char[]> tmp_sendbuf;
678 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, recvcount * comm->size(), recvtype);
680 if (sendbuf == MPI_IN_PLACE) {
681 real_sendcount = recvcount;
682 real_sendtype = recvtype;
685 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoall" : "PMPI_Ialltoall",
686 new simgrid::instr::CollTIData(
687 request == MPI_REQUEST_IGNORED ? "alltoall" : "ialltoall", -1, -1.0,
688 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
689 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
690 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
692 if (request == MPI_REQUEST_IGNORED)
694 simgrid::smpi::colls::alltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, comm);
696 retval = simgrid::smpi::colls::ialltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype,
699 TRACE_smpi_comm_out(rank);
704 int PMPI_Alltoallv(const void* sendbuf, const int* sendcounts, const int* senddispls, MPI_Datatype sendtype, void* recvbuf,
705 const int* recvcounts, const int* recvdispls, MPI_Datatype recvtype, MPI_Comm comm)
707 return PMPI_Ialltoallv(sendbuf, sendcounts, senddispls, sendtype, recvbuf, recvcounts, recvdispls, recvtype, comm, MPI_REQUEST_IGNORED);
710 int PMPI_Ialltoallv(const void* sendbuf, const int* sendcounts, const int* senddispls, MPI_Datatype sendtype, void* recvbuf,
711 const int* recvcounts, const int* recvdispls, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
714 if(sendbuf != MPI_IN_PLACE){
715 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
716 CHECK_NULL(3, MPI_ERR_ARG, senddispls)
717 CHECK_TYPE(4, sendtype)
719 CHECK_TYPE(8, recvtype)
720 CHECK_NULL(6, MPI_ERR_COUNT, recvcounts)
721 CHECK_NULL(7, MPI_ERR_ARG, recvdispls)
724 int rank = simgrid::s4u::this_actor::get_pid();
725 int size = comm->size();
726 for (int i = 0; i < size; i++) {
727 if(sendbuf != MPI_IN_PLACE){
728 CHECK_BUFFER(1, sendbuf, sendcounts[i])
729 CHECK_COUNT(2, sendcounts[i])
731 CHECK_BUFFER(5, recvbuf, recvcounts[i])
732 CHECK_COUNT(6, recvcounts[i])
738 auto* trace_sendcounts = new std::vector<int>();
739 auto* trace_recvcounts = new std::vector<int>();
740 int dt_size_recv = recvtype->size();
742 const int* real_sendcounts = sendcounts;
743 const int* real_senddispls = senddispls;
744 MPI_Datatype real_sendtype = sendtype;
746 for (int i = 0; i < size; i++) { // copy data to avoid bad free
747 recv_size += recvcounts[i] * dt_size_recv;
748 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
749 if (((recvdispls[i] + recvcounts[i]) * dt_size_recv) > maxsize)
750 maxsize = (recvdispls[i] + recvcounts[i]) * dt_size_recv;
753 std::unique_ptr<unsigned char[]> tmp_sendbuf;
754 std::unique_ptr<int[]> tmp_sendcounts;
755 std::unique_ptr<int[]> tmp_senddispls;
756 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
757 if (sendbuf == MPI_IN_PLACE) {
758 tmp_sendcounts = std::make_unique<int[]>(size);
759 std::copy(recvcounts, recvcounts + size, tmp_sendcounts.get());
760 real_sendcounts = tmp_sendcounts.get();
761 tmp_senddispls = std::make_unique<int[]>(size);
762 std::copy(recvdispls, recvdispls + size, tmp_senddispls.get());
763 real_senddispls = tmp_senddispls.get();
764 real_sendtype = recvtype;
767 int dt_size_send = real_sendtype->size();
769 for (int i = 0; i < size; i++) { // copy data to avoid bad free
770 send_size += real_sendcounts[i] * dt_size_send;
771 trace_sendcounts->push_back(real_sendcounts[i] * dt_size_send);
774 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallv" : "PMPI_Ialltoallv",
775 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
776 send_size, trace_sendcounts, recv_size, trace_recvcounts,
777 simgrid::smpi::Datatype::encode(real_sendtype),
778 simgrid::smpi::Datatype::encode(recvtype)));
781 if (request == MPI_REQUEST_IGNORED)
782 retval = simgrid::smpi::colls::alltoallv(real_sendbuf, real_sendcounts, real_senddispls, real_sendtype, recvbuf,
783 recvcounts, recvdispls, recvtype, comm);
785 retval = simgrid::smpi::colls::ialltoallv(real_sendbuf, real_sendcounts, real_senddispls, real_sendtype, recvbuf,
786 recvcounts, recvdispls, recvtype, comm, request);
788 TRACE_smpi_comm_out(rank);
793 int PMPI_Alltoallw(const void* sendbuf, const int* sendcounts, const int* senddispls, const MPI_Datatype* sendtypes, void* recvbuf,
794 const int* recvcounts, const int* recvdispls, const MPI_Datatype* recvtypes, MPI_Comm comm)
796 return PMPI_Ialltoallw(sendbuf, sendcounts, senddispls, sendtypes, recvbuf, recvcounts, recvdispls, recvtypes, comm, MPI_REQUEST_IGNORED);
799 int PMPI_Ialltoallw(const void* sendbuf, const int* sendcounts, const int* senddispls, const MPI_Datatype* sendtypes, void* recvbuf,
800 const int* recvcounts, const int* recvdispls, const MPI_Datatype* recvtypes, MPI_Comm comm, MPI_Request* request)
803 if(sendbuf != MPI_IN_PLACE){
804 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
805 CHECK_NULL(3, MPI_ERR_ARG, senddispls)
806 CHECK_NULL(4, MPI_ERR_TYPE, sendtypes)
808 CHECK_NULL(6, MPI_ERR_COUNT, recvcounts)
809 CHECK_NULL(7, MPI_ERR_ARG, recvdispls)
810 CHECK_NULL(8, MPI_ERR_TYPE, recvtypes)
812 int rank = simgrid::s4u::this_actor::get_pid();
813 int size = comm->size();
814 for (int i = 0; i < size; i++) {
815 if(sendbuf != MPI_IN_PLACE){
816 CHECK_BUFFER(1, sendbuf, sendcounts[i])
817 CHECK_COUNT(2, sendcounts[i])
818 CHECK_TYPE(4, sendtypes[i])
820 CHECK_BUFFER(5, recvbuf, recvcounts[i])
821 CHECK_COUNT(6, recvcounts[i])
822 CHECK_TYPE(8, recvtypes[i])
829 auto* trace_sendcounts = new std::vector<int>();
830 auto* trace_recvcounts = new std::vector<int>();
832 const int* real_sendcounts = sendcounts;
833 const int* real_senddispls = senddispls;
834 const MPI_Datatype* real_sendtypes = sendtypes;
835 unsigned long maxsize = 0;
836 for (int i = 0; i < size; i++) { // copy data to avoid bad free
837 if (recvtypes[i] == MPI_DATATYPE_NULL) {
838 delete trace_recvcounts;
839 delete trace_sendcounts;
842 recv_size += recvcounts[i] * recvtypes[i]->size();
843 trace_recvcounts->push_back(recvcounts[i] * recvtypes[i]->size());
844 if ((recvdispls[i] + (recvcounts[i] * recvtypes[i]->size())) > maxsize)
845 maxsize = recvdispls[i] + (recvcounts[i] * recvtypes[i]->size());
848 std::unique_ptr<unsigned char[]> tmp_sendbuf;
849 std::unique_ptr<int[]> tmp_sendcounts;
850 std::unique_ptr<int[]> tmp_senddispls;
851 std::unique_ptr<MPI_Datatype[]> tmp_sendtypes;
852 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
853 if (sendbuf == MPI_IN_PLACE) {
854 tmp_sendcounts = std::make_unique<int[]>(size);
855 std::copy(recvcounts, recvcounts + size, tmp_sendcounts.get());
856 real_sendcounts = tmp_sendcounts.get();
857 tmp_senddispls = std::make_unique<int[]>(size);
858 std::copy(recvdispls, recvdispls + size, tmp_senddispls.get());
859 real_senddispls = tmp_senddispls.get();
860 tmp_sendtypes = std::make_unique<MPI_Datatype[]>(size);
861 std::copy(recvtypes, recvtypes + size, tmp_sendtypes.get());
862 real_sendtypes = tmp_sendtypes.get();
865 for (int i = 0; i < size; i++) { // copy data to avoid bad free
866 send_size += real_sendcounts[i] * real_sendtypes[i]->size();
867 trace_sendcounts->push_back(real_sendcounts[i] * real_sendtypes[i]->size());
870 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallw" : "PMPI_Ialltoallw",
871 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
872 send_size, trace_sendcounts, recv_size, trace_recvcounts,
873 simgrid::smpi::Datatype::encode(real_sendtypes[0]),
874 simgrid::smpi::Datatype::encode(recvtypes[0])));
877 if (request == MPI_REQUEST_IGNORED)
878 retval = simgrid::smpi::colls::alltoallw(real_sendbuf, real_sendcounts, real_senddispls, real_sendtypes, recvbuf,
879 recvcounts, recvdispls, recvtypes, comm);
881 retval = simgrid::smpi::colls::ialltoallw(real_sendbuf, real_sendcounts, real_senddispls, real_sendtypes, recvbuf,
882 recvcounts, recvdispls, recvtypes, comm, request);
884 TRACE_smpi_comm_out(rank);