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 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){
23 if (inplacebuf == MPI_IN_PLACE) {
24 tmp_sendbuf.reset(new unsigned char[count * datatype->get_extent()]);
25 simgrid::smpi::Datatype::copy(otherbuf, count, datatype, tmp_sendbuf.get(), count, datatype);
26 return tmp_sendbuf.get();
31 /* PMPI User level calls */
33 int PMPI_Barrier(MPI_Comm comm)
35 return PMPI_Ibarrier(comm, MPI_REQUEST_IGNORED);
38 int PMPI_Ibarrier(MPI_Comm comm, MPI_Request *request)
40 if (comm == MPI_COMM_NULL)
42 if (request == nullptr)
46 int rank = simgrid::s4u::this_actor::get_pid();
47 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Barrier" : "PMPI_Ibarrier",
48 new simgrid::instr::NoOpTIData(request == MPI_REQUEST_IGNORED ? "barrier" : "ibarrier"));
49 if (request == MPI_REQUEST_IGNORED) {
50 simgrid::smpi::colls::barrier(comm);
51 // Barrier can be used to synchronize RMA calls. Finish all requests from comm before.
52 comm->finish_rma_calls();
54 simgrid::smpi::colls::ibarrier(comm, request);
56 TRACE_smpi_comm_out(rank);
61 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
63 return PMPI_Ibcast(buf, count, datatype, root, comm, MPI_REQUEST_IGNORED);
66 int PMPI_Ibcast(void *buf, int count, MPI_Datatype datatype,
67 int root, MPI_Comm comm, MPI_Request* request)
69 if (comm == MPI_COMM_NULL)
71 if (buf == nullptr && count > 0)
72 return MPI_ERR_BUFFER;
73 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
77 if (root < 0 || root >= comm->size())
79 if (request == nullptr)
83 int rank = simgrid::s4u::this_actor::get_pid();
84 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Bcast" : "PMPI_Ibcast",
85 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "bcast" : "ibcast", root, -1.0,
86 datatype->is_replayable() ? count : count * datatype->size(), -1,
87 simgrid::smpi::Datatype::encode(datatype), ""));
88 if (comm->size() > 1) {
89 if (request == MPI_REQUEST_IGNORED)
90 simgrid::smpi::colls::bcast(buf, count, datatype, root, comm);
92 simgrid::smpi::colls::ibcast(buf, count, datatype, root, comm, request);
94 if (request != MPI_REQUEST_IGNORED)
95 *request = MPI_REQUEST_NULL;
98 TRACE_smpi_comm_out(rank);
103 int PMPI_Gather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
104 int root, MPI_Comm comm){
105 return PMPI_Igather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
108 int PMPI_Igather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
109 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
111 if (comm == MPI_COMM_NULL)
113 if ((sendbuf == nullptr && sendcount > 0) || ((comm->rank() == root) && recvbuf == nullptr && recvcount > 0))
114 return MPI_ERR_BUFFER;
115 if (((sendbuf != MPI_IN_PLACE && sendcount > 0) && (sendtype == MPI_DATATYPE_NULL)) ||
116 ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL)))
118 if (((sendbuf != MPI_IN_PLACE) && (sendcount < 0)) || ((comm->rank() == root) && (recvcount < 0)))
119 return MPI_ERR_COUNT;
120 if (root < 0 || root >= comm->size())
122 if (request == nullptr)
126 const void* real_sendbuf = sendbuf;
127 int real_sendcount = sendcount;
128 MPI_Datatype real_sendtype = sendtype;
129 if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
131 real_sendtype = recvtype;
133 int rank = simgrid::s4u::this_actor::get_pid();
135 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gather" : "PMPI_Igather",
136 new simgrid::instr::CollTIData(
137 request == MPI_REQUEST_IGNORED ? "gather" : "igather", root, -1.0,
138 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
139 (comm->rank() != root || recvtype->is_replayable()) ? recvcount : recvcount * recvtype->size(),
140 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
141 if (request == MPI_REQUEST_IGNORED)
142 simgrid::smpi::colls::gather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm);
144 simgrid::smpi::colls::igather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm,
147 TRACE_smpi_comm_out(rank);
152 int PMPI_Gatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, const int *recvcounts, const int *displs,
153 MPI_Datatype recvtype, int root, MPI_Comm comm){
154 return PMPI_Igatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, root, comm, MPI_REQUEST_IGNORED);
157 int PMPI_Igatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
158 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
160 if (comm == MPI_COMM_NULL)
162 if ((sendbuf == nullptr && sendcount > 0) || ((comm->rank() == root) && recvbuf == nullptr))
163 return MPI_ERR_BUFFER;
164 if (((sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
165 ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL)))
167 if ((sendbuf != MPI_IN_PLACE) && (sendcount < 0))
168 return MPI_ERR_COUNT;
169 if ((comm->rank() == root) && (recvcounts == nullptr || displs == nullptr))
171 if (root < 0 || root >= comm->size())
173 if (request == nullptr)
176 for (int i = 0; i < comm->size(); i++) {
177 if ((comm->rank() == root) && (recvcounts[i] < 0))
178 return MPI_ERR_COUNT;
182 const void* real_sendbuf = sendbuf;
183 int real_sendcount = sendcount;
184 MPI_Datatype real_sendtype = sendtype;
185 if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
187 real_sendtype = recvtype;
190 int rank = simgrid::s4u::this_actor::get_pid();
191 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
193 std::vector<int>* trace_recvcounts = new std::vector<int>;
194 if (comm->rank() == root) {
195 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
196 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
199 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gatherv" : "PMPI_Igatherv",
200 new simgrid::instr::VarCollTIData(
201 request == MPI_REQUEST_IGNORED ? "gatherv" : "igatherv", root,
202 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
203 nullptr, dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(real_sendtype),
204 simgrid::smpi::Datatype::encode(recvtype)));
205 if (request == MPI_REQUEST_IGNORED)
206 simgrid::smpi::colls::gatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
209 simgrid::smpi::colls::igatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
210 root, comm, request);
212 TRACE_smpi_comm_out(rank);
217 int PMPI_Allgather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
218 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm){
219 return PMPI_Iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
222 int PMPI_Iallgather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
223 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
225 CHECK_ARGS(comm == MPI_COMM_NULL, MPI_ERR_COMM, "Iallgather: the communicator cannot be MPI_COMM_NULL");
226 CHECK_ARGS(recvbuf == nullptr && recvcount > 0, MPI_ERR_BUFFER, "Iallgather: param 4 recvbuf cannot be NULL");
227 CHECK_ARGS(sendbuf == nullptr && sendcount > 0, MPI_ERR_BUFFER,
228 "Iallgather: param 1 sendbuf cannot be NULL when sendcound > 0");
229 CHECK_ARGS((sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL), MPI_ERR_TYPE,
230 "Iallgather: param 3 sendtype cannot be MPI_DATATYPE_NULL when sendbuff is not MPI_IN_PLACE");
231 CHECK_ARGS(recvtype == MPI_DATATYPE_NULL, MPI_ERR_TYPE, "Iallgather: param 6 recvtype cannot be MPI_DATATYPE_NULL");
232 CHECK_ARGS(recvcount < 0, MPI_ERR_COUNT, "Iallgather: param 5 recvcount cannot be negative");
233 CHECK_ARGS((sendbuf != MPI_IN_PLACE) && (sendcount < 0), MPI_ERR_COUNT,
234 "Iallgather: param 2 sendcount cannot be negative when sendbuf is not MPI_IN_PLACE");
235 CHECK_ARGS(request == nullptr, MPI_ERR_ARG, "Iallgather: param 8 request cannot be NULL");
238 if (sendbuf == MPI_IN_PLACE) {
239 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * recvcount * comm->rank();
240 sendcount = recvcount;
243 int rank = simgrid::s4u::this_actor::get_pid();
245 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgather" : "PMPI_Iallggather",
246 new simgrid::instr::CollTIData(
247 request == MPI_REQUEST_IGNORED ? "allgather" : "iallgather", -1, -1.0,
248 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
249 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
250 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
251 if (request == MPI_REQUEST_IGNORED)
252 simgrid::smpi::colls::allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
254 simgrid::smpi::colls::iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, request);
256 TRACE_smpi_comm_out(rank);
261 int PMPI_Allgatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
262 void *recvbuf, const int *recvcounts, const int *displs, MPI_Datatype recvtype, MPI_Comm comm){
263 return PMPI_Iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, MPI_REQUEST_IGNORED);
266 int PMPI_Iallgatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
267 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
269 if (comm == MPI_COMM_NULL)
271 if (sendbuf == nullptr && sendcount > 0)
272 return MPI_ERR_BUFFER;
273 if (((sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) || (recvtype == MPI_DATATYPE_NULL))
275 if ((sendbuf != MPI_IN_PLACE) && (sendcount < 0))
276 return MPI_ERR_COUNT;
277 if (recvcounts == nullptr || displs == nullptr)
279 if (request == nullptr)
282 for (int i = 0; i < comm->size(); i++) {
283 if (recvcounts[i] < 0)
284 return MPI_ERR_COUNT;
285 else if (recvcounts[i] > 0 && recvbuf == nullptr)
286 return MPI_ERR_BUFFER;
290 if (sendbuf == MPI_IN_PLACE) {
291 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * displs[comm->rank()];
292 sendcount = recvcounts[comm->rank()];
295 int rank = simgrid::s4u::this_actor::get_pid();
296 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
298 std::vector<int>* trace_recvcounts = new std::vector<int>;
299 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
300 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
304 rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgatherv" : "PMPI_Iallgatherv",
305 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "allgatherv" : "iallgatherv", -1,
306 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(), nullptr,
307 dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
308 simgrid::smpi::Datatype::encode(recvtype)));
309 if (request == MPI_REQUEST_IGNORED)
310 simgrid::smpi::colls::allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
312 simgrid::smpi::colls::iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm,
315 TRACE_smpi_comm_out(rank);
320 int PMPI_Scatter(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
321 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
322 return PMPI_Iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
325 int PMPI_Iscatter(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
326 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
328 if (comm == MPI_COMM_NULL)
330 if (((comm->rank() == root) && (sendtype == MPI_DATATYPE_NULL || not sendtype->is_valid())) ||
331 ((recvbuf != MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL || not recvtype->is_valid())))
333 if (((comm->rank() == root) && (sendcount < 0)) || ((recvbuf != MPI_IN_PLACE) && (recvcount < 0)))
334 return MPI_ERR_COUNT;
335 if ((sendbuf == recvbuf) || ((comm->rank() == root) && sendcount > 0 && (sendbuf == nullptr)) ||
336 (recvcount > 0 && recvbuf == nullptr))
337 return MPI_ERR_BUFFER;
338 if (root < 0 || root >= comm->size())
340 if (request == nullptr)
344 if (recvbuf == MPI_IN_PLACE) {
346 recvcount = sendcount;
348 int rank = simgrid::s4u::this_actor::get_pid();
350 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatter" : "PMPI_Iscatter",
351 new simgrid::instr::CollTIData(
352 request == MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
353 (comm->rank() != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
354 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
355 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
356 if (request == MPI_REQUEST_IGNORED)
357 simgrid::smpi::colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
359 simgrid::smpi::colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
361 TRACE_smpi_comm_out(rank);
366 int PMPI_Scatterv(const void *sendbuf, const int *sendcounts, const int *displs,
367 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
368 return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
371 int PMPI_Iscatterv(const void* sendbuf, const int* sendcounts, const int* displs, MPI_Datatype sendtype, void* recvbuf, int recvcount,
372 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
374 CHECK_ARGS(comm == MPI_COMM_NULL, MPI_ERR_COMM, "Iscatterv: the communicator cannot be MPI_COMM_NULL");
375 CHECK_ARGS((comm->rank() == root) && (sendcounts == nullptr), MPI_ERR_ARG,
376 "Iscatterv: param 2 sendcounts cannot be NULL on the root rank");
377 CHECK_ARGS((comm->rank() == root) && (displs == nullptr), MPI_ERR_ARG,
378 "Iscatterv: param 3 displs cannot be NULL on the root rank");
379 CHECK_ARGS((comm->rank() == root) && (sendtype == MPI_DATATYPE_NULL), MPI_ERR_TYPE,
380 "Iscatterv: The sendtype cannot be NULL on the root rank");
381 CHECK_ARGS((recvbuf != MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL), MPI_ERR_TYPE,
382 "Iscatterv: the recvtype cannot be NULL when not receiving in place");
383 CHECK_ARGS(request == nullptr, MPI_ERR_ARG, "Iscatterv: param 10 request cannot be NULL");
384 CHECK_ARGS(recvbuf != MPI_IN_PLACE && recvcount < 0, MPI_ERR_COUNT,
385 "Iscatterv: When not receiving in place, the recvcount cannot be negative");
386 CHECK_ARGS(root < 0, MPI_ERR_ROOT, "Iscatterv: root cannot be negative");
387 CHECK_ARGS(root >= comm->size(), MPI_ERR_ROOT, "Iscatterv: root (=%d) is larger than communicator size (=%d)", root,
390 if (comm->rank() == root) {
391 if (recvbuf == MPI_IN_PLACE) {
393 recvcount = sendcounts[comm->rank()];
395 for (int i = 0; i < comm->size(); i++)
396 CHECK_ARGS(sendcounts[i] < 0, MPI_ERR_COUNT, "Iscatterv: sendcounts[%d]=%d but this cannot be negative", i,
402 int rank = simgrid::s4u::this_actor::get_pid();
403 int dt_size_send = sendtype->is_replayable() ? 1 : sendtype->size();
405 std::vector<int>* trace_sendcounts = new std::vector<int>;
406 if (comm->rank() == root) {
407 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
408 trace_sendcounts->push_back(sendcounts[i] * dt_size_send);
412 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatterv" : "PMPI_Iscatterv",
413 new simgrid::instr::VarCollTIData(
414 request == MPI_REQUEST_IGNORED ? "scatterv" : "iscatterv", root, dt_size_send,
415 trace_sendcounts, recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
416 nullptr, simgrid::smpi::Datatype::encode(sendtype),
417 simgrid::smpi::Datatype::encode(recvtype)));
418 if (request == MPI_REQUEST_IGNORED)
419 simgrid::smpi::colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
421 simgrid::smpi::colls::iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm,
424 TRACE_smpi_comm_out(rank);
429 int PMPI_Reduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
431 return PMPI_Ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, MPI_REQUEST_IGNORED);
434 int PMPI_Ireduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm, MPI_Request* request)
436 if (comm == MPI_COMM_NULL)
438 if ((sendbuf == nullptr && count > 0) || ((comm->rank() == root) && recvbuf == nullptr))
439 return MPI_ERR_BUFFER;
440 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
442 if (op == MPI_OP_NULL)
444 if (request == nullptr)
446 if (root < 0 || root >= comm->size())
449 return MPI_ERR_COUNT;
452 int rank = simgrid::s4u::this_actor::get_pid();
454 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce" : "PMPI_Ireduce",
455 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "reduce" : "ireduce", root, 0,
456 datatype->is_replayable() ? count : count * datatype->size(), -1,
457 simgrid::smpi::Datatype::encode(datatype), ""));
458 if (request == MPI_REQUEST_IGNORED)
459 simgrid::smpi::colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
461 simgrid::smpi::colls::ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, request);
463 TRACE_smpi_comm_out(rank);
468 int PMPI_Reduce_local(const void* inbuf, void* inoutbuf, int count, MPI_Datatype datatype, MPI_Op op)
470 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
472 if (op == MPI_OP_NULL)
475 return MPI_ERR_COUNT;
478 op->apply(inbuf, inoutbuf, &count, datatype);
483 int PMPI_Allreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
485 return PMPI_Iallreduce(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
488 int PMPI_Iallreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
490 if (comm == MPI_COMM_NULL)
492 if ((sendbuf == nullptr && count > 0) || (recvbuf == nullptr))
493 return MPI_ERR_BUFFER;
494 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
497 return MPI_ERR_COUNT;
498 if (op == MPI_OP_NULL)
500 if (request == nullptr)
504 std::unique_ptr<unsigned char[]> tmp_sendbuf;
505 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
507 int rank = simgrid::s4u::this_actor::get_pid();
509 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allreduce" : "PMPI_Iallreduce",
510 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "allreduce" : "iallreduce", -1, 0,
511 datatype->is_replayable() ? count : count * datatype->size(), -1,
512 simgrid::smpi::Datatype::encode(datatype), ""));
514 if (request == MPI_REQUEST_IGNORED)
515 simgrid::smpi::colls::allreduce(real_sendbuf, recvbuf, count, datatype, op, comm);
517 simgrid::smpi::colls::iallreduce(real_sendbuf, recvbuf, count, datatype, op, comm, request);
519 TRACE_smpi_comm_out(rank);
524 int PMPI_Scan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
526 return PMPI_Iscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
529 int PMPI_Iscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request)
531 if (comm == MPI_COMM_NULL)
533 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid())
535 if (op == MPI_OP_NULL)
537 if (request == nullptr)
540 return MPI_ERR_COUNT;
541 if (sendbuf == nullptr || recvbuf == nullptr)
542 return MPI_ERR_BUFFER;
545 int rank = simgrid::s4u::this_actor::get_pid();
546 std::unique_ptr<unsigned char[]> tmp_sendbuf;
547 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);;
549 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scan" : "PMPI_Iscan",
550 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "scan" : "iscan", -1,
551 datatype->is_replayable() ? count : count * datatype->size(),
552 simgrid::smpi::Datatype::encode(datatype)));
555 if (request == MPI_REQUEST_IGNORED)
556 retval = simgrid::smpi::colls::scan(real_sendbuf, recvbuf, count, datatype, op, comm);
558 retval = simgrid::smpi::colls::iscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
560 TRACE_smpi_comm_out(rank);
565 int PMPI_Exscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
567 return PMPI_Iexscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
570 int PMPI_Iexscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request){
571 if (comm == MPI_COMM_NULL)
573 if (not datatype->is_valid())
575 if (op == MPI_OP_NULL)
577 if (request == nullptr)
580 return MPI_ERR_COUNT;
581 if (sendbuf == nullptr || recvbuf == nullptr)
582 return MPI_ERR_BUFFER;
585 int rank = simgrid::s4u::this_actor::get_pid();
586 std::unique_ptr<unsigned char[]> tmp_sendbuf;
587 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);;
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];
640 std::unique_ptr<unsigned char[]> tmp_sendbuf;
641 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, totalcount, datatype);;
643 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter" : "PMPI_Ireduce_scatter",
644 new simgrid::instr::VarCollTIData(
645 request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, dt_send_size, nullptr,
646 -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
648 if (request == MPI_REQUEST_IGNORED)
649 simgrid::smpi::colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm);
651 simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm, request);
653 TRACE_smpi_comm_out(rank);
658 int PMPI_Reduce_scatter_block(const void *sendbuf, void *recvbuf, int recvcount,
659 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
661 return PMPI_Ireduce_scatter_block(sendbuf, recvbuf, recvcount, datatype, op, comm, MPI_REQUEST_IGNORED);
664 int PMPI_Ireduce_scatter_block(const void* sendbuf, void* recvbuf, int recvcount, MPI_Datatype datatype, MPI_Op op,
665 MPI_Comm comm, MPI_Request* request)
667 if (comm == MPI_COMM_NULL)
669 if (not datatype->is_valid())
671 if (op == MPI_OP_NULL)
675 if (request == nullptr)
679 int count = comm->size();
681 int rank = simgrid::s4u::this_actor::get_pid();
682 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
683 std::vector<int>* trace_recvcounts = new std::vector<int>(recvcount * dt_send_size); // copy data to avoid bad free
684 std::unique_ptr<unsigned char[]> tmp_sendbuf;
685 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, recvcount * count, datatype);
688 rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter_block" : "PMPI_Ireduce_scatter_block",
689 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, 0,
690 nullptr, -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
692 int* recvcounts = new int[count];
693 for (int i = 0; i < count; i++)
694 recvcounts[i] = recvcount;
695 if (request == MPI_REQUEST_IGNORED)
696 simgrid::smpi::colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm);
698 simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm, request);
701 TRACE_smpi_comm_out(rank);
706 int PMPI_Alltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
707 MPI_Datatype recvtype, MPI_Comm comm){
708 return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
711 int PMPI_Ialltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
712 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
714 if (comm == MPI_COMM_NULL)
716 if ((sendbuf == nullptr && sendcount > 0) || (recvbuf == nullptr && recvcount > 0))
717 return MPI_ERR_BUFFER;
718 if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL)
720 if ((sendbuf != MPI_IN_PLACE && sendcount < 0) || recvcount < 0)
721 return MPI_ERR_COUNT;
722 if (request == nullptr)
726 int rank = simgrid::s4u::this_actor::get_pid();
727 int real_sendcount = sendcount;
728 MPI_Datatype real_sendtype = sendtype;
730 std::unique_ptr<unsigned char[]> tmp_sendbuf;
731 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, recvcount * comm->size(), recvtype);
733 if (sendbuf == MPI_IN_PLACE) {
734 real_sendcount = recvcount;
735 real_sendtype = recvtype;
738 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoall" : "PMPI_Ialltoall",
739 new simgrid::instr::CollTIData(
740 request == MPI_REQUEST_IGNORED ? "alltoall" : "ialltoall", -1, -1.0,
741 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
742 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
743 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
745 if (request == MPI_REQUEST_IGNORED)
747 simgrid::smpi::colls::alltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, comm);
749 retval = simgrid::smpi::colls::ialltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype,
752 TRACE_smpi_comm_out(rank);
757 int PMPI_Alltoallv(const void* sendbuf, const int* sendcounts, const int* senddisps, MPI_Datatype sendtype, void* recvbuf,
758 const int* recvcounts, const int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm)
760 return PMPI_Ialltoallv(sendbuf, sendcounts, senddisps, sendtype, recvbuf, recvcounts, recvdisps, recvtype, comm, MPI_REQUEST_IGNORED);
763 int PMPI_Ialltoallv(const void* sendbuf, const int* sendcounts, const int* senddisps, MPI_Datatype sendtype, void* recvbuf,
764 const int* recvcounts, const int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
766 if (comm == MPI_COMM_NULL)
768 if (sendbuf == nullptr || recvbuf == nullptr)
769 return MPI_ERR_BUFFER;
770 if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL)
772 if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
773 recvdisps == nullptr)
775 if (request == nullptr)
778 int rank = simgrid::s4u::this_actor::get_pid();
779 int size = comm->size();
780 for (int i = 0; i < size; i++) {
781 if (recvcounts[i] < 0 || (sendbuf != MPI_IN_PLACE && sendcounts[i] < 0))
782 return MPI_ERR_COUNT;
788 std::vector<int>* trace_sendcounts = new std::vector<int>;
789 std::vector<int>* trace_recvcounts = new std::vector<int>;
790 int dt_size_recv = recvtype->size();
792 const int* real_sendcounts = sendcounts;
793 const int* real_senddisps = senddisps;
794 MPI_Datatype real_sendtype = sendtype;
796 for (int i = 0; i < size; i++) { // copy data to avoid bad free
797 recv_size += recvcounts[i] * dt_size_recv;
798 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
799 if (((recvdisps[i] + recvcounts[i]) * dt_size_recv) > maxsize)
800 maxsize = (recvdisps[i] + recvcounts[i]) * dt_size_recv;
803 std::unique_ptr<unsigned char[]> tmp_sendbuf;
804 std::unique_ptr<int[]> tmp_sendcounts;
805 std::unique_ptr<int[]> tmp_senddisps;
806 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
807 if (sendbuf == MPI_IN_PLACE) {
808 tmp_sendcounts.reset(new int[size]);
809 std::copy(recvcounts, recvcounts + size, tmp_sendcounts.get());
810 real_sendcounts = tmp_sendcounts.get();
811 tmp_senddisps.reset(new int[size]);
812 std::copy(recvdisps, recvdisps + size, tmp_senddisps.get());
813 real_senddisps = tmp_senddisps.get();
814 real_sendtype = recvtype;
817 int dt_size_send = real_sendtype->size();
819 for (int i = 0; i < size; i++) { // copy data to avoid bad free
820 send_size += real_sendcounts[i] * dt_size_send;
821 trace_sendcounts->push_back(real_sendcounts[i] * dt_size_send);
824 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallv" : "PMPI_Ialltoallv",
825 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
826 send_size, trace_sendcounts, recv_size, trace_recvcounts,
827 simgrid::smpi::Datatype::encode(real_sendtype),
828 simgrid::smpi::Datatype::encode(recvtype)));
831 if (request == MPI_REQUEST_IGNORED)
832 retval = simgrid::smpi::colls::alltoallv(real_sendbuf, real_sendcounts, real_senddisps, real_sendtype, recvbuf,
833 recvcounts, recvdisps, recvtype, comm);
835 retval = simgrid::smpi::colls::ialltoallv(real_sendbuf, real_sendcounts, real_senddisps, real_sendtype, recvbuf,
836 recvcounts, recvdisps, recvtype, comm, request);
838 TRACE_smpi_comm_out(rank);
843 int PMPI_Alltoallw(const void* sendbuf, const int* sendcounts, const int* senddisps, const MPI_Datatype* sendtypes, void* recvbuf,
844 const int* recvcounts, const int* recvdisps, const MPI_Datatype* recvtypes, MPI_Comm comm)
846 return PMPI_Ialltoallw(sendbuf, sendcounts, senddisps, sendtypes, recvbuf, recvcounts, recvdisps, recvtypes, comm, MPI_REQUEST_IGNORED);
849 int PMPI_Ialltoallw(const void* sendbuf, const int* sendcounts, const int* senddisps, const MPI_Datatype* sendtypes, void* recvbuf,
850 const int* recvcounts, const int* recvdisps, const MPI_Datatype* recvtypes, MPI_Comm comm, MPI_Request* request)
852 if (comm == MPI_COMM_NULL)
854 if (sendbuf == nullptr || recvbuf == nullptr)
855 return MPI_ERR_BUFFER;
856 if ((sendbuf != MPI_IN_PLACE && sendtypes == nullptr) || recvtypes == nullptr)
858 if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
859 recvdisps == nullptr)
861 if (request == nullptr)
865 int rank = simgrid::s4u::this_actor::get_pid();
866 int size = comm->size();
867 for (int i = 0; i < size; i++) {
868 if (recvcounts[i] < 0 || (sendbuf != MPI_IN_PLACE && sendcounts[i] < 0))
869 return MPI_ERR_COUNT;
873 std::vector<int>* trace_sendcounts = new std::vector<int>;
874 std::vector<int>* trace_recvcounts = new std::vector<int>;
876 const int* real_sendcounts = sendcounts;
877 const int* real_senddisps = senddisps;
878 const MPI_Datatype* real_sendtypes = sendtypes;
879 unsigned long maxsize = 0;
880 for (int i = 0; i < size; i++) { // copy data to avoid bad free
881 if (recvtypes[i] == MPI_DATATYPE_NULL) {
882 delete trace_recvcounts;
883 delete trace_sendcounts;
886 recv_size += recvcounts[i] * recvtypes[i]->size();
887 trace_recvcounts->push_back(recvcounts[i] * recvtypes[i]->size());
888 if ((recvdisps[i] + (recvcounts[i] * recvtypes[i]->size())) > maxsize)
889 maxsize = recvdisps[i] + (recvcounts[i] * recvtypes[i]->size());
892 std::unique_ptr<unsigned char[]> tmp_sendbuf;
893 std::unique_ptr<int[]> tmp_sendcounts;
894 std::unique_ptr<int[]> tmp_senddisps;
895 std::unique_ptr<MPI_Datatype[]> tmp_sendtypes;
896 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
897 if (sendbuf == MPI_IN_PLACE) {
898 tmp_sendcounts.reset(new int[size]);
899 std::copy(recvcounts, recvcounts + size, tmp_sendcounts.get());
900 real_sendcounts = tmp_sendcounts.get();
901 tmp_senddisps.reset(new int[size]);
902 std::copy(recvdisps, recvdisps + size, tmp_senddisps.get());
903 real_senddisps = tmp_senddisps.get();
904 tmp_sendtypes.reset(new MPI_Datatype[size]);
905 std::copy(recvtypes, recvtypes + size, tmp_sendtypes.get());
906 real_sendtypes = tmp_sendtypes.get();
909 for (int i = 0; i < size; i++) { // copy data to avoid bad free
910 send_size += real_sendcounts[i] * real_sendtypes[i]->size();
911 trace_sendcounts->push_back(real_sendcounts[i] * real_sendtypes[i]->size());
914 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallw" : "PMPI_Ialltoallw",
915 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
916 send_size, trace_sendcounts, recv_size, trace_recvcounts,
917 simgrid::smpi::Datatype::encode(real_sendtypes[0]),
918 simgrid::smpi::Datatype::encode(recvtypes[0])));
921 if (request == MPI_REQUEST_IGNORED)
922 retval = simgrid::smpi::colls::alltoallw(real_sendbuf, real_sendcounts, real_senddisps, real_sendtypes, recvbuf,
923 recvcounts, recvdisps, recvtypes, comm);
925 retval = simgrid::smpi::colls::ialltoallw(real_sendbuf, real_sendcounts, real_senddisps, real_sendtypes, recvbuf,
926 recvcounts, recvdisps, recvtypes, comm, request);
928 TRACE_smpi_comm_out(rank);