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 #define CHECK_ARGS(test, errcode, ...) \
18 XBT_WARN(__VA_ARGS__); \
22 /* PMPI User level calls */
24 int PMPI_Barrier(MPI_Comm comm)
26 return PMPI_Ibarrier(comm, MPI_REQUEST_IGNORED);
29 int PMPI_Ibarrier(MPI_Comm comm, MPI_Request *request)
31 if (comm == MPI_COMM_NULL)
33 if (request == nullptr)
37 int rank = simgrid::s4u::this_actor::get_pid();
38 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Barrier" : "PMPI_Ibarrier",
39 new simgrid::instr::NoOpTIData(request == MPI_REQUEST_IGNORED ? "barrier" : "ibarrier"));
40 if (request == MPI_REQUEST_IGNORED) {
41 simgrid::smpi::Colls::barrier(comm);
42 // Barrier can be used to synchronize RMA calls. Finish all requests from comm before.
43 comm->finish_rma_calls();
45 simgrid::smpi::Colls::ibarrier(comm, request);
47 TRACE_smpi_comm_out(rank);
52 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
54 return PMPI_Ibcast(buf, count, datatype, root, comm, MPI_REQUEST_IGNORED);
57 int PMPI_Ibcast(void *buf, int count, MPI_Datatype datatype,
58 int root, MPI_Comm comm, MPI_Request* request)
60 if (comm == MPI_COMM_NULL)
62 if (buf == nullptr && count > 0)
63 return MPI_ERR_BUFFER;
64 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
68 if (root < 0 || root >= comm->size())
70 if (request == nullptr)
74 int rank = simgrid::s4u::this_actor::get_pid();
75 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Bcast" : "PMPI_Ibcast",
76 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "bcast" : "ibcast", root, -1.0,
77 datatype->is_replayable() ? count : count * datatype->size(), -1,
78 simgrid::smpi::Datatype::encode(datatype), ""));
79 if (comm->size() > 1) {
80 if (request == MPI_REQUEST_IGNORED)
81 simgrid::smpi::Colls::bcast(buf, count, datatype, root, comm);
83 simgrid::smpi::Colls::ibcast(buf, count, datatype, root, comm, request);
85 if (request != MPI_REQUEST_IGNORED)
86 *request = MPI_REQUEST_NULL;
89 TRACE_smpi_comm_out(rank);
94 int PMPI_Gather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
95 int root, MPI_Comm comm){
96 return PMPI_Igather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
99 int PMPI_Igather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
100 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
102 if (comm == MPI_COMM_NULL)
104 if ((sendbuf == nullptr && sendcount > 0) || ((comm->rank() == root) && recvbuf == nullptr && recvcount > 0))
105 return MPI_ERR_BUFFER;
106 if (((sendbuf != MPI_IN_PLACE && sendcount > 0) && (sendtype == MPI_DATATYPE_NULL)) ||
107 ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL)))
109 if (((sendbuf != MPI_IN_PLACE) && (sendcount < 0)) || ((comm->rank() == root) && (recvcount < 0)))
110 return MPI_ERR_COUNT;
111 if (root < 0 || root >= comm->size())
113 if (request == nullptr)
117 const void* real_sendbuf = sendbuf;
118 int real_sendcount = sendcount;
119 MPI_Datatype real_sendtype = sendtype;
120 if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
122 real_sendtype = recvtype;
124 int rank = simgrid::s4u::this_actor::get_pid();
126 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gather" : "PMPI_Igather",
127 new simgrid::instr::CollTIData(
128 request == MPI_REQUEST_IGNORED ? "gather" : "igather", root, -1.0,
129 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
130 (comm->rank() != root || recvtype->is_replayable()) ? recvcount : recvcount * recvtype->size(),
131 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
132 if (request == MPI_REQUEST_IGNORED)
133 simgrid::smpi::Colls::gather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm);
135 simgrid::smpi::Colls::igather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm,
138 TRACE_smpi_comm_out(rank);
143 int PMPI_Gatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, const int *recvcounts, const int *displs,
144 MPI_Datatype recvtype, int root, MPI_Comm comm){
145 return PMPI_Igatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, root, comm, MPI_REQUEST_IGNORED);
148 int PMPI_Igatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
149 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
151 if (comm == MPI_COMM_NULL)
153 if ((sendbuf == nullptr && sendcount > 0) || ((comm->rank() == root) && recvbuf == nullptr))
154 return MPI_ERR_BUFFER;
155 if (((sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
156 ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL)))
158 if ((sendbuf != MPI_IN_PLACE) && (sendcount < 0))
159 return MPI_ERR_COUNT;
160 if ((comm->rank() == root) && (recvcounts == nullptr || displs == nullptr))
162 if (root < 0 || root >= comm->size())
164 if (request == nullptr)
167 for (int i = 0; i < comm->size(); i++) {
168 if ((comm->rank() == root) && (recvcounts[i] < 0))
169 return MPI_ERR_COUNT;
173 const void* real_sendbuf = sendbuf;
174 int real_sendcount = sendcount;
175 MPI_Datatype real_sendtype = sendtype;
176 if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
178 real_sendtype = recvtype;
181 int rank = simgrid::s4u::this_actor::get_pid();
182 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
184 std::vector<int>* trace_recvcounts = new std::vector<int>;
185 if (comm->rank() == root) {
186 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
187 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
190 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gatherv" : "PMPI_Igatherv",
191 new simgrid::instr::VarCollTIData(
192 request == MPI_REQUEST_IGNORED ? "gatherv" : "igatherv", root,
193 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
194 nullptr, dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(real_sendtype),
195 simgrid::smpi::Datatype::encode(recvtype)));
196 if (request == MPI_REQUEST_IGNORED)
197 simgrid::smpi::Colls::gatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
200 simgrid::smpi::Colls::igatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
201 root, comm, request);
203 TRACE_smpi_comm_out(rank);
208 int PMPI_Allgather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
209 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm){
210 return PMPI_Iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
213 int PMPI_Iallgather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
214 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
216 CHECK_ARGS(comm == MPI_COMM_NULL, MPI_ERR_COMM, "Iallgather: the communicator cannot be MPI_COMM_NULL");
217 CHECK_ARGS(recvbuf == nullptr && recvcount > 0, MPI_ERR_BUFFER, "Iallgather: param 4 recvbuf cannot be NULL");
218 CHECK_ARGS(sendbuf == nullptr && sendcount > 0, MPI_ERR_BUFFER,
219 "Iallgather: param 1 sendbuf cannot be NULL when sendcound > 0");
220 CHECK_ARGS((sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL), MPI_ERR_TYPE,
221 "Iallgather: param 3 sendtype cannot be MPI_DATATYPE_NULL when sendbuff is not MPI_IN_PLACE");
222 CHECK_ARGS(recvtype == MPI_DATATYPE_NULL, MPI_ERR_TYPE, "Iallgather: param 6 recvtype cannot be MPI_DATATYPE_NULL");
223 CHECK_ARGS(recvcount < 0, MPI_ERR_COUNT, "Iallgather: param 5 recvcount cannot be negative");
224 CHECK_ARGS((sendbuf != MPI_IN_PLACE) && (sendcount < 0), MPI_ERR_COUNT,
225 "Iallgather: param 2 sendcount cannot be negative when sendbuf is not MPI_IN_PLACE");
226 CHECK_ARGS(request == nullptr, MPI_ERR_ARG, "Iallgather: param 8 request cannot be NULL");
229 if (sendbuf == MPI_IN_PLACE) {
230 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * recvcount * comm->rank();
231 sendcount = recvcount;
234 int rank = simgrid::s4u::this_actor::get_pid();
236 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgather" : "PMPI_Iallggather",
237 new simgrid::instr::CollTIData(
238 request == MPI_REQUEST_IGNORED ? "allgather" : "iallgather", -1, -1.0,
239 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
240 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
241 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
242 if (request == MPI_REQUEST_IGNORED)
243 simgrid::smpi::Colls::allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
245 simgrid::smpi::Colls::iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, request);
247 TRACE_smpi_comm_out(rank);
252 int PMPI_Allgatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
253 void *recvbuf, const int *recvcounts, const int *displs, MPI_Datatype recvtype, MPI_Comm comm){
254 return PMPI_Iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, MPI_REQUEST_IGNORED);
257 int PMPI_Iallgatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
258 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
260 if (comm == MPI_COMM_NULL)
262 if ((sendbuf == nullptr && sendcount > 0) || (recvbuf == nullptr))
263 return MPI_ERR_BUFFER;
264 if (((sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) || (recvtype == MPI_DATATYPE_NULL))
266 if ((sendbuf != MPI_IN_PLACE) && (sendcount < 0))
267 return MPI_ERR_COUNT;
268 if (recvcounts == nullptr || displs == nullptr)
270 if (request == nullptr)
273 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
274 if (recvcounts[i] < 0)
275 return MPI_ERR_COUNT;
279 if (sendbuf == MPI_IN_PLACE) {
280 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * displs[comm->rank()];
281 sendcount = recvcounts[comm->rank()];
284 int rank = simgrid::s4u::this_actor::get_pid();
285 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
287 std::vector<int>* trace_recvcounts = new std::vector<int>;
288 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
289 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
293 rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgatherv" : "PMPI_Iallgatherv",
294 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "allgatherv" : "iallgatherv", -1,
295 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(), nullptr,
296 dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
297 simgrid::smpi::Datatype::encode(recvtype)));
298 if (request == MPI_REQUEST_IGNORED)
299 simgrid::smpi::Colls::allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
301 simgrid::smpi::Colls::iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm,
304 TRACE_smpi_comm_out(rank);
309 int PMPI_Scatter(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
310 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
311 return PMPI_Iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
314 int PMPI_Iscatter(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
315 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
317 if (comm == MPI_COMM_NULL)
319 if (((comm->rank() == root) && (sendtype == MPI_DATATYPE_NULL || not sendtype->is_valid())) ||
320 ((recvbuf != MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL || not recvtype->is_valid())))
322 if (((comm->rank() == root) && (sendcount < 0)) || ((recvbuf != MPI_IN_PLACE) && (recvcount < 0)))
323 return MPI_ERR_COUNT;
324 if ((sendbuf == recvbuf) || ((comm->rank() == root) && sendcount > 0 && (sendbuf == nullptr)) ||
325 (recvcount > 0 && recvbuf == nullptr))
326 return MPI_ERR_BUFFER;
327 if (root < 0 || root >= comm->size())
329 if (request == nullptr)
333 if (recvbuf == MPI_IN_PLACE) {
335 recvcount = sendcount;
337 int rank = simgrid::s4u::this_actor::get_pid();
339 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatter" : "PMPI_Iscatter",
340 new simgrid::instr::CollTIData(
341 request == MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
342 (comm->rank() != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
343 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
344 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
345 if (request == MPI_REQUEST_IGNORED)
346 simgrid::smpi::Colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
348 simgrid::smpi::Colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
350 TRACE_smpi_comm_out(rank);
355 int PMPI_Scatterv(const void *sendbuf, const int *sendcounts, const int *displs,
356 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
357 return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
360 int PMPI_Iscatterv(const void* sendbuf, const int* sendcounts, const int* displs, MPI_Datatype sendtype, void* recvbuf, int recvcount,
361 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
363 CHECK_ARGS(comm == MPI_COMM_NULL, MPI_ERR_COMM, "Iscatterv: the communicator cannot be MPI_COMM_NULL");
364 CHECK_ARGS((comm->rank() == root) && (sendcounts == nullptr), MPI_ERR_ARG,
365 "Iscatterv: param 2 sendcounts cannot be NULL on the root rank");
366 CHECK_ARGS((comm->rank() == root) && (displs == nullptr), MPI_ERR_ARG,
367 "Iscatterv: param 3 displs cannot be NULL on the root rank");
368 CHECK_ARGS((comm->rank() == root) && (sendtype == MPI_DATATYPE_NULL), MPI_ERR_TYPE,
369 "Iscatterv: The sendtype cannot be NULL on the root rank");
370 CHECK_ARGS((recvbuf != MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL), MPI_ERR_TYPE,
371 "Iscatterv: the recvtype cannot be NULL when not receiving in place");
372 CHECK_ARGS(request == nullptr, MPI_ERR_ARG, "Iscatterv: param 10 request cannot be NULL");
373 CHECK_ARGS(recvbuf != MPI_IN_PLACE && recvcount < 0, MPI_ERR_COUNT,
374 "Iscatterv: When not receiving in place, the recvcound cannot be negative");
375 CHECK_ARGS(root < 0, MPI_ERR_ROOT, "Iscatterv: root cannot be negative");
376 CHECK_ARGS(root >= comm->size(), MPI_ERR_ROOT, "Iscatterv: root (=%d) is larger than communicator size (=%d)", root,
379 if (comm->rank() == root) {
380 if (recvbuf == MPI_IN_PLACE) {
382 recvcount = sendcounts[comm->rank()];
384 for (int i = 0; i < comm->size(); i++)
385 CHECK_ARGS(sendcounts[i] < 0, MPI_ERR_COUNT, "Iscatterv: sendcounts[%d]=%d but this cannot be negative", i,
391 int rank = simgrid::s4u::this_actor::get_pid();
392 int dt_size_send = sendtype->is_replayable() ? 1 : sendtype->size();
394 std::vector<int>* trace_sendcounts = new std::vector<int>;
395 if (comm->rank() == root) {
396 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
397 trace_sendcounts->push_back(sendcounts[i] * dt_size_send);
401 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatterv" : "PMPI_Iscatterv",
402 new simgrid::instr::VarCollTIData(
403 request == MPI_REQUEST_IGNORED ? "scatterv" : "iscatterv", root, dt_size_send,
404 trace_sendcounts, recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
405 nullptr, simgrid::smpi::Datatype::encode(sendtype),
406 simgrid::smpi::Datatype::encode(recvtype)));
407 if (request == MPI_REQUEST_IGNORED)
408 simgrid::smpi::Colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
410 simgrid::smpi::Colls::iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm,
413 TRACE_smpi_comm_out(rank);
418 int PMPI_Reduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
420 return PMPI_Ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, MPI_REQUEST_IGNORED);
423 int PMPI_Ireduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm, MPI_Request* request)
425 if (comm == MPI_COMM_NULL)
427 if ((sendbuf == nullptr && count > 0) || ((comm->rank() == root) && recvbuf == nullptr))
428 return MPI_ERR_BUFFER;
429 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
431 if (op == MPI_OP_NULL)
433 if (request == nullptr)
435 if (root < 0 || root >= comm->size())
438 return MPI_ERR_COUNT;
441 int rank = simgrid::s4u::this_actor::get_pid();
443 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce" : "PMPI_Ireduce",
444 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "reduce" : "ireduce", root, 0,
445 datatype->is_replayable() ? count : count * datatype->size(), -1,
446 simgrid::smpi::Datatype::encode(datatype), ""));
447 if (request == MPI_REQUEST_IGNORED)
448 simgrid::smpi::Colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
450 simgrid::smpi::Colls::ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, request);
452 TRACE_smpi_comm_out(rank);
457 int PMPI_Reduce_local(const void* inbuf, void* inoutbuf, int count, MPI_Datatype datatype, MPI_Op op)
459 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
461 if (op == MPI_OP_NULL)
464 return MPI_ERR_COUNT;
467 op->apply(inbuf, inoutbuf, &count, datatype);
472 int PMPI_Allreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
474 return PMPI_Iallreduce(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
477 int PMPI_Iallreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
479 if (comm == MPI_COMM_NULL)
481 if ((sendbuf == nullptr && count > 0) || (recvbuf == nullptr))
482 return MPI_ERR_BUFFER;
483 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
486 return MPI_ERR_COUNT;
487 if (op == MPI_OP_NULL)
489 if (request == nullptr)
493 const void* real_sendbuf = sendbuf;
494 std::unique_ptr<unsigned char[]> tmp_sendbuf;
495 if (sendbuf == MPI_IN_PLACE) {
496 tmp_sendbuf.reset(new unsigned char[count * datatype->get_extent()]);
497 simgrid::smpi::Datatype::copy(recvbuf, count, datatype, tmp_sendbuf.get(), count, datatype);
498 real_sendbuf = tmp_sendbuf.get();
500 int rank = simgrid::s4u::this_actor::get_pid();
502 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allreduce" : "PMPI_Iallreduce",
503 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "allreduce" : "iallreduce", -1, 0,
504 datatype->is_replayable() ? count : count * datatype->size(), -1,
505 simgrid::smpi::Datatype::encode(datatype), ""));
507 if (request == MPI_REQUEST_IGNORED)
508 simgrid::smpi::Colls::allreduce(real_sendbuf, recvbuf, count, datatype, op, comm);
510 simgrid::smpi::Colls::iallreduce(real_sendbuf, recvbuf, count, datatype, op, comm, request);
512 TRACE_smpi_comm_out(rank);
517 int PMPI_Scan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
519 return PMPI_Iscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
522 int PMPI_Iscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request)
524 if (comm == MPI_COMM_NULL)
526 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
528 if (op == MPI_OP_NULL)
530 if (request == nullptr)
533 return MPI_ERR_COUNT;
534 if (sendbuf == nullptr || recvbuf == nullptr)
535 return MPI_ERR_BUFFER;
538 int rank = simgrid::s4u::this_actor::get_pid();
539 const void* real_sendbuf = sendbuf;
540 std::unique_ptr<unsigned char[]> tmp_sendbuf;
541 if (sendbuf == MPI_IN_PLACE) {
542 tmp_sendbuf.reset(new unsigned char[count * datatype->size()]);
543 real_sendbuf = memcpy(tmp_sendbuf.get(), recvbuf, count * datatype->size());
545 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scan" : "PMPI_Iscan",
546 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "scan" : "iscan", -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::scan(real_sendbuf, recvbuf, count, datatype, op, comm);
554 retval = simgrid::smpi::Colls::iscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
556 TRACE_smpi_comm_out(rank);
561 int PMPI_Exscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
563 return PMPI_Iexscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
566 int PMPI_Iexscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request){
567 if (comm == MPI_COMM_NULL)
569 if (not datatype->is_valid())
571 if (op == MPI_OP_NULL)
573 if (request == nullptr)
576 return MPI_ERR_COUNT;
577 if (sendbuf == nullptr || recvbuf == nullptr)
578 return MPI_ERR_BUFFER;
581 int rank = simgrid::s4u::this_actor::get_pid();
582 const void* real_sendbuf = sendbuf;
583 std::unique_ptr<unsigned char[]> tmp_sendbuf;
584 if (sendbuf == MPI_IN_PLACE) {
585 tmp_sendbuf.reset(new unsigned char[count * datatype->size()]);
586 real_sendbuf = memcpy(tmp_sendbuf.get(), recvbuf, count * datatype->size());
589 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Exscan" : "PMPI_Iexscan",
590 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "exscan" : "iexscan", -1,
591 datatype->is_replayable() ? count : count * datatype->size(),
592 simgrid::smpi::Datatype::encode(datatype)));
595 if (request == MPI_REQUEST_IGNORED)
596 retval = simgrid::smpi::Colls::exscan(real_sendbuf, recvbuf, count, datatype, op, comm);
598 retval = simgrid::smpi::Colls::iexscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
600 TRACE_smpi_comm_out(rank);
605 int PMPI_Reduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
607 return PMPI_Ireduce_scatter(sendbuf, recvbuf, recvcounts, datatype, op, comm, MPI_REQUEST_IGNORED);
610 int PMPI_Ireduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
612 if (comm == MPI_COMM_NULL)
614 if ((sendbuf == nullptr) || (recvbuf == nullptr))
615 return MPI_ERR_BUFFER;
616 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
618 if (op == MPI_OP_NULL)
620 if (recvcounts == nullptr)
622 if (request == nullptr)
625 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
626 if (recvcounts[i] < 0)
627 return MPI_ERR_COUNT;
631 int rank = simgrid::s4u::this_actor::get_pid();
632 std::vector<int>* trace_recvcounts = new std::vector<int>;
633 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
636 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
637 trace_recvcounts->push_back(recvcounts[i] * dt_send_size);
638 totalcount += recvcounts[i];
641 const void* real_sendbuf = sendbuf;
642 std::unique_ptr<unsigned char[]> tmp_sendbuf;
643 if (sendbuf == MPI_IN_PLACE) {
644 tmp_sendbuf.reset(new unsigned char[totalcount * datatype->size()]);
645 real_sendbuf = memcpy(tmp_sendbuf.get(), recvbuf, totalcount * datatype->size());
648 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter" : "PMPI_Ireduce_scatter",
649 new simgrid::instr::VarCollTIData(
650 request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, dt_send_size, nullptr,
651 -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
653 if (request == MPI_REQUEST_IGNORED)
654 simgrid::smpi::Colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm);
656 simgrid::smpi::Colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm, request);
658 TRACE_smpi_comm_out(rank);
663 int PMPI_Reduce_scatter_block(const void *sendbuf, void *recvbuf, int recvcount,
664 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
666 return PMPI_Ireduce_scatter_block(sendbuf, recvbuf, recvcount, datatype, op, comm, MPI_REQUEST_IGNORED);
669 int PMPI_Ireduce_scatter_block(const void* sendbuf, void* recvbuf, int recvcount, MPI_Datatype datatype, MPI_Op op,
670 MPI_Comm comm, MPI_Request* request)
672 if (comm == MPI_COMM_NULL)
674 if (not datatype->is_valid())
676 if (op == MPI_OP_NULL)
680 if (request == nullptr)
684 int count = comm->size();
686 int rank = simgrid::s4u::this_actor::get_pid();
687 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
688 std::vector<int>* trace_recvcounts = new std::vector<int>(recvcount * dt_send_size); // copy data to avoid bad free
690 const void* real_sendbuf = sendbuf;
691 std::unique_ptr<unsigned char[]> tmp_sendbuf;
692 if (sendbuf == MPI_IN_PLACE) {
693 tmp_sendbuf.reset(new unsigned char[recvcount * count * datatype->size()]);
694 real_sendbuf = memcpy(tmp_sendbuf.get(), recvbuf, recvcount * count * datatype->size());
698 rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter_block" : "PMPI_Ireduce_scatter_block",
699 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, 0,
700 nullptr, -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
702 int* recvcounts = new int[count];
703 for (int i = 0; i < count; i++)
704 recvcounts[i] = recvcount;
705 if (request == MPI_REQUEST_IGNORED)
706 simgrid::smpi::Colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm);
708 simgrid::smpi::Colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm, request);
711 TRACE_smpi_comm_out(rank);
716 int PMPI_Alltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
717 MPI_Datatype recvtype, MPI_Comm comm){
718 return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
721 int PMPI_Ialltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
722 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
724 if (comm == MPI_COMM_NULL)
726 if ((sendbuf == nullptr && sendcount > 0) || (recvbuf == nullptr && recvcount > 0))
727 return MPI_ERR_BUFFER;
728 if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL)
730 if ((sendbuf != MPI_IN_PLACE && sendcount < 0) || recvcount < 0)
731 return MPI_ERR_COUNT;
732 if (request == nullptr)
736 int rank = simgrid::s4u::this_actor::get_pid();
737 const void* real_sendbuf = sendbuf;
738 int real_sendcount = sendcount;
739 MPI_Datatype real_sendtype = sendtype;
740 std::unique_ptr<unsigned char[]> tmp_sendbuf;
741 if (sendbuf == MPI_IN_PLACE) {
742 tmp_sendbuf.reset(new unsigned char[recvcount * comm->size() * recvtype->size()]);
743 // memcpy(??,nullptr,0) is actually undefined behavor, even if harmless.
744 if (recvbuf != nullptr)
745 memcpy(tmp_sendbuf.get(), recvbuf, recvcount * comm->size() * recvtype->size());
746 real_sendbuf = tmp_sendbuf.get();
747 real_sendcount = recvcount;
748 real_sendtype = recvtype;
751 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoall" : "PMPI_Ialltoall",
752 new simgrid::instr::CollTIData(
753 request == MPI_REQUEST_IGNORED ? "alltoall" : "ialltoall", -1, -1.0,
754 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
755 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
756 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
758 if (request == MPI_REQUEST_IGNORED)
760 simgrid::smpi::Colls::alltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, comm);
762 retval = simgrid::smpi::Colls::ialltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype,
765 TRACE_smpi_comm_out(rank);
770 int PMPI_Alltoallv(const void* sendbuf, const int* sendcounts, const int* senddisps, MPI_Datatype sendtype, void* recvbuf,
771 const int* recvcounts, const int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm)
773 return PMPI_Ialltoallv(sendbuf, sendcounts, senddisps, sendtype, recvbuf, recvcounts, recvdisps, recvtype, comm, MPI_REQUEST_IGNORED);
776 int PMPI_Ialltoallv(const void* sendbuf, const int* sendcounts, const int* senddisps, MPI_Datatype sendtype, void* recvbuf,
777 const int* recvcounts, const int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
779 if (comm == MPI_COMM_NULL)
781 if (sendbuf == nullptr || recvbuf == nullptr)
782 return MPI_ERR_BUFFER;
783 if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL)
785 if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
786 recvdisps == nullptr)
788 if (request == nullptr)
791 int rank = simgrid::s4u::this_actor::get_pid();
792 int size = comm->size();
793 for (int i = 0; i < size; i++) {
794 if (recvcounts[i] < 0 || (sendbuf != MPI_IN_PLACE && sendcounts[i] < 0))
795 return MPI_ERR_COUNT;
801 std::vector<int>* trace_sendcounts = new std::vector<int>;
802 std::vector<int>* trace_recvcounts = new std::vector<int>;
803 int dt_size_recv = recvtype->size();
805 const void* real_sendbuf = sendbuf;
806 const int* real_sendcounts = sendcounts;
807 const int* real_senddisps = senddisps;
808 MPI_Datatype real_sendtype = sendtype;
810 for (int i = 0; i < size; i++) { // copy data to avoid bad free
811 recv_size += recvcounts[i] * dt_size_recv;
812 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
813 if (((recvdisps[i] + recvcounts[i]) * dt_size_recv) > maxsize)
814 maxsize = (recvdisps[i] + recvcounts[i]) * dt_size_recv;
817 std::unique_ptr<unsigned char[]> tmp_sendbuf;
818 std::unique_ptr<int[]> tmp_sendcounts;
819 std::unique_ptr<int[]> tmp_senddisps;
820 if (sendbuf == MPI_IN_PLACE) {
821 tmp_sendbuf.reset(new unsigned char[maxsize]);
822 real_sendbuf = memcpy(tmp_sendbuf.get(), recvbuf, maxsize);
823 tmp_sendcounts.reset(new int[size]);
824 std::copy(recvcounts, recvcounts + size, tmp_sendcounts.get());
825 real_sendcounts = tmp_sendcounts.get();
826 tmp_senddisps.reset(new int[size]);
827 std::copy(recvdisps, recvdisps + size, tmp_senddisps.get());
828 real_senddisps = tmp_senddisps.get();
829 real_sendtype = recvtype;
832 int dt_size_send = real_sendtype->size();
834 for (int i = 0; i < size; i++) { // copy data to avoid bad free
835 send_size += real_sendcounts[i] * dt_size_send;
836 trace_sendcounts->push_back(real_sendcounts[i] * dt_size_send);
839 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallv" : "PMPI_Ialltoallv",
840 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
841 send_size, trace_sendcounts, recv_size, trace_recvcounts,
842 simgrid::smpi::Datatype::encode(real_sendtype),
843 simgrid::smpi::Datatype::encode(recvtype)));
846 if (request == MPI_REQUEST_IGNORED)
847 retval = simgrid::smpi::Colls::alltoallv(real_sendbuf, real_sendcounts, real_senddisps, real_sendtype, recvbuf,
848 recvcounts, recvdisps, recvtype, comm);
850 retval = simgrid::smpi::Colls::ialltoallv(real_sendbuf, real_sendcounts, real_senddisps, real_sendtype, recvbuf,
851 recvcounts, recvdisps, recvtype, comm, request);
853 TRACE_smpi_comm_out(rank);
858 int PMPI_Alltoallw(const void* sendbuf, const int* sendcounts, const int* senddisps, const MPI_Datatype* sendtypes, void* recvbuf,
859 const int* recvcounts, const int* recvdisps, const MPI_Datatype* recvtypes, MPI_Comm comm)
861 return PMPI_Ialltoallw(sendbuf, sendcounts, senddisps, sendtypes, recvbuf, recvcounts, recvdisps, recvtypes, comm, MPI_REQUEST_IGNORED);
864 int PMPI_Ialltoallw(const void* sendbuf, const int* sendcounts, const int* senddisps, const MPI_Datatype* sendtypes, void* recvbuf,
865 const int* recvcounts, const int* recvdisps, const MPI_Datatype* recvtypes, MPI_Comm comm, MPI_Request* request)
867 if (comm == MPI_COMM_NULL)
869 if (sendbuf == nullptr || recvbuf == nullptr)
870 return MPI_ERR_BUFFER;
871 if ((sendbuf != MPI_IN_PLACE && sendtypes == nullptr) || recvtypes == nullptr)
873 if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
874 recvdisps == nullptr)
876 if (request == nullptr)
880 int rank = simgrid::s4u::this_actor::get_pid();
881 int size = comm->size();
882 for (int i = 0; i < size; i++) { // copy data to avoid bad free
883 if (recvcounts[i] < 0 || (sendbuf != MPI_IN_PLACE && sendcounts[i] < 0))
884 return MPI_ERR_COUNT;
888 std::vector<int>* trace_sendcounts = new std::vector<int>;
889 std::vector<int>* trace_recvcounts = new std::vector<int>;
891 const void* real_sendbuf = sendbuf;
892 const int* real_sendcounts = sendcounts;
893 const int* real_senddisps = senddisps;
894 const MPI_Datatype* real_sendtypes = sendtypes;
895 unsigned long maxsize = 0;
896 for (int i = 0; i < size; i++) { // copy data to avoid bad free
897 if (recvtypes[i] == MPI_DATATYPE_NULL) {
898 delete trace_recvcounts;
899 delete trace_sendcounts;
902 recv_size += recvcounts[i] * recvtypes[i]->size();
903 trace_recvcounts->push_back(recvcounts[i] * recvtypes[i]->size());
904 if ((recvdisps[i] + (recvcounts[i] * recvtypes[i]->size())) > maxsize)
905 maxsize = recvdisps[i] + (recvcounts[i] * recvtypes[i]->size());
908 std::unique_ptr<unsigned char[]> tmp_sendbuf;
909 std::unique_ptr<int[]> tmp_sendcounts;
910 std::unique_ptr<int[]> tmp_senddisps;
911 std::unique_ptr<MPI_Datatype[]> tmp_sendtypes;
912 if (sendbuf == MPI_IN_PLACE) {
913 tmp_sendbuf.reset(new unsigned char[maxsize]);
914 real_sendbuf = memcpy(tmp_sendbuf.get(), recvbuf, maxsize);
915 tmp_sendcounts.reset(new int[size]);
916 std::copy(recvcounts, recvcounts + size, tmp_sendcounts.get());
917 real_sendcounts = tmp_sendcounts.get();
918 tmp_senddisps.reset(new int[size]);
919 std::copy(recvdisps, recvdisps + size, tmp_senddisps.get());
920 real_senddisps = tmp_senddisps.get();
921 tmp_sendtypes.reset(new MPI_Datatype[size]);
922 std::copy(recvtypes, recvtypes + size, tmp_sendtypes.get());
923 real_sendtypes = tmp_sendtypes.get();
926 for (int i = 0; i < size; i++) { // copy data to avoid bad free
927 send_size += real_sendcounts[i] * real_sendtypes[i]->size();
928 trace_sendcounts->push_back(real_sendcounts[i] * real_sendtypes[i]->size());
931 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallw" : "PMPI_Ialltoallw",
932 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
933 send_size, trace_sendcounts, recv_size, trace_recvcounts,
934 simgrid::smpi::Datatype::encode(real_sendtypes[0]),
935 simgrid::smpi::Datatype::encode(recvtypes[0])));
938 if (request == MPI_REQUEST_IGNORED)
939 retval = simgrid::smpi::Colls::alltoallw(real_sendbuf, real_sendcounts, real_senddisps, real_sendtypes, recvbuf,
940 recvcounts, recvdisps, recvtypes, comm);
942 retval = simgrid::smpi::Colls::ialltoallw(real_sendbuf, real_sendcounts, real_senddisps, real_sendtypes, recvbuf,
943 recvcounts, recvdisps, recvtypes, comm, request);
945 TRACE_smpi_comm_out(rank);