1 /* Copyright (c) 2007-2020. The SimGrid Team. All rights reserved. */
3 /* This program is free software; you can redistribute it and/or modify it
4 * under the terms of the license (GNU LGPL) which comes with this package. */
7 #include "smpi_coll.hpp"
8 #include "smpi_comm.hpp"
9 #include "smpi_request.hpp"
10 #include "smpi_datatype_derived.hpp"
11 #include "smpi_op.hpp"
12 #include "src/smpi/include/smpi_actor.hpp"
14 XBT_LOG_EXTERNAL_DEFAULT_CATEGORY(smpi_pmpi);
16 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){
17 if (inplacebuf == MPI_IN_PLACE) {
18 tmp_sendbuf.reset(new unsigned char[count * datatype->get_extent()]);
19 simgrid::smpi::Datatype::copy(otherbuf, count, datatype, tmp_sendbuf.get(), count, datatype);
20 return tmp_sendbuf.get();
25 /* PMPI User level calls */
27 int PMPI_Barrier(MPI_Comm comm)
29 return PMPI_Ibarrier(comm, MPI_REQUEST_IGNORED);
32 int PMPI_Ibarrier(MPI_Comm comm, MPI_Request *request)
38 int rank = simgrid::s4u::this_actor::get_pid();
39 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Barrier" : "PMPI_Ibarrier",
40 new simgrid::instr::NoOpTIData(request == MPI_REQUEST_IGNORED ? "barrier" : "ibarrier"));
41 if (request == MPI_REQUEST_IGNORED) {
42 simgrid::smpi::colls::barrier(comm);
43 // Barrier can be used to synchronize RMA calls. Finish all requests from comm before.
44 comm->finish_rma_calls();
46 simgrid::smpi::colls::ibarrier(comm, request);
48 TRACE_smpi_comm_out(rank);
53 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
55 return PMPI_Ibcast(buf, count, datatype, root, comm, MPI_REQUEST_IGNORED);
58 int PMPI_Ibcast(void *buf, int count, MPI_Datatype datatype,
59 int root, MPI_Comm comm, MPI_Request* request)
62 CHECK_BUFFER(1, buf, count)
64 CHECK_TYPE(3, datatype)
69 int rank = simgrid::s4u::this_actor::get_pid();
70 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Bcast" : "PMPI_Ibcast",
71 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "bcast" : "ibcast", root, -1.0,
72 datatype->is_replayable() ? count : count * datatype->size(), -1,
73 simgrid::smpi::Datatype::encode(datatype), ""));
74 if (comm->size() > 1) {
75 if (request == MPI_REQUEST_IGNORED)
76 simgrid::smpi::colls::bcast(buf, count, datatype, root, comm);
78 simgrid::smpi::colls::ibcast(buf, count, datatype, root, comm, request);
80 if (request != MPI_REQUEST_IGNORED)
81 *request = MPI_REQUEST_NULL;
84 TRACE_smpi_comm_out(rank);
89 int PMPI_Gather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
90 int root, MPI_Comm comm){
91 return PMPI_Igather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
94 int PMPI_Igather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
95 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
98 if(sendbuf != MPI_IN_PLACE){
99 CHECK_BUFFER(1,sendbuf, sendcount)
100 CHECK_COUNT(2, sendcount)
101 CHECK_TYPE(3, sendtype)
103 if(comm->rank() == root){
104 CHECK_TYPE(6, recvtype)
105 CHECK_COUNT(5, recvcount)
106 CHECK_BUFFER(4, recvbuf, recvcount)
112 const void* real_sendbuf = sendbuf;
113 int real_sendcount = sendcount;
114 MPI_Datatype real_sendtype = sendtype;
115 if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
117 real_sendtype = recvtype;
119 int rank = simgrid::s4u::this_actor::get_pid();
121 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gather" : "PMPI_Igather",
122 new simgrid::instr::CollTIData(
123 request == MPI_REQUEST_IGNORED ? "gather" : "igather", root, -1.0,
124 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
125 (comm->rank() != root || recvtype->is_replayable()) ? recvcount : recvcount * recvtype->size(),
126 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
127 if (request == MPI_REQUEST_IGNORED)
128 simgrid::smpi::colls::gather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm);
130 simgrid::smpi::colls::igather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm,
133 TRACE_smpi_comm_out(rank);
138 int PMPI_Gatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, const int *recvcounts, const int *displs,
139 MPI_Datatype recvtype, int root, MPI_Comm comm){
140 return PMPI_Igatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, root, comm, MPI_REQUEST_IGNORED);
143 int PMPI_Igatherv(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, MPI_Request* request)
147 CHECK_BUFFER(1, sendbuf, sendcount)
148 if(sendbuf != MPI_IN_PLACE){
149 CHECK_TYPE(3, sendtype)
150 CHECK_COUNT(2, sendcount)
152 if(comm->rank() == root){
153 CHECK_TYPE(6, recvtype)
154 CHECK_NULL(5, MPI_ERR_COUNT, recvcounts)
155 CHECK_NULL(6, MPI_ERR_ARG, displs)
160 if (comm->rank() == root){
161 for (int i = 0; i < comm->size(); i++) {
162 CHECK_COUNT(5, recvcounts[i])
163 CHECK_BUFFER(4,recvbuf,recvcounts[i])
168 const void* real_sendbuf = sendbuf;
169 int real_sendcount = sendcount;
170 MPI_Datatype real_sendtype = sendtype;
171 if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
173 real_sendtype = recvtype;
176 int rank = simgrid::s4u::this_actor::get_pid();
177 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
179 std::vector<int>* trace_recvcounts = new std::vector<int>();
180 if (comm->rank() == root) {
181 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
182 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
185 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gatherv" : "PMPI_Igatherv",
186 new simgrid::instr::VarCollTIData(
187 request == MPI_REQUEST_IGNORED ? "gatherv" : "igatherv", root,
188 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
189 nullptr, dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(real_sendtype),
190 simgrid::smpi::Datatype::encode(recvtype)));
191 if (request == MPI_REQUEST_IGNORED)
192 simgrid::smpi::colls::gatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
195 simgrid::smpi::colls::igatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
196 root, comm, request);
198 TRACE_smpi_comm_out(rank);
203 int PMPI_Allgather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
204 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm){
205 return PMPI_Iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
208 int PMPI_Iallgather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
209 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
212 CHECK_BUFFER(1, sendbuf, sendcount)
213 CHECK_BUFFER(4, recvbuf, recvcount)
214 if(sendbuf != MPI_IN_PLACE){
215 CHECK_COUNT(2, sendcount)
216 CHECK_TYPE(3, sendtype)
218 CHECK_TYPE(6, recvtype)
219 CHECK_COUNT(5, recvcount)
223 if (sendbuf == MPI_IN_PLACE) {
224 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * recvcount * comm->rank();
225 sendcount = recvcount;
228 int rank = simgrid::s4u::this_actor::get_pid();
230 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgather" : "PMPI_Iallggather",
231 new simgrid::instr::CollTIData(
232 request == MPI_REQUEST_IGNORED ? "allgather" : "iallgather", -1, -1.0,
233 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
234 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
235 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
236 if (request == MPI_REQUEST_IGNORED)
237 simgrid::smpi::colls::allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
239 simgrid::smpi::colls::iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, request);
241 TRACE_smpi_comm_out(rank);
246 int PMPI_Allgatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
247 void *recvbuf, const int *recvcounts, const int *displs, MPI_Datatype recvtype, MPI_Comm comm){
248 return PMPI_Iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, MPI_REQUEST_IGNORED);
251 int PMPI_Iallgatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
252 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
255 CHECK_BUFFER(1, sendbuf, sendcount)
256 if(sendbuf != MPI_IN_PLACE)
257 CHECK_TYPE(3, sendtype)
258 CHECK_TYPE(6, recvtype)
259 CHECK_NULL(5, MPI_ERR_COUNT, recvcounts)
260 CHECK_NULL(6, MPI_ERR_ARG, displs)
261 if(sendbuf != MPI_IN_PLACE)
262 CHECK_COUNT(2, sendcount)
264 for (int i = 0; i < comm->size(); i++) {
265 CHECK_COUNT(5, recvcounts[i])
266 CHECK_BUFFER(4, recvbuf, recvcounts[i])
270 if (sendbuf == MPI_IN_PLACE) {
271 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * displs[comm->rank()];
272 sendcount = recvcounts[comm->rank()];
275 int rank = simgrid::s4u::this_actor::get_pid();
276 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
278 std::vector<int>* trace_recvcounts = new std::vector<int>();
279 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
280 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
284 rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgatherv" : "PMPI_Iallgatherv",
285 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "allgatherv" : "iallgatherv", -1,
286 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(), nullptr,
287 dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
288 simgrid::smpi::Datatype::encode(recvtype)));
289 if (request == MPI_REQUEST_IGNORED)
290 simgrid::smpi::colls::allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
292 simgrid::smpi::colls::iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm,
295 TRACE_smpi_comm_out(rank);
300 int PMPI_Scatter(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
301 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
302 return PMPI_Iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
305 int PMPI_Iscatter(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
306 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
309 if(comm->rank() == root){
310 CHECK_BUFFER(1, sendbuf, sendcount)
311 CHECK_COUNT(2, sendcount)
312 CHECK_TYPE(3, sendtype)
314 if(recvbuf != MPI_IN_PLACE){
315 CHECK_BUFFER(4, recvbuf, recvcount)
316 CHECK_COUNT(5, recvcount)
317 CHECK_TYPE(6, recvtype)
323 if (recvbuf == MPI_IN_PLACE) {
325 recvcount = sendcount;
327 int rank = simgrid::s4u::this_actor::get_pid();
329 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatter" : "PMPI_Iscatter",
330 new simgrid::instr::CollTIData(
331 request == MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
332 (comm->rank() != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
333 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
334 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
335 if (request == MPI_REQUEST_IGNORED)
336 simgrid::smpi::colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
338 simgrid::smpi::colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
340 TRACE_smpi_comm_out(rank);
345 int PMPI_Scatterv(const void *sendbuf, const int *sendcounts, const int *displs,
346 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
347 return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
350 int PMPI_Iscatterv(const void* sendbuf, const int* sendcounts, const int* displs, MPI_Datatype sendtype, void* recvbuf, int recvcount,
351 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
354 if(recvbuf != MPI_IN_PLACE){
355 CHECK_BUFFER(4, recvbuf, recvcount)
356 CHECK_COUNT(5, recvcount)
357 CHECK_TYPE(7, recvtype)
361 if (comm->rank() == root) {
362 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
363 CHECK_NULL(3, MPI_ERR_ARG, displs)
364 CHECK_TYPE(4, sendtype)
365 for (int i = 0; i < comm->size(); i++){
366 CHECK_BUFFER(1, sendbuf, sendcounts[i])
367 CHECK_COUNT(2, sendcounts[i])
369 if (recvbuf == MPI_IN_PLACE) {
371 recvcount = sendcounts[comm->rank()];
377 int rank = simgrid::s4u::this_actor::get_pid();
378 int dt_size_send = sendtype->is_replayable() ? 1 : sendtype->size();
380 std::vector<int>* trace_sendcounts = new std::vector<int>();
381 if (comm->rank() == root) {
382 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
383 trace_sendcounts->push_back(sendcounts[i] * dt_size_send);
387 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatterv" : "PMPI_Iscatterv",
388 new simgrid::instr::VarCollTIData(
389 request == MPI_REQUEST_IGNORED ? "scatterv" : "iscatterv", root, dt_size_send,
390 trace_sendcounts, recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
391 nullptr, simgrid::smpi::Datatype::encode(sendtype),
392 simgrid::smpi::Datatype::encode(recvtype)));
393 if (request == MPI_REQUEST_IGNORED)
394 simgrid::smpi::colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
396 simgrid::smpi::colls::iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm,
399 TRACE_smpi_comm_out(rank);
404 int PMPI_Reduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
406 return PMPI_Ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, MPI_REQUEST_IGNORED);
409 int PMPI_Ireduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm, MPI_Request* request)
412 CHECK_BUFFER(1, sendbuf, count)
413 if(comm->rank() == root)
414 CHECK_BUFFER(5, recvbuf, count)
415 CHECK_TYPE(4, datatype)
416 CHECK_COUNT(3, count)
422 int rank = simgrid::s4u::this_actor::get_pid();
424 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce" : "PMPI_Ireduce",
425 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "reduce" : "ireduce", root, 0,
426 datatype->is_replayable() ? count : count * datatype->size(), -1,
427 simgrid::smpi::Datatype::encode(datatype), ""));
428 if (request == MPI_REQUEST_IGNORED)
429 simgrid::smpi::colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
431 simgrid::smpi::colls::ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, request);
433 TRACE_smpi_comm_out(rank);
438 int PMPI_Reduce_local(const void* inbuf, void* inoutbuf, int count, MPI_Datatype datatype, MPI_Op op)
440 CHECK_BUFFER(1, inbuf, count)
441 CHECK_BUFFER(2, inoutbuf, count)
442 CHECK_TYPE(4, datatype)
443 CHECK_COUNT(3, count)
447 op->apply(inbuf, inoutbuf, &count, datatype);
452 int PMPI_Allreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
454 return PMPI_Iallreduce(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
457 int PMPI_Iallreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
460 CHECK_BUFFER(1, sendbuf, count)
461 CHECK_BUFFER(2, recvbuf, count)
462 CHECK_TYPE(4, datatype)
463 CHECK_COUNT(3, count)
468 std::unique_ptr<unsigned char[]> tmp_sendbuf;
469 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
471 int rank = simgrid::s4u::this_actor::get_pid();
473 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allreduce" : "PMPI_Iallreduce",
474 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "allreduce" : "iallreduce", -1, 0,
475 datatype->is_replayable() ? count : count * datatype->size(), -1,
476 simgrid::smpi::Datatype::encode(datatype), ""));
478 if (request == MPI_REQUEST_IGNORED)
479 simgrid::smpi::colls::allreduce(real_sendbuf, recvbuf, count, datatype, op, comm);
481 simgrid::smpi::colls::iallreduce(real_sendbuf, recvbuf, count, datatype, op, comm, request);
483 TRACE_smpi_comm_out(rank);
488 int PMPI_Scan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
490 return PMPI_Iscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
493 int PMPI_Iscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request)
496 CHECK_BUFFER(1,sendbuf,count)
497 CHECK_BUFFER(2,recvbuf,count)
498 CHECK_TYPE(4, datatype)
499 CHECK_COUNT(3, count)
504 int rank = simgrid::s4u::this_actor::get_pid();
505 std::unique_ptr<unsigned char[]> tmp_sendbuf;
506 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
508 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scan" : "PMPI_Iscan",
509 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "scan" : "iscan", -1,
510 datatype->is_replayable() ? count : count * datatype->size(),
511 simgrid::smpi::Datatype::encode(datatype)));
514 if (request == MPI_REQUEST_IGNORED)
515 retval = simgrid::smpi::colls::scan(real_sendbuf, recvbuf, count, datatype, op, comm);
517 retval = simgrid::smpi::colls::iscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
519 TRACE_smpi_comm_out(rank);
524 int PMPI_Exscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
526 return PMPI_Iexscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
529 int PMPI_Iexscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request){
531 CHECK_BUFFER(1, sendbuf, count)
532 CHECK_BUFFER(2, recvbuf, count)
533 CHECK_TYPE(4, datatype)
534 CHECK_COUNT(3, count)
539 int rank = simgrid::s4u::this_actor::get_pid();
540 std::unique_ptr<unsigned char[]> tmp_sendbuf;
541 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
543 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Exscan" : "PMPI_Iexscan",
544 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "exscan" : "iexscan", -1,
545 datatype->is_replayable() ? count : count * datatype->size(),
546 simgrid::smpi::Datatype::encode(datatype)));
549 if (request == MPI_REQUEST_IGNORED)
550 retval = simgrid::smpi::colls::exscan(real_sendbuf, recvbuf, count, datatype, op, comm);
552 retval = simgrid::smpi::colls::iexscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
554 TRACE_smpi_comm_out(rank);
559 int PMPI_Reduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
561 return PMPI_Ireduce_scatter(sendbuf, recvbuf, recvcounts, datatype, op, comm, MPI_REQUEST_IGNORED);
564 int PMPI_Ireduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
567 CHECK_TYPE(4, datatype)
568 CHECK_NULL(3, MPI_ERR_COUNT, recvcounts)
571 for (int i = 0; i < comm->size(); i++) {
572 CHECK_COUNT(3, recvcounts[i])
573 CHECK_BUFFER(1, sendbuf, recvcounts[i])
574 CHECK_BUFFER(2, recvbuf, recvcounts[i])
578 int rank = simgrid::s4u::this_actor::get_pid();
579 std::vector<int>* trace_recvcounts = new std::vector<int>();
580 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
583 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
584 trace_recvcounts->push_back(recvcounts[i] * dt_send_size);
585 totalcount += recvcounts[i];
587 std::unique_ptr<unsigned char[]> tmp_sendbuf;
588 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, totalcount, datatype);
590 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter" : "PMPI_Ireduce_scatter",
591 new simgrid::instr::VarCollTIData(
592 request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, dt_send_size, nullptr,
593 -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
595 if (request == MPI_REQUEST_IGNORED)
596 simgrid::smpi::colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm);
598 simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm, request);
600 TRACE_smpi_comm_out(rank);
605 int PMPI_Reduce_scatter_block(const void *sendbuf, void *recvbuf, int recvcount,
606 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
608 return PMPI_Ireduce_scatter_block(sendbuf, recvbuf, recvcount, datatype, op, comm, MPI_REQUEST_IGNORED);
611 int PMPI_Ireduce_scatter_block(const void* sendbuf, void* recvbuf, int recvcount, MPI_Datatype datatype, MPI_Op op,
612 MPI_Comm comm, MPI_Request* request)
615 CHECK_BUFFER(1, sendbuf, recvcount)
616 CHECK_BUFFER(2, recvbuf, recvcount)
617 CHECK_TYPE(4, datatype)
618 CHECK_COUNT(3, recvcount)
623 int count = comm->size();
625 int rank = simgrid::s4u::this_actor::get_pid();
626 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
627 std::vector<int>* trace_recvcounts = new std::vector<int>(recvcount * dt_send_size); // copy data to avoid bad free
628 std::unique_ptr<unsigned char[]> tmp_sendbuf;
629 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, recvcount * count, datatype);
632 rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter_block" : "PMPI_Ireduce_scatter_block",
633 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, 0,
634 nullptr, -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
636 int* recvcounts = new int[count];
637 for (int i = 0; i < count; i++)
638 recvcounts[i] = recvcount;
639 if (request == MPI_REQUEST_IGNORED)
640 simgrid::smpi::colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm);
642 simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm, request);
645 TRACE_smpi_comm_out(rank);
650 int PMPI_Alltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
651 MPI_Datatype recvtype, MPI_Comm comm){
652 return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
655 int PMPI_Ialltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
656 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
659 CHECK_BUFFER(1, sendbuf, sendcount)
660 CHECK_BUFFER(4, recvbuf, recvcount)
661 if(sendbuf != MPI_IN_PLACE)
662 CHECK_TYPE(3, sendtype)
663 CHECK_TYPE(6, recvtype)
664 CHECK_COUNT(5, recvcount)
665 if(sendbuf != MPI_IN_PLACE)
666 CHECK_COUNT(2, sendcount)
667 CHECK_COUNT(5, recvcount)
671 int rank = simgrid::s4u::this_actor::get_pid();
672 int real_sendcount = sendcount;
673 MPI_Datatype real_sendtype = sendtype;
675 std::unique_ptr<unsigned char[]> tmp_sendbuf;
676 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, recvcount * comm->size(), recvtype);
678 if (sendbuf == MPI_IN_PLACE) {
679 real_sendcount = recvcount;
680 real_sendtype = recvtype;
683 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoall" : "PMPI_Ialltoall",
684 new simgrid::instr::CollTIData(
685 request == MPI_REQUEST_IGNORED ? "alltoall" : "ialltoall", -1, -1.0,
686 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
687 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
688 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
690 if (request == MPI_REQUEST_IGNORED)
692 simgrid::smpi::colls::alltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, comm);
694 retval = simgrid::smpi::colls::ialltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype,
697 TRACE_smpi_comm_out(rank);
702 int PMPI_Alltoallv(const void* sendbuf, const int* sendcounts, const int* senddispls, MPI_Datatype sendtype, void* recvbuf,
703 const int* recvcounts, const int* recvdispls, MPI_Datatype recvtype, MPI_Comm comm)
705 return PMPI_Ialltoallv(sendbuf, sendcounts, senddispls, sendtype, recvbuf, recvcounts, recvdispls, recvtype, comm, MPI_REQUEST_IGNORED);
708 int PMPI_Ialltoallv(const void* sendbuf, const int* sendcounts, const int* senddispls, MPI_Datatype sendtype, void* recvbuf,
709 const int* recvcounts, const int* recvdispls, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
712 if(sendbuf != MPI_IN_PLACE){
713 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
714 CHECK_NULL(3, MPI_ERR_ARG, senddispls)
715 CHECK_TYPE(4, sendtype)
717 CHECK_TYPE(8, recvtype)
718 CHECK_NULL(6, MPI_ERR_COUNT, recvcounts)
719 CHECK_NULL(7, MPI_ERR_ARG, recvdispls)
722 int rank = simgrid::s4u::this_actor::get_pid();
723 int size = comm->size();
724 for (int i = 0; i < size; i++) {
725 if(sendbuf != MPI_IN_PLACE){
726 CHECK_BUFFER(1, sendbuf, sendcounts[i])
727 CHECK_COUNT(2, sendcounts[i])
729 CHECK_BUFFER(5, recvbuf, recvcounts[i])
730 CHECK_COUNT(6, recvcounts[i])
736 std::vector<int>* trace_sendcounts = new std::vector<int>();
737 std::vector<int>* trace_recvcounts = new std::vector<int>();
738 int dt_size_recv = recvtype->size();
740 const int* real_sendcounts = sendcounts;
741 const int* real_senddispls = senddispls;
742 MPI_Datatype real_sendtype = sendtype;
744 for (int i = 0; i < size; i++) { // copy data to avoid bad free
745 recv_size += recvcounts[i] * dt_size_recv;
746 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
747 if (((recvdispls[i] + recvcounts[i]) * dt_size_recv) > maxsize)
748 maxsize = (recvdispls[i] + recvcounts[i]) * dt_size_recv;
751 std::unique_ptr<unsigned char[]> tmp_sendbuf;
752 std::unique_ptr<int[]> tmp_sendcounts;
753 std::unique_ptr<int[]> tmp_senddispls;
754 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
755 if (sendbuf == MPI_IN_PLACE) {
756 tmp_sendcounts.reset(new int[size]);
757 std::copy(recvcounts, recvcounts + size, tmp_sendcounts.get());
758 real_sendcounts = tmp_sendcounts.get();
759 tmp_senddispls.reset(new int[size]);
760 std::copy(recvdispls, recvdispls + size, tmp_senddispls.get());
761 real_senddispls = tmp_senddispls.get();
762 real_sendtype = recvtype;
765 int dt_size_send = real_sendtype->size();
767 for (int i = 0; i < size; i++) { // copy data to avoid bad free
768 send_size += real_sendcounts[i] * dt_size_send;
769 trace_sendcounts->push_back(real_sendcounts[i] * dt_size_send);
772 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallv" : "PMPI_Ialltoallv",
773 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
774 send_size, trace_sendcounts, recv_size, trace_recvcounts,
775 simgrid::smpi::Datatype::encode(real_sendtype),
776 simgrid::smpi::Datatype::encode(recvtype)));
779 if (request == MPI_REQUEST_IGNORED)
780 retval = simgrid::smpi::colls::alltoallv(real_sendbuf, real_sendcounts, real_senddispls, real_sendtype, recvbuf,
781 recvcounts, recvdispls, recvtype, comm);
783 retval = simgrid::smpi::colls::ialltoallv(real_sendbuf, real_sendcounts, real_senddispls, real_sendtype, recvbuf,
784 recvcounts, recvdispls, recvtype, comm, request);
786 TRACE_smpi_comm_out(rank);
791 int PMPI_Alltoallw(const void* sendbuf, const int* sendcounts, const int* senddispls, const MPI_Datatype* sendtypes, void* recvbuf,
792 const int* recvcounts, const int* recvdispls, const MPI_Datatype* recvtypes, MPI_Comm comm)
794 return PMPI_Ialltoallw(sendbuf, sendcounts, senddispls, sendtypes, recvbuf, recvcounts, recvdispls, recvtypes, comm, MPI_REQUEST_IGNORED);
797 int PMPI_Ialltoallw(const void* sendbuf, const int* sendcounts, const int* senddispls, const MPI_Datatype* sendtypes, void* recvbuf,
798 const int* recvcounts, const int* recvdispls, const MPI_Datatype* recvtypes, MPI_Comm comm, MPI_Request* request)
801 if(sendbuf != MPI_IN_PLACE){
802 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
803 CHECK_NULL(3, MPI_ERR_ARG, senddispls)
804 CHECK_NULL(4, MPI_ERR_TYPE, sendtypes)
806 CHECK_NULL(6, MPI_ERR_COUNT, recvcounts)
807 CHECK_NULL(7, MPI_ERR_ARG, recvdispls)
808 CHECK_NULL(8, MPI_ERR_TYPE, recvtypes)
810 int rank = simgrid::s4u::this_actor::get_pid();
811 int size = comm->size();
812 for (int i = 0; i < size; i++) {
813 if(sendbuf != MPI_IN_PLACE){
814 CHECK_BUFFER(1, sendbuf, sendcounts[i])
815 CHECK_COUNT(2, sendcounts[i])
816 CHECK_TYPE(4, sendtypes[i])
818 CHECK_BUFFER(5, recvbuf, recvcounts[i])
819 CHECK_COUNT(6, recvcounts[i])
820 CHECK_TYPE(8, recvtypes[i])
827 std::vector<int>* trace_sendcounts = new std::vector<int>();
828 std::vector<int>* trace_recvcounts = new std::vector<int>();
830 const int* real_sendcounts = sendcounts;
831 const int* real_senddispls = senddispls;
832 const MPI_Datatype* real_sendtypes = sendtypes;
833 unsigned long maxsize = 0;
834 for (int i = 0; i < size; i++) { // copy data to avoid bad free
835 if (recvtypes[i] == MPI_DATATYPE_NULL) {
836 delete trace_recvcounts;
837 delete trace_sendcounts;
840 recv_size += recvcounts[i] * recvtypes[i]->size();
841 trace_recvcounts->push_back(recvcounts[i] * recvtypes[i]->size());
842 if ((recvdispls[i] + (recvcounts[i] * recvtypes[i]->size())) > maxsize)
843 maxsize = recvdispls[i] + (recvcounts[i] * recvtypes[i]->size());
846 std::unique_ptr<unsigned char[]> tmp_sendbuf;
847 std::unique_ptr<int[]> tmp_sendcounts;
848 std::unique_ptr<int[]> tmp_senddispls;
849 std::unique_ptr<MPI_Datatype[]> tmp_sendtypes;
850 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
851 if (sendbuf == MPI_IN_PLACE) {
852 tmp_sendcounts.reset(new int[size]);
853 std::copy(recvcounts, recvcounts + size, tmp_sendcounts.get());
854 real_sendcounts = tmp_sendcounts.get();
855 tmp_senddispls.reset(new int[size]);
856 std::copy(recvdispls, recvdispls + size, tmp_senddispls.get());
857 real_senddispls = tmp_senddispls.get();
858 tmp_sendtypes.reset(new MPI_Datatype[size]);
859 std::copy(recvtypes, recvtypes + size, tmp_sendtypes.get());
860 real_sendtypes = tmp_sendtypes.get();
863 for (int i = 0; i < size; i++) { // copy data to avoid bad free
864 send_size += real_sendcounts[i] * real_sendtypes[i]->size();
865 trace_sendcounts->push_back(real_sendcounts[i] * real_sendtypes[i]->size());
868 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallw" : "PMPI_Ialltoallw",
869 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
870 send_size, trace_sendcounts, recv_size, trace_recvcounts,
871 simgrid::smpi::Datatype::encode(real_sendtypes[0]),
872 simgrid::smpi::Datatype::encode(recvtypes[0])));
875 if (request == MPI_REQUEST_IGNORED)
876 retval = simgrid::smpi::colls::alltoallw(real_sendbuf, real_sendcounts, real_senddispls, real_sendtypes, recvbuf,
877 recvcounts, recvdispls, recvtypes, comm);
879 retval = simgrid::smpi::colls::ialltoallw(real_sendbuf, real_sendcounts, real_senddispls, real_sendtypes, recvbuf,
880 recvcounts, recvdispls, recvtypes, comm, request);
882 TRACE_smpi_comm_out(rank);