1 /* Copyright (c) 2007-2021. 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"
16 XBT_LOG_EXTERNAL_DEFAULT_CATEGORY(smpi_pmpi);
18 static const void* smpi_get_in_place_buf(const void* inplacebuf, const void* otherbuf,
19 std::vector<unsigned char>& tmp_sendbuf, int count, MPI_Datatype datatype)
21 if (inplacebuf == MPI_IN_PLACE) {
22 tmp_sendbuf.resize(count * datatype->get_extent());
23 simgrid::smpi::Datatype::copy(otherbuf, count, datatype, tmp_sendbuf.data(), count, datatype);
24 return tmp_sendbuf.data();
29 /* PMPI User level calls */
31 int PMPI_Barrier(MPI_Comm comm)
33 return PMPI_Ibarrier(comm, MPI_REQUEST_IGNORED);
36 int PMPI_Ibarrier(MPI_Comm comm, MPI_Request *request)
42 int rank = simgrid::s4u::this_actor::get_pid();
43 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Barrier" : "PMPI_Ibarrier",
44 new simgrid::instr::NoOpTIData(request == MPI_REQUEST_IGNORED ? "barrier" : "ibarrier"));
45 if (request == MPI_REQUEST_IGNORED) {
46 simgrid::smpi::colls::barrier(comm);
47 // Barrier can be used to synchronize RMA calls. Finish all requests from comm before.
48 comm->finish_rma_calls();
50 simgrid::smpi::colls::ibarrier(comm, request);
52 TRACE_smpi_comm_out(rank);
57 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
59 return PMPI_Ibcast(buf, count, datatype, root, comm, MPI_REQUEST_IGNORED);
62 int PMPI_Ibcast(void *buf, int count, MPI_Datatype datatype,
63 int root, MPI_Comm comm, MPI_Request* request)
66 CHECK_BUFFER(1, buf, count)
68 CHECK_TYPE(3, datatype)
73 int rank = simgrid::s4u::this_actor::get_pid();
74 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Bcast" : "PMPI_Ibcast",
75 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "bcast" : "ibcast", root, -1.0,
76 datatype->is_replayable() ? count : count * datatype->size(), -1,
77 simgrid::smpi::Datatype::encode(datatype), ""));
78 if (comm->size() > 1) {
79 if (request == MPI_REQUEST_IGNORED)
80 simgrid::smpi::colls::bcast(buf, count, datatype, root, comm);
82 simgrid::smpi::colls::ibcast(buf, count, datatype, root, comm, request);
84 if (request != MPI_REQUEST_IGNORED)
85 *request = MPI_REQUEST_NULL;
88 TRACE_smpi_comm_out(rank);
93 int PMPI_Gather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
94 int root, MPI_Comm comm){
95 return PMPI_Igather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
98 int PMPI_Igather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
99 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
102 if(sendbuf != MPI_IN_PLACE){
103 CHECK_BUFFER(1,sendbuf, sendcount)
104 CHECK_COUNT(2, sendcount)
105 CHECK_TYPE(3, sendtype)
107 if(comm->rank() == root){
108 CHECK_TYPE(6, recvtype)
109 CHECK_COUNT(5, recvcount)
110 CHECK_BUFFER(4, recvbuf, recvcount)
115 const void* real_sendbuf = sendbuf;
116 int real_sendcount = sendcount;
117 MPI_Datatype real_sendtype = sendtype;
118 if (comm->rank() == root){
119 if (sendbuf == MPI_IN_PLACE) {
121 real_sendtype = recvtype;
122 } else if(recvtype->size() * recvcount != sendtype->size() * sendcount){
123 XBT_WARN("MPI_(I)Gather : received size at root differs from sent size %lu %lu", recvtype->size() * recvcount , sendtype->size() * sendcount);
124 return MPI_ERR_TRUNCATE;
130 int rank = simgrid::s4u::this_actor::get_pid();
132 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gather" : "PMPI_Igather",
133 new simgrid::instr::CollTIData(
134 request == MPI_REQUEST_IGNORED ? "gather" : "igather", root, -1.0,
135 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
136 (comm->rank() != root || recvtype->is_replayable()) ? recvcount : recvcount * recvtype->size(),
137 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
138 if (request == MPI_REQUEST_IGNORED)
139 simgrid::smpi::colls::gather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm);
141 simgrid::smpi::colls::igather(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, root, comm,
144 TRACE_smpi_comm_out(rank);
149 int PMPI_Gatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, const int *recvcounts, const int *displs,
150 MPI_Datatype recvtype, int root, MPI_Comm comm){
151 return PMPI_Igatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, root, comm, MPI_REQUEST_IGNORED);
154 int PMPI_Igatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
155 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
158 CHECK_BUFFER(1, sendbuf, sendcount)
159 if(sendbuf != MPI_IN_PLACE){
160 CHECK_TYPE(3, sendtype)
161 CHECK_COUNT(2, sendcount)
163 if(comm->rank() == root){
164 CHECK_TYPE(6, recvtype)
165 CHECK_NULL(5, MPI_ERR_COUNT, recvcounts)
166 CHECK_NULL(6, MPI_ERR_ARG, displs)
171 if (comm->rank() == root){
172 for (int i = 0; i < comm->size(); i++) {
173 CHECK_COUNT(5, recvcounts[i])
174 CHECK_BUFFER(4,recvbuf,recvcounts[i])
179 const void* real_sendbuf = sendbuf;
180 int real_sendcount = sendcount;
181 MPI_Datatype real_sendtype = sendtype;
182 if ((comm->rank() == root) && (sendbuf == MPI_IN_PLACE)) {
184 real_sendtype = recvtype;
187 int rank = simgrid::s4u::this_actor::get_pid();
188 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
190 auto* trace_recvcounts = new std::vector<int>();
191 if (comm->rank() == root) {
192 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
193 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
196 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Gatherv" : "PMPI_Igatherv",
197 new simgrid::instr::VarCollTIData(
198 request == MPI_REQUEST_IGNORED ? "gatherv" : "igatherv", root,
199 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
200 nullptr, dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(real_sendtype),
201 simgrid::smpi::Datatype::encode(recvtype)));
202 if (request == MPI_REQUEST_IGNORED)
203 simgrid::smpi::colls::gatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
206 simgrid::smpi::colls::igatherv(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcounts, displs, recvtype,
207 root, comm, request);
209 TRACE_smpi_comm_out(rank);
214 int PMPI_Allgather(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
215 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm){
216 return PMPI_Iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
219 int PMPI_Iallgather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
220 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
223 CHECK_BUFFER(1, sendbuf, sendcount)
224 CHECK_BUFFER(4, recvbuf, recvcount)
225 if(sendbuf != MPI_IN_PLACE){
226 CHECK_COUNT(2, sendcount)
227 CHECK_TYPE(3, sendtype)
229 CHECK_TYPE(6, recvtype)
230 CHECK_COUNT(5, recvcount)
233 if (sendbuf == MPI_IN_PLACE) {
234 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * recvcount * comm->rank();
235 sendcount = recvcount;
239 if(recvtype->size() * recvcount != sendtype->size() * sendcount){
240 XBT_WARN("MPI_(I)Allgather : received size from each process differs from sent size");
241 return MPI_ERR_TRUNCATE;
246 int rank = simgrid::s4u::this_actor::get_pid();
248 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgather" : "PMPI_Iallggather",
249 new simgrid::instr::CollTIData(
250 request == MPI_REQUEST_IGNORED ? "allgather" : "iallgather", -1, -1.0,
251 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
252 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
253 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
254 if (request == MPI_REQUEST_IGNORED)
255 simgrid::smpi::colls::allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
257 simgrid::smpi::colls::iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, request);
259 TRACE_smpi_comm_out(rank);
264 int PMPI_Allgatherv(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
265 void *recvbuf, const int *recvcounts, const int *displs, MPI_Datatype recvtype, MPI_Comm comm){
266 return PMPI_Iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, MPI_REQUEST_IGNORED);
269 int PMPI_Iallgatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int* recvcounts, const int* displs,
270 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
273 CHECK_BUFFER(1, sendbuf, sendcount)
274 if(sendbuf != MPI_IN_PLACE)
275 CHECK_TYPE(3, sendtype)
276 CHECK_TYPE(6, recvtype)
277 CHECK_NULL(5, MPI_ERR_COUNT, recvcounts)
278 CHECK_NULL(6, MPI_ERR_ARG, displs)
279 if(sendbuf != MPI_IN_PLACE)
280 CHECK_COUNT(2, sendcount)
282 for (int i = 0; i < comm->size(); i++) {
283 CHECK_COUNT(5, recvcounts[i])
284 CHECK_BUFFER(4, recvbuf, recvcounts[i])
288 if (sendbuf == MPI_IN_PLACE) {
289 sendbuf = static_cast<char*>(recvbuf) + recvtype->get_extent() * displs[comm->rank()];
290 sendcount = recvcounts[comm->rank()];
293 int rank = simgrid::s4u::this_actor::get_pid();
294 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
296 auto* trace_recvcounts = new std::vector<int>();
297 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
298 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
302 rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allgatherv" : "PMPI_Iallgatherv",
303 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "allgatherv" : "iallgatherv", -1,
304 sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(), nullptr,
305 dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
306 simgrid::smpi::Datatype::encode(recvtype)));
307 if (request == MPI_REQUEST_IGNORED)
308 simgrid::smpi::colls::allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
310 simgrid::smpi::colls::iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm,
313 TRACE_smpi_comm_out(rank);
318 int PMPI_Scatter(const void *sendbuf, int sendcount, MPI_Datatype sendtype,
319 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
320 return PMPI_Iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
323 int PMPI_Iscatter(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
324 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
327 if(comm->rank() == root){
328 CHECK_BUFFER(1, sendbuf, sendcount)
329 CHECK_COUNT(2, sendcount)
330 CHECK_TYPE(3, sendtype)
332 if(recvbuf != MPI_IN_PLACE){
333 CHECK_BUFFER(4, recvbuf, recvcount)
334 CHECK_COUNT(5, recvcount)
335 CHECK_TYPE(6, recvtype)
340 if (recvbuf == MPI_IN_PLACE) {
342 recvcount = sendcount;
345 if((comm->rank() == root) && (recvtype->size() * recvcount != sendtype->size() * sendcount)){
346 XBT_WARN("MPI_(I)Scatter : sent size to each process differs from receive size");
347 return MPI_ERR_TRUNCATE;
352 int rank = simgrid::s4u::this_actor::get_pid();
354 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scatter" : "PMPI_Iscatter",
355 new simgrid::instr::CollTIData(
356 request == MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
357 (comm->rank() != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
358 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
359 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
360 if (request == MPI_REQUEST_IGNORED)
361 simgrid::smpi::colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
363 simgrid::smpi::colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
365 TRACE_smpi_comm_out(rank);
370 int PMPI_Scatterv(const void *sendbuf, const int *sendcounts, const int *displs,
371 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
372 return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
375 int PMPI_Iscatterv(const void* sendbuf, const int* sendcounts, const int* displs, MPI_Datatype sendtype, void* recvbuf, int recvcount,
376 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
379 if(recvbuf != MPI_IN_PLACE){
380 CHECK_BUFFER(4, recvbuf, recvcount)
381 CHECK_COUNT(5, recvcount)
382 CHECK_TYPE(7, recvtype)
386 if (comm->rank() == root) {
387 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
388 CHECK_NULL(3, MPI_ERR_ARG, displs)
389 CHECK_TYPE(4, sendtype)
390 for (int i = 0; i < comm->size(); i++){
391 CHECK_BUFFER(1, sendbuf, sendcounts[i])
392 CHECK_COUNT(2, sendcounts[i])
394 if (recvbuf == MPI_IN_PLACE) {
396 recvcount = sendcounts[comm->rank()];
402 int rank = simgrid::s4u::this_actor::get_pid();
403 int dt_size_send = sendtype->is_replayable() ? 1 : sendtype->size();
405 auto* 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)
437 CHECK_BUFFER(1, sendbuf, count)
438 if(comm->rank() == root)
439 CHECK_BUFFER(5, recvbuf, count)
440 CHECK_TYPE(4, datatype)
441 CHECK_COUNT(3, count)
447 int rank = simgrid::s4u::this_actor::get_pid();
449 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce" : "PMPI_Ireduce",
450 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "reduce" : "ireduce", root, 0,
451 datatype->is_replayable() ? count : count * datatype->size(), -1,
452 simgrid::smpi::Datatype::encode(datatype), ""));
453 if (request == MPI_REQUEST_IGNORED)
454 simgrid::smpi::colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
456 simgrid::smpi::colls::ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, request);
458 TRACE_smpi_comm_out(rank);
463 int PMPI_Reduce_local(const void* inbuf, void* inoutbuf, int count, MPI_Datatype datatype, MPI_Op op)
465 CHECK_BUFFER(1, inbuf, count)
466 CHECK_BUFFER(2, inoutbuf, count)
467 CHECK_TYPE(4, datatype)
468 CHECK_COUNT(3, count)
472 op->apply(inbuf, inoutbuf, &count, datatype);
477 int PMPI_Allreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
479 return PMPI_Iallreduce(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
482 int PMPI_Iallreduce(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
485 CHECK_BUFFER(1, sendbuf, count)
486 CHECK_BUFFER(2, recvbuf, count)
487 CHECK_TYPE(4, datatype)
488 CHECK_COUNT(3, count)
493 std::vector<unsigned char> tmp_sendbuf;
494 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
496 int rank = simgrid::s4u::this_actor::get_pid();
498 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Allreduce" : "PMPI_Iallreduce",
499 new simgrid::instr::CollTIData(request == MPI_REQUEST_IGNORED ? "allreduce" : "iallreduce", -1, 0,
500 datatype->is_replayable() ? count : count * datatype->size(), -1,
501 simgrid::smpi::Datatype::encode(datatype), ""));
503 if (request == MPI_REQUEST_IGNORED)
504 simgrid::smpi::colls::allreduce(real_sendbuf, recvbuf, count, datatype, op, comm);
506 simgrid::smpi::colls::iallreduce(real_sendbuf, recvbuf, count, datatype, op, comm, request);
508 TRACE_smpi_comm_out(rank);
513 int PMPI_Scan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
515 return PMPI_Iscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
518 int PMPI_Iscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request)
521 CHECK_BUFFER(1,sendbuf,count)
522 CHECK_BUFFER(2,recvbuf,count)
523 CHECK_TYPE(4, datatype)
524 CHECK_COUNT(3, count)
529 int rank = simgrid::s4u::this_actor::get_pid();
530 std::vector<unsigned char> tmp_sendbuf;
531 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
533 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Scan" : "PMPI_Iscan",
534 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "scan" : "iscan", -1,
535 datatype->is_replayable() ? count : count * datatype->size(),
536 simgrid::smpi::Datatype::encode(datatype)));
539 if (request == MPI_REQUEST_IGNORED)
540 retval = simgrid::smpi::colls::scan(real_sendbuf, recvbuf, count, datatype, op, comm);
542 retval = simgrid::smpi::colls::iscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
544 TRACE_smpi_comm_out(rank);
549 int PMPI_Exscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
551 return PMPI_Iexscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
554 int PMPI_Iexscan(const void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request){
556 CHECK_BUFFER(1, sendbuf, count)
557 CHECK_BUFFER(2, recvbuf, count)
558 CHECK_TYPE(4, datatype)
559 CHECK_COUNT(3, count)
564 int rank = simgrid::s4u::this_actor::get_pid();
565 std::vector<unsigned char> tmp_sendbuf;
566 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, count, datatype);
568 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Exscan" : "PMPI_Iexscan",
569 new simgrid::instr::Pt2PtTIData(request == MPI_REQUEST_IGNORED ? "exscan" : "iexscan", -1,
570 datatype->is_replayable() ? count : count * datatype->size(),
571 simgrid::smpi::Datatype::encode(datatype)));
574 if (request == MPI_REQUEST_IGNORED)
575 retval = simgrid::smpi::colls::exscan(real_sendbuf, recvbuf, count, datatype, op, comm);
577 retval = simgrid::smpi::colls::iexscan(real_sendbuf, recvbuf, count, datatype, op, comm, request);
579 TRACE_smpi_comm_out(rank);
584 int PMPI_Reduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
586 return PMPI_Ireduce_scatter(sendbuf, recvbuf, recvcounts, datatype, op, comm, MPI_REQUEST_IGNORED);
589 int PMPI_Ireduce_scatter(const void *sendbuf, void *recvbuf, const int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
592 CHECK_TYPE(4, datatype)
593 CHECK_NULL(3, MPI_ERR_COUNT, recvcounts)
596 for (int i = 0; i < comm->size(); i++) {
597 CHECK_COUNT(3, recvcounts[i])
598 CHECK_BUFFER(1, sendbuf, recvcounts[i])
599 CHECK_BUFFER(2, recvbuf, recvcounts[i])
603 int rank = simgrid::s4u::this_actor::get_pid();
604 auto* trace_recvcounts = new std::vector<int>();
605 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
608 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
609 trace_recvcounts->push_back(recvcounts[i] * dt_send_size);
610 totalcount += recvcounts[i];
612 std::vector<unsigned char> tmp_sendbuf;
613 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, totalcount, datatype);
615 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter" : "PMPI_Ireduce_scatter",
616 new simgrid::instr::VarCollTIData(
617 request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, dt_send_size, nullptr,
618 -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
620 if (request == MPI_REQUEST_IGNORED)
621 simgrid::smpi::colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm);
623 simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts, datatype, op, comm, request);
625 TRACE_smpi_comm_out(rank);
630 int PMPI_Reduce_scatter_block(const void *sendbuf, void *recvbuf, int recvcount,
631 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
633 return PMPI_Ireduce_scatter_block(sendbuf, recvbuf, recvcount, datatype, op, comm, MPI_REQUEST_IGNORED);
636 int PMPI_Ireduce_scatter_block(const void* sendbuf, void* recvbuf, int recvcount, MPI_Datatype datatype, MPI_Op op,
637 MPI_Comm comm, MPI_Request* request)
640 CHECK_BUFFER(1, sendbuf, recvcount)
641 CHECK_BUFFER(2, recvbuf, recvcount)
642 CHECK_TYPE(4, datatype)
643 CHECK_COUNT(3, recvcount)
648 int count = comm->size();
650 int rank = simgrid::s4u::this_actor::get_pid();
651 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
652 auto* trace_recvcounts = new std::vector<int>(recvcount * dt_send_size); // copy data to avoid bad free
653 std::vector<unsigned char> tmp_sendbuf;
654 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, recvcount * count, datatype);
657 rank, request == MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter_block" : "PMPI_Ireduce_scatter_block",
658 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "reducescatter" : "ireducescatter", -1, 0,
659 nullptr, -1, trace_recvcounts, simgrid::smpi::Datatype::encode(datatype), ""));
661 std::vector<int> recvcounts(count);
662 for (int i = 0; i < count; i++)
663 recvcounts[i] = recvcount;
664 if (request == MPI_REQUEST_IGNORED)
665 simgrid::smpi::colls::reduce_scatter(real_sendbuf, recvbuf, recvcounts.data(), datatype, op, comm);
667 simgrid::smpi::colls::ireduce_scatter(real_sendbuf, recvbuf, recvcounts.data(), datatype, op, comm, request);
669 TRACE_smpi_comm_out(rank);
674 int PMPI_Alltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
675 MPI_Datatype recvtype, MPI_Comm comm){
676 return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
679 int PMPI_Ialltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
680 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
683 CHECK_BUFFER(1, sendbuf, sendcount)
684 CHECK_BUFFER(4, recvbuf, recvcount)
685 if(sendbuf != MPI_IN_PLACE)
686 CHECK_TYPE(3, sendtype)
687 CHECK_TYPE(6, recvtype)
688 CHECK_COUNT(5, recvcount)
689 if(sendbuf != MPI_IN_PLACE)
690 CHECK_COUNT(2, sendcount)
691 CHECK_COUNT(5, recvcount)
694 int rank = simgrid::s4u::this_actor::get_pid();
695 int real_sendcount = sendcount;
696 MPI_Datatype real_sendtype = sendtype;
698 std::vector<unsigned char> tmp_sendbuf;
699 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, recvcount * comm->size(), recvtype);
701 if (sendbuf == MPI_IN_PLACE) {
702 real_sendcount = recvcount;
703 real_sendtype = recvtype;
706 if(recvtype->size() * recvcount != real_sendtype->size() * real_sendcount){
707 XBT_WARN("MPI_(I)Alltoall : receive size from each process differs from sent size");
708 return MPI_ERR_TRUNCATE;
713 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoall" : "PMPI_Ialltoall",
714 new simgrid::instr::CollTIData(
715 request == MPI_REQUEST_IGNORED ? "alltoall" : "ialltoall", -1, -1.0,
716 real_sendtype->is_replayable() ? real_sendcount : real_sendcount * real_sendtype->size(),
717 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
718 simgrid::smpi::Datatype::encode(real_sendtype), simgrid::smpi::Datatype::encode(recvtype)));
720 if (request == MPI_REQUEST_IGNORED)
722 simgrid::smpi::colls::alltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype, comm);
724 retval = simgrid::smpi::colls::ialltoall(real_sendbuf, real_sendcount, real_sendtype, recvbuf, recvcount, recvtype,
727 TRACE_smpi_comm_out(rank);
732 int PMPI_Alltoallv(const void* sendbuf, const int* sendcounts, const int* senddispls, MPI_Datatype sendtype, void* recvbuf,
733 const int* recvcounts, const int* recvdispls, MPI_Datatype recvtype, MPI_Comm comm)
735 return PMPI_Ialltoallv(sendbuf, sendcounts, senddispls, sendtype, recvbuf, recvcounts, recvdispls, recvtype, comm, MPI_REQUEST_IGNORED);
738 int PMPI_Ialltoallv(const void* sendbuf, const int* sendcounts, const int* senddispls, MPI_Datatype sendtype, void* recvbuf,
739 const int* recvcounts, const int* recvdispls, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
742 if(sendbuf != MPI_IN_PLACE){
743 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
744 CHECK_NULL(3, MPI_ERR_ARG, senddispls)
745 CHECK_TYPE(4, sendtype)
747 CHECK_TYPE(8, recvtype)
748 CHECK_NULL(6, MPI_ERR_COUNT, recvcounts)
749 CHECK_NULL(7, MPI_ERR_ARG, recvdispls)
752 int rank = simgrid::s4u::this_actor::get_pid();
753 int size = comm->size();
754 for (int i = 0; i < size; i++) {
755 if(sendbuf != MPI_IN_PLACE){
756 CHECK_BUFFER(1, sendbuf, sendcounts[i])
757 CHECK_COUNT(2, sendcounts[i])
759 CHECK_BUFFER(5, recvbuf, recvcounts[i])
760 CHECK_COUNT(6, recvcounts[i])
766 auto* trace_sendcounts = new std::vector<int>();
767 auto* trace_recvcounts = new std::vector<int>();
768 int dt_size_recv = recvtype->size();
770 const int* real_sendcounts = sendcounts;
771 const int* real_senddispls = senddispls;
772 MPI_Datatype real_sendtype = sendtype;
774 for (int i = 0; i < size; i++) { // copy data to avoid bad free
775 recv_size += recvcounts[i] * dt_size_recv;
776 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
777 if (((recvdispls[i] + recvcounts[i]) * dt_size_recv) > maxsize)
778 maxsize = (recvdispls[i] + recvcounts[i]) * dt_size_recv;
781 std::vector<unsigned char> tmp_sendbuf;
782 std::vector<int> tmp_sendcounts;
783 std::vector<int> tmp_senddispls;
784 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
785 if (sendbuf == MPI_IN_PLACE) {
786 tmp_sendcounts.assign(recvcounts, recvcounts + size);
787 real_sendcounts = tmp_sendcounts.data();
788 tmp_senddispls.assign(recvdispls, recvdispls + size);
789 real_senddispls = tmp_senddispls.data();
790 real_sendtype = recvtype;
793 if(recvtype->size() * recvcounts[comm->rank()] != real_sendtype->size() * real_sendcounts[comm->rank()]){
794 XBT_WARN("MPI_(I)Alltoallv : receive size from me differs from sent size to me");
796 return MPI_ERR_TRUNCATE;
799 int dt_size_send = real_sendtype->size();
801 for (int i = 0; i < size; i++) { // copy data to avoid bad free
802 send_size += real_sendcounts[i] * dt_size_send;
803 trace_sendcounts->push_back(real_sendcounts[i] * dt_size_send);
806 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallv" : "PMPI_Ialltoallv",
807 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
808 send_size, trace_sendcounts, recv_size, trace_recvcounts,
809 simgrid::smpi::Datatype::encode(real_sendtype),
810 simgrid::smpi::Datatype::encode(recvtype)));
813 if (request == MPI_REQUEST_IGNORED)
814 retval = simgrid::smpi::colls::alltoallv(real_sendbuf, real_sendcounts, real_senddispls, real_sendtype, recvbuf,
815 recvcounts, recvdispls, recvtype, comm);
817 retval = simgrid::smpi::colls::ialltoallv(real_sendbuf, real_sendcounts, real_senddispls, real_sendtype, recvbuf,
818 recvcounts, recvdispls, recvtype, comm, request);
820 TRACE_smpi_comm_out(rank);
825 int PMPI_Alltoallw(const void* sendbuf, const int* sendcounts, const int* senddispls, const MPI_Datatype* sendtypes, void* recvbuf,
826 const int* recvcounts, const int* recvdispls, const MPI_Datatype* recvtypes, MPI_Comm comm)
828 return PMPI_Ialltoallw(sendbuf, sendcounts, senddispls, sendtypes, recvbuf, recvcounts, recvdispls, recvtypes, comm, MPI_REQUEST_IGNORED);
831 int PMPI_Ialltoallw(const void* sendbuf, const int* sendcounts, const int* senddispls, const MPI_Datatype* sendtypes, void* recvbuf,
832 const int* recvcounts, const int* recvdispls, const MPI_Datatype* recvtypes, MPI_Comm comm, MPI_Request* request)
835 if(sendbuf != MPI_IN_PLACE){
836 CHECK_NULL(2, MPI_ERR_COUNT, sendcounts)
837 CHECK_NULL(3, MPI_ERR_ARG, senddispls)
838 CHECK_NULL(4, MPI_ERR_TYPE, sendtypes)
840 CHECK_NULL(6, MPI_ERR_COUNT, recvcounts)
841 CHECK_NULL(7, MPI_ERR_ARG, recvdispls)
842 CHECK_NULL(8, MPI_ERR_TYPE, recvtypes)
844 int rank = simgrid::s4u::this_actor::get_pid();
845 int size = comm->size();
846 for (int i = 0; i < size; i++) {
847 if(sendbuf != MPI_IN_PLACE){
848 CHECK_BUFFER(1, sendbuf, sendcounts[i])
849 CHECK_COUNT(2, sendcounts[i])
850 CHECK_TYPE(4, sendtypes[i])
852 CHECK_BUFFER(5, recvbuf, recvcounts[i])
853 CHECK_COUNT(6, recvcounts[i])
854 CHECK_TYPE(8, recvtypes[i])
861 auto* trace_sendcounts = new std::vector<int>();
862 auto* trace_recvcounts = new std::vector<int>();
864 const int* real_sendcounts = sendcounts;
865 const int* real_senddispls = senddispls;
866 const MPI_Datatype* real_sendtypes = sendtypes;
867 unsigned long maxsize = 0;
868 for (int i = 0; i < size; i++) { // copy data to avoid bad free
869 if (recvtypes[i] == MPI_DATATYPE_NULL) {
870 delete trace_recvcounts;
871 delete trace_sendcounts;
874 recv_size += recvcounts[i] * recvtypes[i]->size();
875 trace_recvcounts->push_back(recvcounts[i] * recvtypes[i]->size());
876 if ((recvdispls[i] + (recvcounts[i] * recvtypes[i]->size())) > maxsize)
877 maxsize = recvdispls[i] + (recvcounts[i] * recvtypes[i]->size());
880 std::vector<unsigned char> tmp_sendbuf;
881 std::vector<int> tmp_sendcounts;
882 std::vector<int> tmp_senddispls;
883 std::vector<MPI_Datatype> tmp_sendtypes;
884 const void* real_sendbuf = smpi_get_in_place_buf(sendbuf, recvbuf, tmp_sendbuf, maxsize, MPI_CHAR);
885 if (sendbuf == MPI_IN_PLACE) {
886 tmp_sendcounts.assign(recvcounts, recvcounts + size);
887 real_sendcounts = tmp_sendcounts.data();
888 tmp_senddispls.assign(recvdispls, recvdispls + size);
889 real_senddispls = tmp_senddispls.data();
890 tmp_sendtypes.assign(recvtypes, recvtypes + size);
891 real_sendtypes = tmp_sendtypes.data();
895 if(recvtypes[comm->rank()]->size() * recvcounts[comm->rank()] != real_sendtypes[comm->rank()]->size() * real_sendcounts[comm->rank()]){
896 XBT_WARN("MPI_(I)Alltoallw : receive size from me differs from sent size to me");
898 return MPI_ERR_TRUNCATE;
901 for (int i = 0; i < size; i++) { // copy data to avoid bad free
902 send_size += real_sendcounts[i] * real_sendtypes[i]->size();
903 trace_sendcounts->push_back(real_sendcounts[i] * real_sendtypes[i]->size());
906 TRACE_smpi_comm_in(rank, request == MPI_REQUEST_IGNORED ? "PMPI_Alltoallw" : "PMPI_Ialltoallw",
907 new simgrid::instr::VarCollTIData(request == MPI_REQUEST_IGNORED ? "alltoallv" : "ialltoallv", -1,
908 send_size, trace_sendcounts, recv_size, trace_recvcounts,
909 simgrid::smpi::Datatype::encode(real_sendtypes[0]),
910 simgrid::smpi::Datatype::encode(recvtypes[0])));
913 if (request == MPI_REQUEST_IGNORED)
914 retval = simgrid::smpi::colls::alltoallw(real_sendbuf, real_sendcounts, real_senddispls, real_sendtypes, recvbuf,
915 recvcounts, recvdispls, recvtypes, comm);
917 retval = simgrid::smpi::colls::ialltoallw(real_sendbuf, real_sendcounts, real_senddispls, real_sendtypes, recvbuf,
918 recvcounts, recvdispls, recvtypes, comm, request);
920 TRACE_smpi_comm_out(rank);