1 /* Copyright (c) 2007-2019. The SimGrid Team. All rights reserved. */
3 /* This program is free software; you can redistribute it and/or modify it
4 * under the terms of the license (GNU LGPL) which comes with this package. */
7 #include "smpi_coll.hpp"
8 #include "smpi_comm.hpp"
9 #include "smpi_request.hpp"
10 #include "smpi_datatype_derived.hpp"
11 #include "smpi_op.hpp"
12 #include "src/smpi/include/smpi_actor.hpp"
14 XBT_LOG_EXTERNAL_DEFAULT_CATEGORY(smpi_pmpi);
16 /* PMPI User level calls */
18 int PMPI_Barrier(MPI_Comm comm)
20 return PMPI_Ibarrier(comm, MPI_REQUEST_IGNORED);
23 int PMPI_Ibarrier(MPI_Comm comm, MPI_Request *request)
25 if (comm == MPI_COMM_NULL)
27 if (request == nullptr)
31 int rank = simgrid::s4u::this_actor::get_pid();
32 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Barrier" : "PMPI_Ibarrier",
33 new simgrid::instr::NoOpTIData(request == MPI_REQUEST_IGNORED ? "barrier" : "ibarrier"));
34 if (request == MPI_REQUEST_IGNORED) {
35 simgrid::smpi::Colls::barrier(comm);
36 // Barrier can be used to synchronize RMA calls. Finish all requests from comm before.
37 comm->finish_rma_calls();
39 simgrid::smpi::Colls::ibarrier(comm, request);
41 TRACE_smpi_comm_out(rank);
46 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
48 return PMPI_Ibcast(buf, count, datatype, root, comm, MPI_REQUEST_IGNORED);
51 int PMPI_Ibcast(void *buf, int count, MPI_Datatype datatype,
52 int root, MPI_Comm comm, MPI_Request* request)
54 if (comm == MPI_COMM_NULL)
56 if (buf == nullptr && count > 0)
57 return MPI_ERR_BUFFER;
58 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
62 if (root < 0 || root >= comm->size())
64 if (request == nullptr)
68 int rank = simgrid::s4u::this_actor::get_pid();
69 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Bcast" : "PMPI_Ibcast",
70 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "bcast" : "ibcast", root, -1.0,
71 datatype->is_replayable() ? count : count * datatype->size(), -1,
72 simgrid::smpi::Datatype::encode(datatype), ""));
73 if (comm->size() > 1) {
74 if (request == MPI_REQUEST_IGNORED)
75 simgrid::smpi::Colls::bcast(buf, count, datatype, root, comm);
77 simgrid::smpi::Colls::ibcast(buf, count, datatype, root, comm, request);
79 if (request != MPI_REQUEST_IGNORED)
80 *request = MPI_REQUEST_NULL;
83 TRACE_smpi_comm_out(rank);
88 int PMPI_Gather(void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
89 int root, MPI_Comm comm){
90 return PMPI_Igather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
93 int PMPI_Igather(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
94 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
96 if (comm == MPI_COMM_NULL)
98 if ((sendbuf == nullptr) || ((comm->rank() == root) && recvbuf == nullptr))
99 return MPI_ERR_BUFFER;
100 if (((sendbuf != MPI_IN_PLACE && sendcount > 0) && (sendtype == MPI_DATATYPE_NULL)) ||
101 ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL)))
103 if (((sendbuf != MPI_IN_PLACE) && (sendcount < 0)) || ((comm->rank() == root) && (recvcount < 0)))
104 return MPI_ERR_COUNT;
105 if (root < 0 || root >= comm->size())
107 if (request == nullptr)
111 char* sendtmpbuf = static_cast<char*>(sendbuf);
112 int sendtmpcount = sendcount;
113 MPI_Datatype sendtmptype = sendtype;
114 if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
116 sendtmptype = recvtype;
118 int rank = simgrid::s4u::this_actor::get_pid();
120 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gather" : "PMPI_Igather",
121 new simgrid::instr::CollTIData(
122 request == MPI_REQUEST_IGNORED ? "gather" : "igather", root, -1.0,
123 sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(),
124 (comm->rank() != root || recvtype->is_replayable()) ? recvcount : recvcount * recvtype->size(),
125 simgrid::smpi::Datatype::encode(sendtmptype), simgrid::smpi::Datatype::encode(recvtype)));
126 if (request == MPI_REQUEST_IGNORED)
127 simgrid::smpi::Colls::gather(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, root, comm);
129 simgrid::smpi::Colls::igather(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, root, comm,
132 TRACE_smpi_comm_out(rank);
137 int PMPI_Gatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, int *recvcounts, int *displs,
138 MPI_Datatype recvtype, int root, MPI_Comm comm){
139 return PMPI_Igatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, root, comm, MPI_REQUEST_IGNORED);
142 int PMPI_Igatherv(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int* recvcounts, int* displs,
143 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
145 if (comm == MPI_COMM_NULL)
147 if ((sendbuf == nullptr && sendcount > 0) || ((comm->rank() == root) && recvbuf == nullptr))
148 return MPI_ERR_BUFFER;
149 if (((sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
150 ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL)))
152 if ((sendbuf != MPI_IN_PLACE) && (sendcount < 0))
153 return MPI_ERR_COUNT;
154 if ((comm->rank() == root) && (recvcounts == nullptr || displs == nullptr))
156 if (root < 0 || root >= comm->size())
158 if (request == nullptr)
161 for (int i = 0; i < comm->size(); i++) {
162 if ((comm->rank() == root) && (recvcounts[i] < 0))
163 return MPI_ERR_COUNT;
167 char* sendtmpbuf = static_cast<char*>(sendbuf);
168 int sendtmpcount = sendcount;
169 MPI_Datatype sendtmptype = sendtype;
170 if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
172 sendtmptype = recvtype;
175 int rank = simgrid::s4u::this_actor::get_pid();
176 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
178 std::vector<int>* trace_recvcounts = new std::vector<int>;
179 if (comm->rank() == root) {
180 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
181 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
184 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gatherv" : "PMPI_Igatherv",
185 new simgrid::instr::VarCollTIData(
186 request == MPI_REQUEST_IGNORED ? "gatherv" : "igatherv", root,
187 sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(), nullptr,
188 dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtmptype),
189 simgrid::smpi::Datatype::encode(recvtype)));
190 if (request == MPI_REQUEST_IGNORED)
191 simgrid::smpi::Colls::gatherv(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcounts, displs, recvtype, root,
194 simgrid::smpi::Colls::igatherv(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcounts, displs, recvtype, root,
197 TRACE_smpi_comm_out(rank);
202 int PMPI_Allgather(void *sendbuf, int sendcount, MPI_Datatype sendtype,
203 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm){
204 return PMPI_Iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
207 int PMPI_Iallgather(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
208 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
210 if (comm == MPI_COMM_NULL)
212 if ((sendbuf == nullptr && sendcount > 0) || (recvbuf == nullptr))
213 return MPI_ERR_BUFFER;
214 if (((sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) || (recvtype == MPI_DATATYPE_NULL))
216 if (((sendbuf != MPI_IN_PLACE) && (sendcount < 0)) || (recvcount < 0))
217 return MPI_ERR_COUNT;
218 if (request == nullptr)
222 if (sendbuf == MPI_IN_PLACE) {
223 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * recvcount * comm->rank();
224 sendcount = recvcount;
227 int rank = simgrid::s4u::this_actor::get_pid();
229 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgather" : "PMPI_Iallggather",
230 new simgrid::instr::CollTIData(
231 request == MPI_REQUEST_IGNORED ? "allgather" : "iallgather", -1, -1.0,
232 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
233 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
234 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
235 if (request == MPI_REQUEST_IGNORED)
236 simgrid::smpi::Colls::allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
238 simgrid::smpi::Colls::iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, request);
240 TRACE_smpi_comm_out(rank);
245 int PMPI_Allgatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype,
246 void *recvbuf, int *recvcounts, int *displs, MPI_Datatype recvtype, MPI_Comm comm){
247 return PMPI_Iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, MPI_REQUEST_IGNORED);
250 int PMPI_Iallgatherv(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int* recvcounts, int* displs,
251 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
253 if (comm == MPI_COMM_NULL)
255 if ((sendbuf == nullptr && sendcount > 0) || (recvbuf == nullptr))
256 return MPI_ERR_BUFFER;
257 if (((sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) || (recvtype == MPI_DATATYPE_NULL))
259 if ((sendbuf != MPI_IN_PLACE) && (sendcount < 0))
260 return MPI_ERR_COUNT;
261 if (recvcounts == nullptr || displs == nullptr)
263 if (request == nullptr)
266 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
267 if (recvcounts[i] < 0)
268 return MPI_ERR_COUNT;
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 std::vector<int>* 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(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(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
308 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
310 if (comm == MPI_COMM_NULL)
312 if (((comm->rank() == root) && (sendtype == MPI_DATATYPE_NULL || not sendtype->is_valid())) ||
313 ((recvbuf != MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL || not recvtype->is_valid())))
315 if (((comm->rank() == root) && (sendcount < 0)) || ((recvbuf != MPI_IN_PLACE) && (recvcount < 0)))
316 return MPI_ERR_COUNT;
317 if ((sendbuf == recvbuf) || ((comm->rank() == root) && sendcount > 0 && (sendbuf == nullptr)) ||
318 (recvcount > 0 && recvbuf == nullptr))
319 return MPI_ERR_BUFFER;
320 if (root < 0 || root >= comm->size())
322 if (request == nullptr)
326 if (recvbuf == MPI_IN_PLACE) {
328 recvcount = sendcount;
330 int rank = simgrid::s4u::this_actor::get_pid();
332 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatter" : "PMPI_Iscatter",
333 new simgrid::instr::CollTIData(
334 request == MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
335 (comm->rank() != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
336 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
337 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
338 if (request == MPI_REQUEST_IGNORED)
339 simgrid::smpi::Colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
341 simgrid::smpi::Colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
343 TRACE_smpi_comm_out(rank);
348 int PMPI_Scatterv(void *sendbuf, int *sendcounts, int *displs,
349 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
350 return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
353 int PMPI_Iscatterv(void* sendbuf, int* sendcounts, int* displs, MPI_Datatype sendtype, void* recvbuf, int recvcount,
354 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
356 if (comm == MPI_COMM_NULL)
358 if (sendcounts == nullptr || displs == nullptr)
360 if (((comm->rank() == root) && (sendtype == MPI_DATATYPE_NULL)) ||
361 ((recvbuf != MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL)))
363 if (request == nullptr)
365 if (recvbuf != MPI_IN_PLACE && recvcount < 0)
366 return MPI_ERR_COUNT;
367 if (root < 0 || root >= comm->size())
370 if (comm->rank() == root) {
371 if (recvbuf == MPI_IN_PLACE) {
373 recvcount = sendcounts[comm->rank()];
375 for (int i = 0; i < comm->size(); i++) {
376 if (sendcounts[i] < 0)
377 return MPI_ERR_COUNT;
383 int rank = simgrid::s4u::this_actor::get_pid();
384 int dt_size_send = sendtype->is_replayable() ? 1 : sendtype->size();
386 std::vector<int>* trace_sendcounts = new std::vector<int>;
387 if (comm->rank() == root) {
388 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
389 trace_sendcounts->push_back(sendcounts[i] * dt_size_send);
393 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatterv" : "PMPI_Iscatterv",
394 new simgrid::instr::VarCollTIData(
395 request == MPI_REQUEST_IGNORED ? "scatterv" : "iscatterv", root, dt_size_send,
396 trace_sendcounts, recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
397 nullptr, simgrid::smpi::Datatype::encode(sendtype),
398 simgrid::smpi::Datatype::encode(recvtype)));
399 if (request == MPI_REQUEST_IGNORED)
400 simgrid::smpi::Colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
402 simgrid::smpi::Colls::iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm,
405 TRACE_smpi_comm_out(rank);
410 int PMPI_Reduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
412 return PMPI_Ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, MPI_REQUEST_IGNORED);
415 int PMPI_Ireduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm, MPI_Request* request)
417 if (comm == MPI_COMM_NULL)
419 if ((sendbuf == nullptr && count > 0) || ((comm->rank() == root) && recvbuf == nullptr))
420 return MPI_ERR_BUFFER;
421 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
423 if (op == MPI_OP_NULL)
425 if (request == nullptr)
427 if (root < 0 || root >= comm->size())
430 return MPI_ERR_COUNT;
433 int rank = simgrid::s4u::this_actor::get_pid();
435 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce" : "PMPI_Ireduce",
436 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "reduce" : "ireduce", root, 0,
437 datatype->is_replayable() ? count : count * datatype->size(), -1,
438 simgrid::smpi::Datatype::encode(datatype), ""));
439 if (request == MPI_REQUEST_IGNORED)
440 simgrid::smpi::Colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
442 simgrid::smpi::Colls::ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, request);
444 TRACE_smpi_comm_out(rank);
449 int PMPI_Reduce_local(void* inbuf, void* inoutbuf, int count, MPI_Datatype datatype, MPI_Op op)
451 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
453 if (op == MPI_OP_NULL)
456 return MPI_ERR_COUNT;
459 op->apply(inbuf, inoutbuf, &count, datatype);
464 int PMPI_Allreduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
466 return PMPI_Iallreduce(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
469 int PMPI_Iallreduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
471 if (comm == MPI_COMM_NULL)
473 if ((sendbuf == nullptr && count > 0) || (recvbuf == nullptr))
474 return MPI_ERR_BUFFER;
475 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
478 return MPI_ERR_COUNT;
479 if (op == MPI_OP_NULL)
481 if (request == nullptr)
485 char* sendtmpbuf = static_cast<char*>(sendbuf);
486 if (sendbuf == MPI_IN_PLACE) {
487 sendtmpbuf = static_cast<char*>(xbt_malloc(count * datatype->get_extent()));
488 simgrid::smpi::Datatype::copy(recvbuf, count, datatype, sendtmpbuf, count, datatype);
490 int rank = simgrid::s4u::this_actor::get_pid();
492 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allreduce" : "PMPI_Iallreduce",
493 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "allreduce" : "iallreduce", -1, 0,
494 datatype->is_replayable() ? count : count * datatype->size(), -1,
495 simgrid::smpi::Datatype::encode(datatype), ""));
497 if (request == MPI_REQUEST_IGNORED)
498 simgrid::smpi::Colls::allreduce(sendtmpbuf, recvbuf, count, datatype, op, comm);
500 simgrid::smpi::Colls::iallreduce(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
502 if (sendbuf == MPI_IN_PLACE)
503 xbt_free(sendtmpbuf);
505 TRACE_smpi_comm_out(rank);
510 int PMPI_Scan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
512 return PMPI_Iscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
515 int PMPI_Iscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request)
517 if (comm == MPI_COMM_NULL)
519 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
521 if (op == MPI_OP_NULL)
523 if (request == nullptr)
526 return MPI_ERR_COUNT;
527 if (sendbuf == nullptr || recvbuf == nullptr)
528 return MPI_ERR_BUFFER;
531 int rank = simgrid::s4u::this_actor::get_pid();
532 void* sendtmpbuf = sendbuf;
533 if (sendbuf == MPI_IN_PLACE) {
534 sendtmpbuf = static_cast<void*>(xbt_malloc(count * datatype->size()));
535 memcpy(sendtmpbuf, recvbuf, count * datatype->size());
537 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scan" : "PMPI_Iscan",
538 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "scan" : "iscan", -1,
539 datatype->is_replayable() ? count : count * datatype->size(),
540 simgrid::smpi::Datatype::encode(datatype)));
543 if (request == MPI_REQUEST_IGNORED)
544 retval = simgrid::smpi::Colls::scan(sendtmpbuf, recvbuf, count, datatype, op, comm);
546 retval = simgrid::smpi::Colls::iscan(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
548 TRACE_smpi_comm_out(rank);
549 if (sendbuf == MPI_IN_PLACE)
550 xbt_free(sendtmpbuf);
555 int PMPI_Exscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
557 return PMPI_Iexscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
560 int PMPI_Iexscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request){
561 if (comm == MPI_COMM_NULL)
563 if (not datatype->is_valid())
565 if (op == MPI_OP_NULL)
567 if (request == nullptr)
570 return MPI_ERR_COUNT;
571 if (sendbuf == nullptr || recvbuf == nullptr)
572 return MPI_ERR_BUFFER;
575 int rank = simgrid::s4u::this_actor::get_pid();
576 void* sendtmpbuf = sendbuf;
577 if (sendbuf == MPI_IN_PLACE) {
578 sendtmpbuf = static_cast<void*>(xbt_malloc(count * datatype->size()));
579 memcpy(sendtmpbuf, recvbuf, count * datatype->size());
582 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Exscan" : "PMPI_Iexscan",
583 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "exscan" : "iexscan", -1,
584 datatype->is_replayable() ? count : count * datatype->size(),
585 simgrid::smpi::Datatype::encode(datatype)));
588 if (request == MPI_REQUEST_IGNORED)
589 retval = simgrid::smpi::Colls::exscan(sendtmpbuf, recvbuf, count, datatype, op, comm);
591 retval = simgrid::smpi::Colls::iexscan(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
593 TRACE_smpi_comm_out(rank);
594 if (sendbuf == MPI_IN_PLACE)
595 xbt_free(sendtmpbuf);
600 int PMPI_Reduce_scatter(void *sendbuf, void *recvbuf, int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
602 return PMPI_Ireduce_scatter(sendbuf, recvbuf, recvcounts, datatype, op, comm, MPI_REQUEST_IGNORED);
605 int PMPI_Ireduce_scatter(void *sendbuf, void *recvbuf, int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
607 if (comm == MPI_COMM_NULL)
609 if ((sendbuf == nullptr) || (recvbuf == nullptr))
610 return MPI_ERR_BUFFER;
611 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
613 if (op == MPI_OP_NULL)
615 if (recvcounts == nullptr)
617 if (request == nullptr)
620 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
621 if (recvcounts[i] < 0)
622 return MPI_ERR_COUNT;
626 int rank = simgrid::s4u::this_actor::get_pid();
627 std::vector<int>* trace_recvcounts = new std::vector<int>;
628 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
631 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
632 trace_recvcounts->push_back(recvcounts[i] * dt_send_size);
633 totalcount += recvcounts[i];
636 void* sendtmpbuf = sendbuf;
637 if (sendbuf == MPI_IN_PLACE) {
638 sendtmpbuf = static_cast<void*>(xbt_malloc(totalcount * datatype->size()));
639 memcpy(sendtmpbuf, recvbuf, totalcount * datatype->size());
642 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter" : "PMPI_Ireduce_scatter",
643 new simgrid::instr::VarCollTIData(
644 request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, dt_send_size, nullptr,
645 -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
647 if (request == MPI_REQUEST_IGNORED)
648 simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
650 simgrid::smpi::Colls::ireduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm, request);
652 TRACE_smpi_comm_out(rank);
653 if (sendbuf == MPI_IN_PLACE)
654 xbt_free(sendtmpbuf);
659 int PMPI_Reduce_scatter_block(void *sendbuf, void *recvbuf, int recvcount,
660 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
662 return PMPI_Ireduce_scatter_block(sendbuf, recvbuf, recvcount, datatype, op, comm, MPI_REQUEST_IGNORED);
665 int PMPI_Ireduce_scatter_block(void* sendbuf, void* recvbuf, int recvcount, MPI_Datatype datatype, MPI_Op op,
666 MPI_Comm comm, MPI_Request* request)
668 if (comm == MPI_COMM_NULL)
670 if (not datatype->is_valid())
672 if (op == MPI_OP_NULL)
676 if (request == nullptr)
680 int count = comm->size();
682 int rank = simgrid::s4u::this_actor::get_pid();
683 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
684 std::vector<int>* trace_recvcounts = new std::vector<int>(recvcount * dt_send_size); // copy data to avoid bad free
686 void* sendtmpbuf = sendbuf;
687 if (sendbuf == MPI_IN_PLACE) {
688 sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * count * datatype->size()));
689 memcpy(sendtmpbuf, recvbuf, recvcount * count * datatype->size());
693 rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter_block" : "PMPI_Ireduce_scatter_block",
694 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, 0,
695 nullptr, -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
697 int* recvcounts = new int[count];
698 for (int i = 0; i < count; i++)
699 recvcounts[i] = recvcount;
700 if (request == MPI_REQUEST_IGNORED)
701 simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
703 simgrid::smpi::Colls::ireduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm, request);
706 TRACE_smpi_comm_out(rank);
707 if (sendbuf == MPI_IN_PLACE)
708 xbt_free(sendtmpbuf);
713 int PMPI_Alltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
714 MPI_Datatype recvtype, MPI_Comm comm){
715 return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
718 int PMPI_Ialltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
719 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
721 if (comm == MPI_COMM_NULL)
723 if ((sendbuf == nullptr && sendcount > 0) || (recvbuf == nullptr && recvcount > 0))
724 return MPI_ERR_BUFFER;
725 if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL)
727 if ((sendbuf != MPI_IN_PLACE && sendcount < 0) || recvcount < 0)
728 return MPI_ERR_COUNT;
729 if (request == nullptr)
733 int rank = simgrid::s4u::this_actor::get_pid();
734 void* sendtmpbuf = static_cast<char*>(sendbuf);
735 int sendtmpcount = sendcount;
736 MPI_Datatype sendtmptype = sendtype;
737 if (sendbuf == MPI_IN_PLACE) {
738 sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * comm->size() * recvtype->size()));
739 // memcpy(??,nullptr,0) is actually undefined behavor, even if harmless.
740 if (recvbuf != nullptr)
741 memcpy(sendtmpbuf, recvbuf, recvcount * comm->size() * recvtype->size());
742 sendtmpcount = recvcount;
743 sendtmptype = recvtype;
746 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoall" : "PMPI_Ialltoall",
747 new simgrid::instr::CollTIData(
748 request == MPI_REQUEST_IGNORED ? "alltoall" : "ialltoall", -1, -1.0,
749 sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(),
750 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
751 simgrid::smpi::Datatype::encode(sendtmptype), simgrid::smpi::Datatype::encode(recvtype)));
753 if (request == MPI_REQUEST_IGNORED)
754 retval = simgrid::smpi::Colls::alltoall(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, comm);
756 retval = simgrid::smpi::Colls::ialltoall(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, comm,
759 TRACE_smpi_comm_out(rank);
760 if (sendbuf == MPI_IN_PLACE)
761 xbt_free(sendtmpbuf);
766 int PMPI_Alltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf,
767 int* recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm)
769 return PMPI_Ialltoallv(sendbuf, sendcounts, senddisps, sendtype, recvbuf, recvcounts, recvdisps, recvtype, comm, MPI_REQUEST_IGNORED);
772 int PMPI_Ialltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf,
773 int* recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
775 if (comm == MPI_COMM_NULL)
777 if (sendbuf == nullptr || recvbuf == nullptr)
778 return MPI_ERR_BUFFER;
779 if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL)
781 if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
782 recvdisps == nullptr)
784 if (request == nullptr)
787 int rank = simgrid::s4u::this_actor::get_pid();
788 int size = comm->size();
789 for (int i = 0; i < size; i++) {
790 if (recvcounts[i] < 0 || (sendbuf != MPI_IN_PLACE && sendcounts[i] < 0))
791 return MPI_ERR_COUNT;
797 std::vector<int>* trace_sendcounts = new std::vector<int>;
798 std::vector<int>* trace_recvcounts = new std::vector<int>;
799 int dt_size_recv = recvtype->size();
801 void* sendtmpbuf = static_cast<char*>(sendbuf);
802 int* sendtmpcounts = sendcounts;
803 int* sendtmpdisps = senddisps;
804 MPI_Datatype sendtmptype = sendtype;
806 for (int i = 0; i < size; i++) { // copy data to avoid bad free
807 recv_size += recvcounts[i] * dt_size_recv;
808 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
809 if (((recvdisps[i] + recvcounts[i]) * dt_size_recv) > maxsize)
810 maxsize = (recvdisps[i] + recvcounts[i]) * dt_size_recv;
813 if (sendbuf == MPI_IN_PLACE) {
814 sendtmpbuf = static_cast<void*>(xbt_malloc(maxsize));
815 memcpy(sendtmpbuf, recvbuf, maxsize);
816 sendtmpcounts = static_cast<int*>(xbt_malloc(size * sizeof(int)));
817 memcpy(sendtmpcounts, recvcounts, size * sizeof(int));
818 sendtmpdisps = static_cast<int*>(xbt_malloc(size * sizeof(int)));
819 memcpy(sendtmpdisps, recvdisps, size * sizeof(int));
820 sendtmptype = recvtype;
823 int dt_size_send = sendtmptype->size();
825 for (int i = 0; i < size; i++) { // copy data to avoid bad free
826 send_size += sendtmpcounts[i] * dt_size_send;
827 trace_sendcounts->push_back(sendtmpcounts[i] * dt_size_send);
830 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallv" : "PMPI_Ialltoallv",
831 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
832 send_size, trace_sendcounts, recv_size, trace_recvcounts,
833 simgrid::smpi::Datatype::encode(sendtmptype),
834 simgrid::smpi::Datatype::encode(recvtype)));
837 if (request == MPI_REQUEST_IGNORED)
838 retval = simgrid::smpi::Colls::alltoallv(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptype, recvbuf, recvcounts,
839 recvdisps, recvtype, comm);
841 retval = simgrid::smpi::Colls::ialltoallv(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptype, recvbuf, recvcounts,
842 recvdisps, recvtype, comm, request);
844 TRACE_smpi_comm_out(rank);
845 if (sendbuf == MPI_IN_PLACE) {
846 xbt_free(sendtmpbuf);
847 xbt_free(sendtmpcounts);
848 xbt_free(sendtmpdisps);
854 int PMPI_Alltoallw(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype* sendtypes, void* recvbuf,
855 int* recvcounts, int* recvdisps, MPI_Datatype* recvtypes, MPI_Comm comm)
857 return PMPI_Ialltoallw(sendbuf, sendcounts, senddisps, sendtypes, recvbuf, recvcounts, recvdisps, recvtypes, comm, MPI_REQUEST_IGNORED);
860 int PMPI_Ialltoallw(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype* sendtypes, void* recvbuf,
861 int* recvcounts, int* recvdisps, MPI_Datatype* recvtypes, MPI_Comm comm, MPI_Request* request)
863 if (comm == MPI_COMM_NULL)
865 if (sendbuf == nullptr || recvbuf == nullptr)
866 return MPI_ERR_BUFFER;
867 if ((sendbuf != MPI_IN_PLACE && sendtypes == nullptr) || recvtypes == nullptr)
869 if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
870 recvdisps == nullptr)
872 if (request == nullptr)
876 int rank = simgrid::s4u::this_actor::get_pid();
877 int size = comm->size();
878 for (int i = 0; i < size; i++) { // copy data to avoid bad free
879 if (recvcounts[i] < 0 || (sendbuf != MPI_IN_PLACE && sendcounts[i] < 0))
880 return MPI_ERR_COUNT;
884 std::vector<int>* trace_sendcounts = new std::vector<int>;
885 std::vector<int>* trace_recvcounts = new std::vector<int>;
887 void* sendtmpbuf = static_cast<char*>(sendbuf);
888 int* sendtmpcounts = sendcounts;
889 int* sendtmpdisps = senddisps;
890 MPI_Datatype* sendtmptypes = sendtypes;
891 unsigned long maxsize = 0;
892 for (int i = 0; i < size; i++) { // copy data to avoid bad free
893 if (recvtypes[i] == MPI_DATATYPE_NULL) {
894 delete trace_recvcounts;
895 delete trace_sendcounts;
898 recv_size += recvcounts[i] * recvtypes[i]->size();
899 trace_recvcounts->push_back(recvcounts[i] * recvtypes[i]->size());
900 if ((recvdisps[i] + (recvcounts[i] * recvtypes[i]->size())) > maxsize)
901 maxsize = recvdisps[i] + (recvcounts[i] * recvtypes[i]->size());
904 if (sendbuf == MPI_IN_PLACE) {
905 sendtmpbuf = static_cast<void*>(xbt_malloc(maxsize));
906 memcpy(sendtmpbuf, recvbuf, maxsize);
907 sendtmpcounts = static_cast<int*>(xbt_malloc(size * sizeof(int)));
908 memcpy(sendtmpcounts, recvcounts, size * sizeof(int));
909 sendtmpdisps = static_cast<int*>(xbt_malloc(size * sizeof(int)));
910 memcpy(sendtmpdisps, recvdisps, size * sizeof(int));
911 sendtmptypes = static_cast<MPI_Datatype*>(xbt_malloc(size * sizeof(MPI_Datatype)));
912 memcpy(sendtmptypes, recvtypes, size * sizeof(MPI_Datatype));
915 for (int i = 0; i < size; i++) { // copy data to avoid bad free
916 send_size += sendtmpcounts[i] * sendtmptypes[i]->size();
917 trace_sendcounts->push_back(sendtmpcounts[i] * sendtmptypes[i]->size());
920 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallw" : "PMPI_Ialltoallw",
921 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
922 send_size, trace_sendcounts, recv_size, trace_recvcounts,
923 simgrid::smpi::Datatype::encode(sendtmptypes[0]),
924 simgrid::smpi::Datatype::encode(recvtypes[0])));
927 if (request == MPI_REQUEST_IGNORED)
928 retval = simgrid::smpi::Colls::alltoallw(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptypes, recvbuf, recvcounts,
929 recvdisps, recvtypes, comm);
931 retval = simgrid::smpi::Colls::ialltoallw(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptypes, recvbuf,
932 recvcounts, recvdisps, recvtypes, comm, request);
934 TRACE_smpi_comm_out(rank);
935 if (sendbuf == MPI_IN_PLACE) {
936 xbt_free(sendtmpbuf);
937 xbt_free(sendtmpcounts);
938 xbt_free(sendtmpdisps);
939 xbt_free(sendtmptypes);