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 /* PMPI User level calls */
18 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
20 return PMPI_Ibcast(buf, count, datatype, root, comm, MPI_REQUEST_IGNORED);
23 int PMPI_Barrier(MPI_Comm comm)
25 return PMPI_Ibarrier(comm, MPI_REQUEST_IGNORED);
28 int PMPI_Ibarrier(MPI_Comm comm, MPI_Request *request)
32 if (comm == MPI_COMM_NULL) {
33 retval = MPI_ERR_COMM;
34 } else if(request == nullptr){
37 int rank = simgrid::s4u::this_actor::get_pid();
38 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED? "PMPI_Barrier" : "PMPI_Ibarrier", new simgrid::instr::NoOpTIData(request==MPI_REQUEST_IGNORED? "barrier" : "ibarrier"));
39 if(request==MPI_REQUEST_IGNORED){
40 simgrid::smpi::Colls::barrier(comm);
41 //Barrier can be used to synchronize RMA calls. Finish all requests from comm before.
42 comm->finish_rma_calls();
44 simgrid::smpi::Colls::ibarrier(comm, request);
45 TRACE_smpi_comm_out(rank);
51 int PMPI_Ibcast(void *buf, int count, MPI_Datatype datatype,
52 int root, MPI_Comm comm, MPI_Request* request)
56 if (comm == MPI_COMM_NULL) {
57 retval = MPI_ERR_COMM;
58 } else if (not datatype->is_valid()) {
60 } else if(request == nullptr){
63 int rank = simgrid::s4u::this_actor::get_pid();
64 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Bcast":"PMPI_Ibcast",
65 new simgrid::instr::CollTIData(request==MPI_REQUEST_IGNORED?"bcast":"ibcast", root, -1.0,
66 datatype->is_replayable() ? count : count * datatype->size(), -1,
67 simgrid::smpi::Datatype::encode(datatype), ""));
68 if (comm->size() > 1){
69 if(request==MPI_REQUEST_IGNORED)
70 simgrid::smpi::Colls::bcast(buf, count, datatype, root, comm);
72 simgrid::smpi::Colls::ibcast(buf, count, datatype, root, comm, request);
74 if(request!=MPI_REQUEST_IGNORED)
75 *request = MPI_REQUEST_NULL;
79 TRACE_smpi_comm_out(rank);
85 int PMPI_Gather(void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
86 int root, MPI_Comm comm){
87 return PMPI_Igather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
90 int PMPI_Igather(void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
91 int root, MPI_Comm comm, MPI_Request *request)
97 if (comm == MPI_COMM_NULL) {
98 retval = MPI_ERR_COMM;
99 } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
100 ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL))){
101 retval = MPI_ERR_TYPE;
102 } else if ((( sendbuf != MPI_IN_PLACE) && (sendcount <0)) || ((comm->rank() == root) && (recvcount <0))){
103 retval = MPI_ERR_COUNT;
104 } else if (request == nullptr){
105 retval = MPI_ERR_ARG;
108 char* sendtmpbuf = static_cast<char*>(sendbuf);
109 int sendtmpcount = sendcount;
110 MPI_Datatype sendtmptype = sendtype;
111 if( (comm->rank() == root) && (sendbuf == MPI_IN_PLACE )) {
113 sendtmptype=recvtype;
115 int rank = simgrid::s4u::this_actor::get_pid();
118 rank, request==MPI_REQUEST_IGNORED?"PMPI_Gather":"PMPI_Igather",
119 new simgrid::instr::CollTIData(
120 request==MPI_REQUEST_IGNORED ? "gather":"igather", root, -1.0, sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(),
121 (comm->rank() != root || recvtype->is_replayable()) ? recvcount : recvcount * recvtype->size(),
122 simgrid::smpi::Datatype::encode(sendtmptype), simgrid::smpi::Datatype::encode(recvtype)));
123 if(request == MPI_REQUEST_IGNORED)
124 simgrid::smpi::Colls::gather(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, root, comm);
126 simgrid::smpi::Colls::igather(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, root, comm, request);
128 retval = MPI_SUCCESS;
129 TRACE_smpi_comm_out(rank);
136 int PMPI_Gatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, int *recvcounts, int *displs,
137 MPI_Datatype recvtype, int root, MPI_Comm comm){
138 return PMPI_Igatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, root, comm, MPI_REQUEST_IGNORED);
141 int PMPI_Igatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, int *recvcounts, int *displs,
142 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request *request)
148 if (comm == MPI_COMM_NULL) {
149 retval = MPI_ERR_COMM;
150 } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
151 ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL))){
152 retval = MPI_ERR_TYPE;
153 } else if (( sendbuf != MPI_IN_PLACE) && (sendcount <0)){
154 retval = MPI_ERR_COUNT;
155 } else if ((comm->rank() == root) && (recvcounts == nullptr || displs == nullptr)) {
156 retval = MPI_ERR_ARG;
157 } else if (request == nullptr){
158 retval = MPI_ERR_ARG;
160 char* sendtmpbuf = static_cast<char*>(sendbuf);
161 int sendtmpcount = sendcount;
162 MPI_Datatype sendtmptype = sendtype;
163 if( (comm->rank() == root) && (sendbuf == MPI_IN_PLACE )) {
165 sendtmptype=recvtype;
168 int rank = simgrid::s4u::this_actor::get_pid();
169 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
171 std::vector<int>* trace_recvcounts = new std::vector<int>;
172 if (comm->rank() == root) {
173 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
174 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
177 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Gatherv":"PMPI_Igatherv",
178 new simgrid::instr::VarCollTIData(
179 request==MPI_REQUEST_IGNORED ? "gatherv":"igatherv", root,
180 sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(), nullptr,
181 dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtmptype),
182 simgrid::smpi::Datatype::encode(recvtype)));
183 if(request == MPI_REQUEST_IGNORED)
184 retval = simgrid::smpi::Colls::gatherv(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcounts, displs, recvtype, root, comm);
186 retval = simgrid::smpi::Colls::igatherv(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcounts, displs, recvtype, root, comm, request);
187 TRACE_smpi_comm_out(rank);
194 int PMPI_Allgather(void *sendbuf, int sendcount, MPI_Datatype sendtype,
195 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm){
196 return PMPI_Iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
199 int PMPI_Iallgather(void *sendbuf, int sendcount, MPI_Datatype sendtype,
200 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
202 int retval = MPI_SUCCESS;
206 if (comm == MPI_COMM_NULL) {
207 retval = MPI_ERR_COMM;
208 } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
209 (recvtype == MPI_DATATYPE_NULL)){
210 retval = MPI_ERR_TYPE;
211 } else if ((( sendbuf != MPI_IN_PLACE) && (sendcount <0)) ||
213 retval = MPI_ERR_COUNT;
214 } else if (request == nullptr){
215 retval = MPI_ERR_ARG;
217 if(sendbuf == MPI_IN_PLACE) {
218 sendbuf=static_cast<char*>(recvbuf)+recvtype->get_extent()*recvcount*comm->rank();
222 int rank = simgrid::s4u::this_actor::get_pid();
224 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Allgather":"PMPI_Iallggather",
225 new simgrid::instr::CollTIData(
226 request==MPI_REQUEST_IGNORED ? "allgather" : "iallgather", -1, -1.0, sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
227 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
228 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
229 if(request == MPI_REQUEST_IGNORED)
230 simgrid::smpi::Colls::allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
232 simgrid::smpi::Colls::iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, request);
233 TRACE_smpi_comm_out(rank);
239 int PMPI_Allgatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype,
240 void *recvbuf, int *recvcounts, int *displs, MPI_Datatype recvtype, MPI_Comm comm){
241 return PMPI_Iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, MPI_REQUEST_IGNORED);
244 int PMPI_Iallgatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype,
245 void *recvbuf, int *recvcounts, int *displs, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
251 if (comm == MPI_COMM_NULL) {
252 retval = MPI_ERR_COMM;
253 } else if (((sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) || (recvtype == MPI_DATATYPE_NULL)) {
254 retval = MPI_ERR_TYPE;
255 } else if (( sendbuf != MPI_IN_PLACE) && (sendcount <0)){
256 retval = MPI_ERR_COUNT;
257 } else if (recvcounts == nullptr || displs == nullptr) {
258 retval = MPI_ERR_ARG;
259 } else if (request == nullptr){
260 retval = MPI_ERR_ARG;
263 if(sendbuf == MPI_IN_PLACE) {
264 sendbuf=static_cast<char*>(recvbuf)+recvtype->get_extent()*displs[comm->rank()];
265 sendcount=recvcounts[comm->rank()];
268 int rank = simgrid::s4u::this_actor::get_pid();
269 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
271 std::vector<int>* trace_recvcounts = new std::vector<int>;
272 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
273 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
275 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Allgatherv":"PMPI_Iallgatherv",
276 new simgrid::instr::VarCollTIData(
277 request==MPI_REQUEST_IGNORED ? "allgatherv" : "iallgatherv", -1, sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
278 nullptr, dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
279 simgrid::smpi::Datatype::encode(recvtype)));
280 if(request == MPI_REQUEST_IGNORED)
281 simgrid::smpi::Colls::allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
283 simgrid::smpi::Colls::iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, request);
284 retval = MPI_SUCCESS;
285 TRACE_smpi_comm_out(rank);
292 int PMPI_Scatter(void *sendbuf, int sendcount, MPI_Datatype sendtype,
293 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
294 return PMPI_Iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
297 int PMPI_Iscatter(void *sendbuf, int sendcount, MPI_Datatype sendtype,
298 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
304 if (comm == MPI_COMM_NULL) {
305 retval = MPI_ERR_COMM;
306 } else if (((comm->rank() == root) && (not sendtype->is_valid())) ||
307 ((recvbuf != MPI_IN_PLACE) && (not recvtype->is_valid()))) {
308 retval = MPI_ERR_TYPE;
309 } else if ((sendbuf == recvbuf) ||
310 ((comm->rank()==root) && sendcount>0 && (sendbuf == nullptr))){
311 retval = MPI_ERR_BUFFER;
312 }else if (request == nullptr){
313 retval = MPI_ERR_ARG;
316 if (recvbuf == MPI_IN_PLACE) {
318 recvcount = sendcount;
320 int rank = simgrid::s4u::this_actor::get_pid();
323 rank, request==MPI_REQUEST_IGNORED?"PMPI_Scatter":"PMPI_Iscatter",
324 new simgrid::instr::CollTIData(
325 request==MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
326 (comm->rank() != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
327 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
328 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
329 if(request == MPI_REQUEST_IGNORED)
330 simgrid::smpi::Colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
332 simgrid::smpi::Colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
333 retval = MPI_SUCCESS;
334 TRACE_smpi_comm_out(rank);
341 int PMPI_Scatterv(void *sendbuf, int *sendcounts, int *displs,
342 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
343 return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
346 int PMPI_Iscatterv(void *sendbuf, int *sendcounts, int *displs,
347 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request *request)
353 if (comm == MPI_COMM_NULL) {
354 retval = MPI_ERR_COMM;
355 } else if (sendcounts == nullptr || displs == nullptr) {
356 retval = MPI_ERR_ARG;
357 } else if (((comm->rank() == root) && (sendtype == MPI_DATATYPE_NULL)) ||
358 ((recvbuf != MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL))) {
359 retval = MPI_ERR_TYPE;
360 } else if (request == nullptr){
361 retval = MPI_ERR_ARG;
363 if (recvbuf == MPI_IN_PLACE) {
365 recvcount = sendcounts[comm->rank()];
367 int rank = simgrid::s4u::this_actor::get_pid();
368 int dt_size_send = sendtype->is_replayable() ? 1 : sendtype->size();
370 std::vector<int>* trace_sendcounts = new std::vector<int>;
371 if (comm->rank() == root) {
372 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
373 trace_sendcounts->push_back(sendcounts[i] * dt_size_send);
376 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Scatterv":"PMPI_Iscatterv",
377 new simgrid::instr::VarCollTIData(
378 request==MPI_REQUEST_IGNORED ? "scatterv":"iscatterv", root, dt_size_send, trace_sendcounts,
379 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(), nullptr,
380 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
381 if(request == MPI_REQUEST_IGNORED)
382 retval = simgrid::smpi::Colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
384 retval = simgrid::smpi::Colls::iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
386 TRACE_smpi_comm_out(rank);
393 int PMPI_Reduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
395 return PMPI_Ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, MPI_REQUEST_IGNORED);
398 int PMPI_Ireduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm, MPI_Request* request)
404 if (comm == MPI_COMM_NULL) {
405 retval = MPI_ERR_COMM;
406 } else if (not datatype->is_valid() || op == MPI_OP_NULL) {
407 retval = MPI_ERR_ARG;
408 } else if (request == nullptr){
409 retval = MPI_ERR_ARG;
411 int rank = simgrid::s4u::this_actor::get_pid();
413 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Reduce":"PMPI_Ireduce",
414 new simgrid::instr::CollTIData(request==MPI_REQUEST_IGNORED ? "reduce":"ireduce", root, 0,
415 datatype->is_replayable() ? count : count * datatype->size(), -1,
416 simgrid::smpi::Datatype::encode(datatype), ""));
417 if(request == MPI_REQUEST_IGNORED)
418 simgrid::smpi::Colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
420 simgrid::smpi::Colls::ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, request);
423 retval = MPI_SUCCESS;
424 TRACE_smpi_comm_out(rank);
431 int PMPI_Reduce_local(void *inbuf, void *inoutbuf, int count, MPI_Datatype datatype, MPI_Op op){
435 if (not datatype->is_valid() || op == MPI_OP_NULL) {
436 retval = MPI_ERR_ARG;
438 op->apply(inbuf, inoutbuf, &count, datatype);
439 retval = MPI_SUCCESS;
445 int PMPI_Allreduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
447 return PMPI_Iallreduce(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
450 int PMPI_Iallreduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
456 if (comm == MPI_COMM_NULL) {
457 retval = MPI_ERR_COMM;
458 } else if (not datatype->is_valid()) {
459 retval = MPI_ERR_TYPE;
460 } else if (op == MPI_OP_NULL) {
462 } else if (request == nullptr){
463 retval = MPI_ERR_ARG;
465 char* sendtmpbuf = static_cast<char*>(sendbuf);
466 if( sendbuf == MPI_IN_PLACE ) {
467 sendtmpbuf = static_cast<char*>(xbt_malloc(count*datatype->get_extent()));
468 simgrid::smpi::Datatype::copy(recvbuf, count, datatype,sendtmpbuf, count, datatype);
470 int rank = simgrid::s4u::this_actor::get_pid();
472 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ?"PMPI_Allreduce":"PMPI_Iallreduce",
473 new simgrid::instr::CollTIData(request==MPI_REQUEST_IGNORED ? "allreduce":"iallreduce", -1, 0,
474 datatype->is_replayable() ? count : count * datatype->size(), -1,
475 simgrid::smpi::Datatype::encode(datatype), ""));
477 if(request == MPI_REQUEST_IGNORED)
478 simgrid::smpi::Colls::allreduce(sendtmpbuf, recvbuf, count, datatype, op, comm);
480 simgrid::smpi::Colls::iallreduce(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
482 if( sendbuf == MPI_IN_PLACE )
483 xbt_free(sendtmpbuf);
485 retval = MPI_SUCCESS;
486 TRACE_smpi_comm_out(rank);
493 int PMPI_Scan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
495 return PMPI_Iscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
498 int PMPI_Iscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request)
504 if (comm == MPI_COMM_NULL) {
505 retval = MPI_ERR_COMM;
506 } else if (not datatype->is_valid()) {
507 retval = MPI_ERR_TYPE;
508 } else if (op == MPI_OP_NULL) {
510 } else if (request == nullptr){
511 retval = MPI_ERR_ARG;
513 int rank = simgrid::s4u::this_actor::get_pid();
514 void* sendtmpbuf = sendbuf;
515 if (sendbuf == MPI_IN_PLACE) {
516 sendtmpbuf = static_cast<void*>(xbt_malloc(count * datatype->size()));
517 memcpy(sendtmpbuf, recvbuf, count * datatype->size());
519 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Scan" : "PMPI_Iscan", new simgrid::instr::Pt2PtTIData(
520 request==MPI_REQUEST_IGNORED ? "scan":"iscan", -1,
521 datatype->is_replayable() ? count : count * datatype->size(),
522 simgrid::smpi::Datatype::encode(datatype)));
524 if(request == MPI_REQUEST_IGNORED)
525 retval = simgrid::smpi::Colls::scan(sendtmpbuf, recvbuf, count, datatype, op, comm);
527 retval = simgrid::smpi::Colls::iscan(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
529 TRACE_smpi_comm_out(rank);
530 if (sendbuf == MPI_IN_PLACE)
531 xbt_free(sendtmpbuf);
538 int PMPI_Exscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
540 return PMPI_Iexscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
543 int PMPI_Iexscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request){
548 if (comm == MPI_COMM_NULL) {
549 retval = MPI_ERR_COMM;
550 } else if (not datatype->is_valid()) {
551 retval = MPI_ERR_TYPE;
552 } else if (op == MPI_OP_NULL) {
554 } else if (request == nullptr){
555 retval = MPI_ERR_ARG;
557 int rank = simgrid::s4u::this_actor::get_pid();
558 void* sendtmpbuf = sendbuf;
559 if (sendbuf == MPI_IN_PLACE) {
560 sendtmpbuf = static_cast<void*>(xbt_malloc(count * datatype->size()));
561 memcpy(sendtmpbuf, recvbuf, count * datatype->size());
564 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Exscan" : "PMPI_Iexscan", new simgrid::instr::Pt2PtTIData(
565 request==MPI_REQUEST_IGNORED ? "exscan":"iexscan", -1, datatype->is_replayable() ? count : count * datatype->size(),
566 simgrid::smpi::Datatype::encode(datatype)));
568 if(request == MPI_REQUEST_IGNORED)
569 retval = simgrid::smpi::Colls::exscan(sendtmpbuf, recvbuf, count, datatype, op, comm);
571 retval = simgrid::smpi::Colls::iexscan(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
573 TRACE_smpi_comm_out(rank);
574 if (sendbuf == MPI_IN_PLACE)
575 xbt_free(sendtmpbuf);
582 int PMPI_Reduce_scatter(void *sendbuf, void *recvbuf, int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
584 return PMPI_Ireduce_scatter(sendbuf, recvbuf, recvcounts, datatype, op, comm, MPI_REQUEST_IGNORED);
587 int PMPI_Ireduce_scatter(void *sendbuf, void *recvbuf, int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
592 if (comm == MPI_COMM_NULL) {
593 retval = MPI_ERR_COMM;
594 } else if (not datatype->is_valid()) {
595 retval = MPI_ERR_TYPE;
596 } else if (op == MPI_OP_NULL) {
598 } else if (recvcounts == nullptr) {
599 retval = MPI_ERR_ARG;
600 } else if (request == nullptr){
601 retval = MPI_ERR_ARG;
603 int rank = simgrid::s4u::this_actor::get_pid();
604 std::vector<int>* 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];
613 void* sendtmpbuf = sendbuf;
614 if (sendbuf == MPI_IN_PLACE) {
615 sendtmpbuf = static_cast<void*>(xbt_malloc(totalcount * datatype->size()));
616 memcpy(sendtmpbuf, recvbuf, totalcount * datatype->size());
619 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter": "PMPI_Ireduce_scatter", new simgrid::instr::VarCollTIData(
620 request==MPI_REQUEST_IGNORED ? "reducescatter":"ireducescatter", -1, dt_send_size, nullptr, -1, trace_recvcounts,
621 simgrid::smpi::Datatype::encode(datatype), ""));
623 if(request == MPI_REQUEST_IGNORED)
624 simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
626 simgrid::smpi::Colls::ireduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm, request);
628 retval = MPI_SUCCESS;
629 TRACE_smpi_comm_out(rank);
631 if (sendbuf == MPI_IN_PLACE)
632 xbt_free(sendtmpbuf);
639 int PMPI_Reduce_scatter_block(void *sendbuf, void *recvbuf, int recvcount,
640 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
642 return PMPI_Ireduce_scatter_block(sendbuf, recvbuf, recvcount, datatype, op, comm, MPI_REQUEST_IGNORED);
645 int PMPI_Ireduce_scatter_block(void *sendbuf, void *recvbuf, int recvcount,
646 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
651 if (comm == MPI_COMM_NULL) {
652 retval = MPI_ERR_COMM;
653 } else if (not datatype->is_valid()) {
654 retval = MPI_ERR_TYPE;
655 } else if (op == MPI_OP_NULL) {
657 } else if (recvcount < 0) {
658 retval = MPI_ERR_ARG;
659 } else if (request == nullptr){
660 retval = MPI_ERR_ARG;
662 int count = comm->size();
664 int rank = simgrid::s4u::this_actor::get_pid();
665 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
666 std::vector<int>* trace_recvcounts = new std::vector<int>(recvcount * dt_send_size); // copy data to avoid bad free
668 void* sendtmpbuf = sendbuf;
669 if (sendbuf == MPI_IN_PLACE) {
670 sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * count * datatype->size()));
671 memcpy(sendtmpbuf, recvbuf, recvcount * count * datatype->size());
674 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter_block":"PMPI_Ireduce_scatter_block",
675 new simgrid::instr::VarCollTIData(request==MPI_REQUEST_IGNORED ? "reducescatter":"ireducescatter", -1, 0, nullptr, -1, trace_recvcounts,
676 simgrid::smpi::Datatype::encode(datatype), ""));
678 int* recvcounts = new int[count];
679 for (int i = 0; i < count; i++)
680 recvcounts[i] = recvcount;
681 if(request == MPI_REQUEST_IGNORED)
682 simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
684 simgrid::smpi::Colls::ireduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm, request);
686 retval = MPI_SUCCESS;
688 TRACE_smpi_comm_out(rank);
690 if (sendbuf == MPI_IN_PLACE)
691 xbt_free(sendtmpbuf);
697 int PMPI_Alltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
698 MPI_Datatype recvtype, MPI_Comm comm){
699 return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
702 int PMPI_Ialltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
703 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request *request)
708 if (comm == MPI_COMM_NULL) {
709 retval = MPI_ERR_COMM;
710 } else if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL) {
711 retval = MPI_ERR_TYPE;
712 } else if (request == nullptr){
713 retval = MPI_ERR_ARG;
715 int rank = simgrid::s4u::this_actor::get_pid();
716 void* sendtmpbuf = static_cast<char*>(sendbuf);
717 int sendtmpcount = sendcount;
718 MPI_Datatype sendtmptype = sendtype;
719 if (sendbuf == MPI_IN_PLACE) {
720 sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * comm->size() * recvtype->size()));
721 memcpy(sendtmpbuf, recvbuf, recvcount * comm->size() * recvtype->size());
722 sendtmpcount = recvcount;
723 sendtmptype = recvtype;
726 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Alltoall":"PMPI_Ialltoall",
727 new simgrid::instr::CollTIData(
728 request==MPI_REQUEST_IGNORED ? "alltoall" : "ialltoall", -1, -1.0,
729 sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(),
730 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
731 simgrid::smpi::Datatype::encode(sendtmptype), simgrid::smpi::Datatype::encode(recvtype)));
732 if(request == MPI_REQUEST_IGNORED)
733 retval = simgrid::smpi::Colls::alltoall(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, comm);
735 retval = simgrid::smpi::Colls::ialltoall(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, comm, request);
737 TRACE_smpi_comm_out(rank);
739 if (sendbuf == MPI_IN_PLACE)
740 xbt_free(sendtmpbuf);
747 int PMPI_Alltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf,
748 int* recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm)
750 return PMPI_Ialltoallv(sendbuf, sendcounts, senddisps, sendtype, recvbuf, recvcounts, recvdisps, recvtype, comm, MPI_REQUEST_IGNORED);
753 int PMPI_Ialltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf,
754 int* recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request *request)
760 if (comm == MPI_COMM_NULL) {
761 retval = MPI_ERR_COMM;
762 } else if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL) {
763 retval = MPI_ERR_TYPE;
764 } else if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
765 recvdisps == nullptr) {
766 retval = MPI_ERR_ARG;
767 } else if (request == nullptr){
768 retval = MPI_ERR_ARG;
770 int rank = simgrid::s4u::this_actor::get_pid();
771 int size = comm->size();
774 std::vector<int>* trace_sendcounts = new std::vector<int>;
775 std::vector<int>* trace_recvcounts = new std::vector<int>;
776 int dt_size_recv = recvtype->size();
778 void* sendtmpbuf = static_cast<char*>(sendbuf);
779 int* sendtmpcounts = sendcounts;
780 int* sendtmpdisps = senddisps;
781 MPI_Datatype sendtmptype = sendtype;
783 for (int i = 0; i < size; i++) { // copy data to avoid bad free
784 recv_size += recvcounts[i] * dt_size_recv;
785 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
786 if (((recvdisps[i] + recvcounts[i]) * dt_size_recv) > maxsize)
787 maxsize = (recvdisps[i] + recvcounts[i]) * dt_size_recv;
790 if (sendbuf == MPI_IN_PLACE) {
791 sendtmpbuf = static_cast<void*>(xbt_malloc(maxsize));
792 memcpy(sendtmpbuf, recvbuf, maxsize);
793 sendtmpcounts = static_cast<int*>(xbt_malloc(size * sizeof(int)));
794 memcpy(sendtmpcounts, recvcounts, size * sizeof(int));
795 sendtmpdisps = static_cast<int*>(xbt_malloc(size * sizeof(int)));
796 memcpy(sendtmpdisps, recvdisps, size * sizeof(int));
797 sendtmptype = recvtype;
800 int dt_size_send = sendtmptype->size();
802 for (int i = 0; i < size; i++) { // copy data to avoid bad free
803 send_size += sendtmpcounts[i] * dt_size_send;
804 trace_sendcounts->push_back(sendtmpcounts[i] * dt_size_send);
807 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Alltoallv":"PMPI_Ialltoallv",
808 new simgrid::instr::VarCollTIData(request==MPI_REQUEST_IGNORED ? "alltoallv":"ialltoallv", -1, send_size, trace_sendcounts, recv_size,
809 trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
810 simgrid::smpi::Datatype::encode(recvtype)));
812 if(request == MPI_REQUEST_IGNORED)
813 retval = simgrid::smpi::Colls::alltoallv(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptype, recvbuf, recvcounts,
814 recvdisps, recvtype, comm);
816 retval = simgrid::smpi::Colls::ialltoallv(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptype, recvbuf, recvcounts,
817 recvdisps, recvtype, comm, request);
818 TRACE_smpi_comm_out(rank);
820 if (sendbuf == MPI_IN_PLACE) {
821 xbt_free(sendtmpbuf);
822 xbt_free(sendtmpcounts);
823 xbt_free(sendtmpdisps);
831 int PMPI_Alltoallw(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype* sendtypes, void* recvbuf,
832 int* recvcounts, int* recvdisps, MPI_Datatype* recvtypes, MPI_Comm comm)
834 return PMPI_Ialltoallw(sendbuf, sendcounts, senddisps, sendtypes, recvbuf, recvcounts, recvdisps, recvtypes, comm, MPI_REQUEST_IGNORED);
837 int PMPI_Ialltoallw(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype* sendtypes, void* recvbuf,
838 int* recvcounts, int* recvdisps, MPI_Datatype* recvtypes, MPI_Comm comm, MPI_Request *request)
844 if (comm == MPI_COMM_NULL) {
845 retval = MPI_ERR_COMM;
846 } else if ((sendbuf != MPI_IN_PLACE && sendtypes == nullptr) || recvtypes == nullptr) {
847 retval = MPI_ERR_TYPE;
848 } else if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
849 recvdisps == nullptr) {
850 retval = MPI_ERR_ARG;
851 } else if (request == nullptr){
852 retval = MPI_ERR_ARG;
854 int rank = simgrid::s4u::this_actor::get_pid();
855 int size = comm->size();
858 std::vector<int>* trace_sendcounts = new std::vector<int>;
859 std::vector<int>* trace_recvcounts = new std::vector<int>;
861 void* sendtmpbuf = static_cast<char*>(sendbuf);
862 int* sendtmpcounts = sendcounts;
863 int* sendtmpdisps = senddisps;
864 MPI_Datatype* sendtmptypes = sendtypes;
865 unsigned long maxsize = 0;
866 for (int i = 0; i < size; i++) { // copy data to avoid bad free
867 if(recvtypes[i]==MPI_DATATYPE_NULL){
868 delete trace_recvcounts;
869 delete trace_sendcounts;
872 recv_size += recvcounts[i] * recvtypes[i]->size();
873 trace_recvcounts->push_back(recvcounts[i] * recvtypes[i]->size());
874 if ((recvdisps[i] + (recvcounts[i] * recvtypes[i]->size())) > maxsize)
875 maxsize = recvdisps[i] + (recvcounts[i] * recvtypes[i]->size());
878 if (sendbuf == MPI_IN_PLACE) {
879 sendtmpbuf = static_cast<void*>(xbt_malloc(maxsize));
880 memcpy(sendtmpbuf, recvbuf, maxsize);
881 sendtmpcounts = static_cast<int*>(xbt_malloc(size * sizeof(int)));
882 memcpy(sendtmpcounts, recvcounts, size * sizeof(int));
883 sendtmpdisps = static_cast<int*>(xbt_malloc(size * sizeof(int)));
884 memcpy(sendtmpdisps, recvdisps, size * sizeof(int));
885 sendtmptypes = static_cast<MPI_Datatype*>(xbt_malloc(size * sizeof(MPI_Datatype)));
886 memcpy(sendtmptypes, recvtypes, size * sizeof(MPI_Datatype));
889 for (int i = 0; i < size; i++) { // copy data to avoid bad free
890 send_size += sendtmpcounts[i] * sendtmptypes[i]->size();
891 trace_sendcounts->push_back(sendtmpcounts[i] * sendtmptypes[i]->size());
894 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Alltoallw":"PMPI_Ialltoallw",
895 new simgrid::instr::VarCollTIData(request==MPI_REQUEST_IGNORED ? "alltoallv":"ialltoallv", -1, send_size, trace_sendcounts, recv_size,
896 trace_recvcounts, simgrid::smpi::Datatype::encode(sendtmptypes[0]),
897 simgrid::smpi::Datatype::encode(recvtypes[0])));
899 if(request == MPI_REQUEST_IGNORED)
900 retval = simgrid::smpi::Colls::alltoallw(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptypes, recvbuf, recvcounts,
901 recvdisps, recvtypes, comm);
903 retval = simgrid::smpi::Colls::ialltoallw(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptypes, recvbuf, recvcounts,
904 recvdisps, recvtypes, comm, request);
905 TRACE_smpi_comm_out(rank);
907 if (sendbuf == MPI_IN_PLACE) {
908 xbt_free(sendtmpbuf);
909 xbt_free(sendtmpcounts);
910 xbt_free(sendtmpdisps);
911 xbt_free(sendtmptypes);