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 #define CHECK_COMM(num)\
23 CHECK_ARGS(comm == MPI_COMM_NULL, MPI_ERR_COMM,\
24 "%s: param %d communicator cannot be MPI_COMM_NULL", __func__, num);
25 #define CHECK_REQUEST(num)\
26 CHECK_ARGS(request == nullptr, MPI_ERR_ARG,\
27 "%s: param %d request cannot be NULL",__func__, num);
28 #define CHECK_BUFFER(num,buf,count)\
29 CHECK_ARGS(buf == nullptr && count > 0, MPI_ERR_BUFFER,\
30 "%s: param %d %s cannot be NULL if %s > 0",__func__, num, #buf, #count);
31 #define CHECK_COUNT(num,count)\
32 CHECK_ARGS(count < 0, MPI_ERR_COUNT,\
33 "%s: param %d %s cannot be negative", __func__, num, #count);
34 #define CHECK_TYPE(num, datatype)\
35 CHECK_ARGS((datatype == MPI_DATATYPE_NULL|| not datatype->is_valid()), MPI_ERR_TYPE,\
36 "%s: param %d %s cannot be MPI_DATATYPE_NULL or invalid", __func__, num, #datatype);
37 #define CHECK_OP(num)\
38 CHECK_ARGS(op == MPI_OP_NULL, MPI_ERR_OP,\
39 "%s: param %d op cannot be MPI_OP_NULL or invalid", __func__, num);
40 #define CHECK_ROOT(num)\
41 CHECK_ARGS((root < 0 || root >= comm->size()), MPI_ERR_ROOT,\
42 "%s: param %d root (=%d) cannot be negative or larger than communicator size (=%d)", __func__, num, root,\
44 #define CHECK_NULL(num,err,buf)\
45 CHECK_ARGS(buf == nullptr, err,\
46 "%s: param %d %s cannot be NULL", __func__, num, #buf);
48 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){
49 if (inplacebuf == MPI_IN_PLACE) {
50 tmp_sendbuf.reset(new unsigned char[count * datatype->get_extent()]);
51 simgrid::smpi::Datatype::copy(otherbuf, count, datatype, tmp_sendbuf.get(), count, datatype);
52 return tmp_sendbuf.get();
57 /* PMPI User level calls */
59 int PMPI_Barrier(MPI_Comm comm)
61 return PMPI_Ibarrier(comm, MPI_REQUEST_IGNORED);
64 int PMPI_Ibarrier(MPI_Comm comm, MPI_Request *request)
70 int rank = simgrid::s4u::this_actor::get_pid();
71 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Barrier" : "PMPI_Ibarrier",
72 new simgrid::instr::NoOpTIData(request == MPI_REQUEST_IGNORED ? "barrier" : "ibarrier"));
73 if (request == MPI_REQUEST_IGNORED) {
74 simgrid::smpi::colls::barrier(comm);
75 // Barrier can be used to synchronize RMA calls. Finish all requests from comm before.
76 comm->finish_rma_calls();
78 simgrid::smpi::colls::ibarrier(comm, request);
80 TRACE_smpi_comm_out(rank);
85 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
87 return PMPI_Ibcast(buf, count, datatype, root, comm, MPI_REQUEST_IGNORED);
90 int PMPI_Ibcast(void *buf, int count, MPI_Datatype datatype,
91 int root, MPI_Comm comm, MPI_Request* request)
94 CHECK_BUFFER(1, buf, count)
96 CHECK_TYPE(3, datatype)
101 int rank = simgrid::s4u::this_actor::get_pid();
102 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Bcast" : "PMPI_Ibcast",
103 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "bcast" : "ibcast", root, -1.0,
104 datatype->is_replayable() ? count : count * datatype->size(), -1,
105 simgrid::smpi::Datatype::encode(datatype), ""));
106 if (comm->size() > 1) {
107 if (request == MPI_REQUEST_IGNORED)
108 simgrid::smpi::colls::bcast(buf, count, datatype, root, comm);
110 simgrid::smpi::colls::ibcast(buf, count, datatype, root, comm, request);
112 if (request != MPI_REQUEST_IGNORED)
113 *request = MPI_REQUEST_NULL;
116 TRACE_smpi_comm_out(rank);
121 int PMPI_Gather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
122 int root, MPI_Comm comm){
123 return PMPI_Igather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
126 int PMPI_Igather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
127 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
130 if(sendbuf != MPI_IN_PLACE){
131 CHECK_BUFFER(1,sendbuf, sendcount)
132 CHECK_COUNT(2, sendcount)
133 CHECK_TYPE(3, sendtype)
135 if(comm->rank() == root){
136 CHECK_TYPE(6, recvtype)
137 CHECK_COUNT(5, recvcount)
138 CHECK_BUFFER(4, recvbuf, recvcount)
144 const void* real_sendbuf = sendbuf;
145 int real_sendcount = sendcount;
146 MPI_Datatype real_sendtype = sendtype;
147 if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
149 real_sendtype = recvtype;
151 int rank = simgrid::s4u::this_actor::get_pid();
153 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gather" : "PMPI_Igather",
154 new simgrid::instr::CollTIData(
155 request == MPI_REQUEST_IGNORED ? "gather" : "igather", root, -1.0,
156 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
157 (comm->rank() != root || recvtype->is_replayable()) ? recvcount : recvcount * recvtype->size(),
158 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
159 if (request == MPI_REQUEST_IGNORED)
160 simgrid::smpi::colls::gather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm);
162 simgrid::smpi::colls::igather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm,
165 TRACE_smpi_comm_out(rank);
170 int PMPI_Gatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, const int *recvcounts, const int *displs,
171 MPI_Datatype recvtype, int root, MPI_Comm comm){
172 return PMPI_Igatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, root, comm, MPI_REQUEST_IGNORED);
175 int PMPI_Igatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
176 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
179 CHECK_BUFFER(1, sendbuf, sendcount)
180 if(sendbuf != MPI_IN_PLACE){
181 CHECK_TYPE(3, sendtype)
182 CHECK_COUNT(2, sendcount)
184 if(comm->rank() == root){
185 CHECK_TYPE(6, recvtype)
186 CHECK_NULL(5, MPI_ERR_COUNT, recvcounts)
187 CHECK_NULL(6, MPI_ERR_ARG, displs)
192 if (comm->rank() == root){
193 for (int i = 0; i < comm->size(); i++) {
194 CHECK_COUNT(5, recvcounts[i])
195 CHECK_BUFFER(4,recvbuf,recvcounts[i])
200 const void* real_sendbuf = sendbuf;
201 int real_sendcount = sendcount;
202 MPI_Datatype real_sendtype = sendtype;
203 if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
205 real_sendtype = recvtype;
208 int rank = simgrid::s4u::this_actor::get_pid();
209 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
211 std::vector<int>* trace_recvcounts = new std::vector<int>;
212 if (comm->rank() == root) {
213 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
214 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
217 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gatherv" : "PMPI_Igatherv",
218 new simgrid::instr::VarCollTIData(
219 request == MPI_REQUEST_IGNORED ? "gatherv" : "igatherv", root,
220 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
221 nullptr, dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(real_sendtype),
222 simgrid::smpi::Datatype::encode(recvtype)));
223 if (request == MPI_REQUEST_IGNORED)
224 simgrid::smpi::colls::gatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
227 simgrid::smpi::colls::igatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
228 root, comm, request);
230 TRACE_smpi_comm_out(rank);
235 int PMPI_Allgather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
236 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm){
237 return PMPI_Iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
240 int PMPI_Iallgather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
241 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
244 CHECK_BUFFER(1, sendbuf, sendcount)
245 CHECK_BUFFER(4, recvbuf, recvcount)
246 if(sendbuf != MPI_IN_PLACE){
247 CHECK_COUNT(2, sendcount)
248 CHECK_TYPE(3, sendtype)
250 CHECK_TYPE(6, recvtype)
251 CHECK_COUNT(5, recvcount)
255 if (sendbuf == MPI_IN_PLACE) {
256 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * recvcount * comm->rank();
257 sendcount = recvcount;
260 int rank = simgrid::s4u::this_actor::get_pid();
262 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgather" : "PMPI_Iallggather",
263 new simgrid::instr::CollTIData(
264 request == MPI_REQUEST_IGNORED ? "allgather" : "iallgather", -1, -1.0,
265 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
266 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
267 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
268 if (request == MPI_REQUEST_IGNORED)
269 simgrid::smpi::colls::allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
271 simgrid::smpi::colls::iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, request);
273 TRACE_smpi_comm_out(rank);
278 int PMPI_Allgatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
279 void *recvbuf, const int *recvcounts, const int *displs, MPI_Datatype recvtype, MPI_Comm comm){
280 return PMPI_Iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, MPI_REQUEST_IGNORED);
283 int PMPI_Iallgatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
284 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
287 CHECK_BUFFER(1, sendbuf, sendcount)
288 if(sendbuf != MPI_IN_PLACE)
289 CHECK_TYPE(3, sendtype)
290 CHECK_TYPE(6, recvtype)
291 CHECK_NULL(5, MPI_ERR_COUNT, recvcounts)
292 CHECK_NULL(6, MPI_ERR_ARG, displs)
293 if(sendbuf != MPI_IN_PLACE)
294 CHECK_COUNT(2, sendcount)
296 for (int i = 0; i < comm->size(); i++) {
297 CHECK_COUNT(5, recvcounts[i])
298 CHECK_BUFFER(4, recvbuf, recvcounts[i])
302 if (sendbuf == MPI_IN_PLACE) {
303 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * displs[comm->rank()];
304 sendcount = recvcounts[comm->rank()];
307 int rank = simgrid::s4u::this_actor::get_pid();
308 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
310 std::vector<int>* trace_recvcounts = new std::vector<int>;
311 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
312 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
316 rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgatherv" : "PMPI_Iallgatherv",
317 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "allgatherv" : "iallgatherv", -1,
318 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(), nullptr,
319 dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
320 simgrid::smpi::Datatype::encode(recvtype)));
321 if (request == MPI_REQUEST_IGNORED)
322 simgrid::smpi::colls::allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
324 simgrid::smpi::colls::iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm,
327 TRACE_smpi_comm_out(rank);
332 int PMPI_Scatter(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
333 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
334 return PMPI_Iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
337 int PMPI_Iscatter(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
338 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
341 if(comm->rank() == root){
342 CHECK_BUFFER(1, sendbuf, sendcount)
343 CHECK_COUNT(2, sendcount)
344 CHECK_TYPE(3, sendtype)
346 if(recvbuf != MPI_IN_PLACE){
347 CHECK_BUFFER(4, recvbuf, recvcount)
348 CHECK_COUNT(5, recvcount)
349 CHECK_TYPE(6, recvtype)
355 if (recvbuf == MPI_IN_PLACE) {
357 recvcount = sendcount;
359 int rank = simgrid::s4u::this_actor::get_pid();
361 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatter" : "PMPI_Iscatter",
362 new simgrid::instr::CollTIData(
363 request == MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
364 (comm->rank() != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
365 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
366 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
367 if (request == MPI_REQUEST_IGNORED)
368 simgrid::smpi::colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
370 simgrid::smpi::colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
372 TRACE_smpi_comm_out(rank);
377 int PMPI_Scatterv(const void *sendbuf, const int *sendcounts, const int *displs,
378 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
379 return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
382 int PMPI_Iscatterv(const void* sendbuf, const int* sendcounts, const int* displs, MPI_Datatype sendtype, void* recvbuf, int recvcount,
383 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
386 if(recvbuf != MPI_IN_PLACE){
387 CHECK_BUFFER(4, recvbuf, recvcount)
388 CHECK_COUNT(5, recvcount)
389 CHECK_TYPE(7, recvtype)
393 if (comm->rank() == root) {
394 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
395 CHECK_NULL(3, MPI_ERR_ARG, displs)
396 CHECK_TYPE(4, sendtype)
397 for (int i = 0; i < comm->size(); i++){
398 CHECK_BUFFER(1, sendbuf, sendcounts[i])
399 CHECK_COUNT(2, sendcounts[i])
401 if (recvbuf == MPI_IN_PLACE) {
403 recvcount = sendcounts[comm->rank()];
409 int rank = simgrid::s4u::this_actor::get_pid();
410 int dt_size_send = sendtype->is_replayable() ? 1 : sendtype->size();
412 std::vector<int>* trace_sendcounts = new std::vector<int>;
413 if (comm->rank() == root) {
414 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
415 trace_sendcounts->push_back(sendcounts[i] * dt_size_send);
419 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatterv" : "PMPI_Iscatterv",
420 new simgrid::instr::VarCollTIData(
421 request == MPI_REQUEST_IGNORED ? "scatterv" : "iscatterv", root, dt_size_send,
422 trace_sendcounts, recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
423 nullptr, simgrid::smpi::Datatype::encode(sendtype),
424 simgrid::smpi::Datatype::encode(recvtype)));
425 if (request == MPI_REQUEST_IGNORED)
426 simgrid::smpi::colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
428 simgrid::smpi::colls::iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm,
431 TRACE_smpi_comm_out(rank);
436 int PMPI_Reduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
438 return PMPI_Ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, MPI_REQUEST_IGNORED);
441 int PMPI_Ireduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm, MPI_Request* request)
444 CHECK_BUFFER(1, sendbuf, count)
445 if(comm->rank() == root)
446 CHECK_BUFFER(5, recvbuf, count)
447 CHECK_TYPE(4, datatype)
448 CHECK_COUNT(3, count)
454 int rank = simgrid::s4u::this_actor::get_pid();
456 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce" : "PMPI_Ireduce",
457 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "reduce" : "ireduce", root, 0,
458 datatype->is_replayable() ? count : count * datatype->size(), -1,
459 simgrid::smpi::Datatype::encode(datatype), ""));
460 if (request == MPI_REQUEST_IGNORED)
461 simgrid::smpi::colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
463 simgrid::smpi::colls::ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, request);
465 TRACE_smpi_comm_out(rank);
470 int PMPI_Reduce_local(const void* inbuf, void* inoutbuf, int count, MPI_Datatype datatype, MPI_Op op)
472 CHECK_BUFFER(1, inbuf, count)
473 CHECK_BUFFER(2, inoutbuf, count)
474 CHECK_TYPE(4, datatype)
475 CHECK_COUNT(3, count)
479 op->apply(inbuf, inoutbuf, &count, datatype);
484 int PMPI_Allreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
486 return PMPI_Iallreduce(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
489 int PMPI_Iallreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
493 CHECK_BUFFER(1, sendbuf, count)
494 CHECK_BUFFER(2, recvbuf, count)
495 CHECK_TYPE(4, datatype)
496 CHECK_COUNT(3, count)
501 std::unique_ptr<unsigned char[]> tmp_sendbuf;
502 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
504 int rank = simgrid::s4u::this_actor::get_pid();
506 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allreduce" : "PMPI_Iallreduce",
507 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "allreduce" : "iallreduce", -1, 0,
508 datatype->is_replayable() ? count : count * datatype->size(), -1,
509 simgrid::smpi::Datatype::encode(datatype), ""));
511 if (request == MPI_REQUEST_IGNORED)
512 simgrid::smpi::colls::allreduce(real_sendbuf, recvbuf, count, datatype, op, comm);
514 simgrid::smpi::colls::iallreduce(real_sendbuf, recvbuf, count, datatype, op, comm, request);
516 TRACE_smpi_comm_out(rank);
521 int PMPI_Scan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
523 return PMPI_Iscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
526 int PMPI_Iscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request)
529 CHECK_BUFFER(1,sendbuf,count)
530 CHECK_BUFFER(2,recvbuf,count)
531 CHECK_TYPE(4, datatype)
532 CHECK_COUNT(3, count)
537 int rank = simgrid::s4u::this_actor::get_pid();
538 std::unique_ptr<unsigned char[]> tmp_sendbuf;
539 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
541 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scan" : "PMPI_Iscan",
542 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "scan" : "iscan", -1,
543 datatype->is_replayable() ? count : count * datatype->size(),
544 simgrid::smpi::Datatype::encode(datatype)));
547 if (request == MPI_REQUEST_IGNORED)
548 retval = simgrid::smpi::colls::scan(real_sendbuf, recvbuf, count, datatype, op, comm);
550 retval = simgrid::smpi::colls::iscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
552 TRACE_smpi_comm_out(rank);
557 int PMPI_Exscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
559 return PMPI_Iexscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
562 int PMPI_Iexscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request){
564 CHECK_BUFFER(1, sendbuf, count)
565 CHECK_BUFFER(2, recvbuf, count)
566 CHECK_TYPE(4, datatype)
567 CHECK_COUNT(3, count)
572 int rank = simgrid::s4u::this_actor::get_pid();
573 std::unique_ptr<unsigned char[]> tmp_sendbuf;
574 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
576 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Exscan" : "PMPI_Iexscan",
577 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "exscan" : "iexscan", -1,
578 datatype->is_replayable() ? count : count * datatype->size(),
579 simgrid::smpi::Datatype::encode(datatype)));
582 if (request == MPI_REQUEST_IGNORED)
583 retval = simgrid::smpi::colls::exscan(real_sendbuf, recvbuf, count, datatype, op, comm);
585 retval = simgrid::smpi::colls::iexscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
587 TRACE_smpi_comm_out(rank);
592 int PMPI_Reduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
594 return PMPI_Ireduce_scatter(sendbuf, recvbuf, recvcounts, datatype, op, comm, MPI_REQUEST_IGNORED);
597 int PMPI_Ireduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
600 CHECK_TYPE(4, datatype)
601 CHECK_NULL(3, MPI_ERR_COUNT, recvcounts)
604 for (int i = 0; i < comm->size(); i++) {
605 CHECK_COUNT(3, recvcounts[i])
606 CHECK_BUFFER(1, sendbuf, recvcounts[i])
607 CHECK_BUFFER(2, recvbuf, recvcounts[i])
611 int rank = simgrid::s4u::this_actor::get_pid();
612 std::vector<int>* trace_recvcounts = new std::vector<int>;
613 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
616 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
617 trace_recvcounts->push_back(recvcounts[i] * dt_send_size);
618 totalcount += recvcounts[i];
620 std::unique_ptr<unsigned char[]> tmp_sendbuf;
621 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, totalcount, datatype);
623 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter" : "PMPI_Ireduce_scatter",
624 new simgrid::instr::VarCollTIData(
625 request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, dt_send_size, nullptr,
626 -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
628 if (request == MPI_REQUEST_IGNORED)
629 simgrid::smpi::colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm);
631 simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm, request);
633 TRACE_smpi_comm_out(rank);
638 int PMPI_Reduce_scatter_block(const void *sendbuf, void *recvbuf, int recvcount,
639 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
641 return PMPI_Ireduce_scatter_block(sendbuf, recvbuf, recvcount, datatype, op, comm, MPI_REQUEST_IGNORED);
644 int PMPI_Ireduce_scatter_block(const void* sendbuf, void* recvbuf, int recvcount, MPI_Datatype datatype, MPI_Op op,
645 MPI_Comm comm, MPI_Request* request)
648 CHECK_BUFFER(1, sendbuf, recvcount)
649 CHECK_BUFFER(2, recvbuf, recvcount)
650 CHECK_TYPE(4, datatype)
651 CHECK_COUNT(3, recvcount)
656 int count = comm->size();
658 int rank = simgrid::s4u::this_actor::get_pid();
659 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
660 std::vector<int>* trace_recvcounts = new std::vector<int>(recvcount * dt_send_size); // copy data to avoid bad free
661 std::unique_ptr<unsigned char[]> tmp_sendbuf;
662 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, recvcount * count, datatype);
665 rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter_block" : "PMPI_Ireduce_scatter_block",
666 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, 0,
667 nullptr, -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
669 int* recvcounts = new int[count];
670 for (int i = 0; i < count; i++)
671 recvcounts[i] = recvcount;
672 if (request == MPI_REQUEST_IGNORED)
673 simgrid::smpi::colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm);
675 simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm, request);
678 TRACE_smpi_comm_out(rank);
683 int PMPI_Alltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
684 MPI_Datatype recvtype, MPI_Comm comm){
685 return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
688 int PMPI_Ialltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
689 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
692 CHECK_BUFFER(1, sendbuf, sendcount)
693 CHECK_BUFFER(4, recvbuf, recvcount)
694 if(sendbuf != MPI_IN_PLACE)
695 CHECK_TYPE(3, sendtype)
696 CHECK_TYPE(6, recvtype)
697 CHECK_COUNT(5, recvcount)
698 if(sendbuf != MPI_IN_PLACE)
699 CHECK_COUNT(2, sendcount)
700 CHECK_COUNT(5, recvcount)
704 int rank = simgrid::s4u::this_actor::get_pid();
705 int real_sendcount = sendcount;
706 MPI_Datatype real_sendtype = sendtype;
708 std::unique_ptr<unsigned char[]> tmp_sendbuf;
709 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, recvcount * comm->size(), recvtype);
711 if (sendbuf == MPI_IN_PLACE) {
712 real_sendcount = recvcount;
713 real_sendtype = recvtype;
716 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoall" : "PMPI_Ialltoall",
717 new simgrid::instr::CollTIData(
718 request == MPI_REQUEST_IGNORED ? "alltoall" : "ialltoall", -1, -1.0,
719 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
720 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
721 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
723 if (request == MPI_REQUEST_IGNORED)
725 simgrid::smpi::colls::alltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, comm);
727 retval = simgrid::smpi::colls::ialltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype,
730 TRACE_smpi_comm_out(rank);
735 int PMPI_Alltoallv(const void* sendbuf, const int* sendcounts, const int* senddispls, MPI_Datatype sendtype, void* recvbuf,
736 const int* recvcounts, const int* recvdispls, MPI_Datatype recvtype, MPI_Comm comm)
738 return PMPI_Ialltoallv(sendbuf, sendcounts, senddispls, sendtype, recvbuf, recvcounts, recvdispls, recvtype, comm, MPI_REQUEST_IGNORED);
741 int PMPI_Ialltoallv(const void* sendbuf, const int* sendcounts, const int* senddispls, MPI_Datatype sendtype, void* recvbuf,
742 const int* recvcounts, const int* recvdispls, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
745 if(sendbuf != MPI_IN_PLACE){
746 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
747 CHECK_NULL(3, MPI_ERR_ARG, senddispls)
748 CHECK_TYPE(4, sendtype)
750 CHECK_TYPE(8, recvtype)
751 CHECK_NULL(6, MPI_ERR_COUNT, recvcounts)
752 CHECK_NULL(7, MPI_ERR_ARG, recvdispls)
755 int rank = simgrid::s4u::this_actor::get_pid();
756 int size = comm->size();
757 for (int i = 0; i < size; i++) {
758 if(sendbuf != MPI_IN_PLACE){
759 CHECK_BUFFER(1, sendbuf, sendcounts[i])
760 CHECK_COUNT(2, sendcounts[i])
762 CHECK_BUFFER(5, recvbuf, recvcounts[i])
763 CHECK_COUNT(6, recvcounts[i])
769 std::vector<int>* trace_sendcounts = new std::vector<int>;
770 std::vector<int>* trace_recvcounts = new std::vector<int>;
771 int dt_size_recv = recvtype->size();
773 const int* real_sendcounts = sendcounts;
774 const int* real_senddispls = senddispls;
775 MPI_Datatype real_sendtype = sendtype;
777 for (int i = 0; i < size; i++) { // copy data to avoid bad free
778 recv_size += recvcounts[i] * dt_size_recv;
779 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
780 if (((recvdispls[i] + recvcounts[i]) * dt_size_recv) > maxsize)
781 maxsize = (recvdispls[i] + recvcounts[i]) * dt_size_recv;
784 std::unique_ptr<unsigned char[]> tmp_sendbuf;
785 std::unique_ptr<int[]> tmp_sendcounts;
786 std::unique_ptr<int[]> tmp_senddispls;
787 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
788 if (sendbuf == MPI_IN_PLACE) {
789 tmp_sendcounts.reset(new int[size]);
790 std::copy(recvcounts, recvcounts + size, tmp_sendcounts.get());
791 real_sendcounts = tmp_sendcounts.get();
792 tmp_senddispls.reset(new int[size]);
793 std::copy(recvdispls, recvdispls + size, tmp_senddispls.get());
794 real_senddispls = tmp_senddispls.get();
795 real_sendtype = recvtype;
798 int dt_size_send = real_sendtype->size();
800 for (int i = 0; i < size; i++) { // copy data to avoid bad free
801 send_size += real_sendcounts[i] * dt_size_send;
802 trace_sendcounts->push_back(real_sendcounts[i] * dt_size_send);
805 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallv" : "PMPI_Ialltoallv",
806 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
807 send_size, trace_sendcounts, recv_size, trace_recvcounts,
808 simgrid::smpi::Datatype::encode(real_sendtype),
809 simgrid::smpi::Datatype::encode(recvtype)));
812 if (request == MPI_REQUEST_IGNORED)
813 retval = simgrid::smpi::colls::alltoallv(real_sendbuf, real_sendcounts, real_senddispls, real_sendtype, recvbuf,
814 recvcounts, recvdispls, recvtype, comm);
816 retval = simgrid::smpi::colls::ialltoallv(real_sendbuf, real_sendcounts, real_senddispls, real_sendtype, recvbuf,
817 recvcounts, recvdispls, recvtype, comm, request);
819 TRACE_smpi_comm_out(rank);
824 int PMPI_Alltoallw(const void* sendbuf, const int* sendcounts, const int* senddispls, const MPI_Datatype* sendtypes, void* recvbuf,
825 const int* recvcounts, const int* recvdispls, const MPI_Datatype* recvtypes, MPI_Comm comm)
827 return PMPI_Ialltoallw(sendbuf, sendcounts, senddispls, sendtypes, recvbuf, recvcounts, recvdispls, recvtypes, comm, MPI_REQUEST_IGNORED);
830 int PMPI_Ialltoallw(const void* sendbuf, const int* sendcounts, const int* senddispls, const MPI_Datatype* sendtypes, void* recvbuf,
831 const int* recvcounts, const int* recvdispls, const MPI_Datatype* recvtypes, MPI_Comm comm, MPI_Request* request)
834 if(sendbuf != MPI_IN_PLACE){
835 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
836 CHECK_NULL(3, MPI_ERR_ARG, senddispls)
837 CHECK_NULL(4, MPI_ERR_TYPE, sendtypes)
839 CHECK_NULL(6, MPI_ERR_COUNT, recvcounts)
840 CHECK_NULL(7, MPI_ERR_ARG, recvdispls)
841 CHECK_NULL(8, MPI_ERR_TYPE, recvtypes)
843 int rank = simgrid::s4u::this_actor::get_pid();
844 int size = comm->size();
845 for (int i = 0; i < size; i++) {
846 if(sendbuf != MPI_IN_PLACE){
847 CHECK_BUFFER(1, sendbuf, sendcounts[i])
848 CHECK_COUNT(2, sendcounts[i])
849 CHECK_TYPE(4, sendtypes[i])
851 CHECK_BUFFER(5, recvbuf, recvcounts[i])
852 CHECK_COUNT(6, recvcounts[i])
853 CHECK_TYPE(8, recvtypes[i])
860 std::vector<int>* trace_sendcounts = new std::vector<int>;
861 std::vector<int>* trace_recvcounts = new std::vector<int>;
863 const int* real_sendcounts = sendcounts;
864 const int* real_senddispls = senddispls;
865 const MPI_Datatype* real_sendtypes = sendtypes;
866 unsigned long maxsize = 0;
867 for (int i = 0; i < size; i++) { // copy data to avoid bad free
868 if (recvtypes[i] == MPI_DATATYPE_NULL) {
869 delete trace_recvcounts;
870 delete trace_sendcounts;
873 recv_size += recvcounts[i] * recvtypes[i]->size();
874 trace_recvcounts->push_back(recvcounts[i] * recvtypes[i]->size());
875 if ((recvdispls[i] + (recvcounts[i] * recvtypes[i]->size())) > maxsize)
876 maxsize = recvdispls[i] + (recvcounts[i] * recvtypes[i]->size());
879 std::unique_ptr<unsigned char[]> tmp_sendbuf;
880 std::unique_ptr<int[]> tmp_sendcounts;
881 std::unique_ptr<int[]> tmp_senddispls;
882 std::unique_ptr<MPI_Datatype[]> tmp_sendtypes;
883 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
884 if (sendbuf == MPI_IN_PLACE) {
885 tmp_sendcounts.reset(new int[size]);
886 std::copy(recvcounts, recvcounts + size, tmp_sendcounts.get());
887 real_sendcounts = tmp_sendcounts.get();
888 tmp_senddispls.reset(new int[size]);
889 std::copy(recvdispls, recvdispls + size, tmp_senddispls.get());
890 real_senddispls = tmp_senddispls.get();
891 tmp_sendtypes.reset(new MPI_Datatype[size]);
892 std::copy(recvtypes, recvtypes + size, tmp_sendtypes.get());
893 real_sendtypes = tmp_sendtypes.get();
896 for (int i = 0; i < size; i++) { // copy data to avoid bad free
897 send_size += real_sendcounts[i] * real_sendtypes[i]->size();
898 trace_sendcounts->push_back(real_sendcounts[i] * real_sendtypes[i]->size());
901 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallw" : "PMPI_Ialltoallw",
902 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
903 send_size, trace_sendcounts, recv_size, trace_recvcounts,
904 simgrid::smpi::Datatype::encode(real_sendtypes[0]),
905 simgrid::smpi::Datatype::encode(recvtypes[0])));
908 if (request == MPI_REQUEST_IGNORED)
909 retval = simgrid::smpi::colls::alltoallw(real_sendbuf, real_sendcounts, real_senddispls, real_sendtypes, recvbuf,
910 recvcounts, recvdispls, recvtypes, comm);
912 retval = simgrid::smpi::colls::ialltoallw(real_sendbuf, real_sendcounts, real_senddispls, real_sendtypes, recvbuf,
913 recvcounts, recvdispls, recvtypes, comm, request);
915 TRACE_smpi_comm_out(rank);