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;
314 if (recvbuf == MPI_IN_PLACE) {
316 recvcount = sendcount;
318 int rank = simgrid::s4u::this_actor::get_pid();
321 rank, request==MPI_REQUEST_IGNORED?"PMPI_Scatter":"PMPI_Iscatter",
322 new simgrid::instr::CollTIData(
323 request==MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
324 (comm->rank() != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
325 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
326 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
327 if(request == MPI_REQUEST_IGNORED)
328 simgrid::smpi::Colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
330 simgrid::smpi::Colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
331 retval = MPI_SUCCESS;
332 TRACE_smpi_comm_out(rank);
339 int PMPI_Scatterv(void *sendbuf, int *sendcounts, int *displs,
340 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
341 return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
344 int PMPI_Iscatterv(void *sendbuf, int *sendcounts, int *displs,
345 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request *request)
351 if (comm == MPI_COMM_NULL) {
352 retval = MPI_ERR_COMM;
353 } else if (sendcounts == nullptr || displs == nullptr) {
354 retval = MPI_ERR_ARG;
355 } else if (((comm->rank() == root) && (sendtype == MPI_DATATYPE_NULL)) ||
356 ((recvbuf != MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL))) {
357 retval = MPI_ERR_TYPE;
358 } else if (request == nullptr){
359 retval = MPI_ERR_ARG;
361 if (recvbuf == MPI_IN_PLACE) {
363 recvcount = sendcounts[comm->rank()];
365 int rank = simgrid::s4u::this_actor::get_pid();
366 int dt_size_send = sendtype->is_replayable() ? 1 : sendtype->size();
368 std::vector<int>* trace_sendcounts = new std::vector<int>;
369 if (comm->rank() == root) {
370 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
371 trace_sendcounts->push_back(sendcounts[i] * dt_size_send);
374 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Scatterv":"PMPI_Iscatterv",
375 new simgrid::instr::VarCollTIData(
376 request==MPI_REQUEST_IGNORED ? "scatterv":"iscatterv", root, dt_size_send, trace_sendcounts,
377 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(), nullptr,
378 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
379 if(request == MPI_REQUEST_IGNORED)
380 retval = simgrid::smpi::Colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
382 retval = simgrid::smpi::Colls::iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
384 TRACE_smpi_comm_out(rank);
391 int PMPI_Reduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
393 return PMPI_Ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, MPI_REQUEST_IGNORED);
396 int PMPI_Ireduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm, MPI_Request* request)
402 if (comm == MPI_COMM_NULL) {
403 retval = MPI_ERR_COMM;
404 } else if (not datatype->is_valid() || op == MPI_OP_NULL) {
405 retval = MPI_ERR_ARG;
407 int rank = simgrid::s4u::this_actor::get_pid();
409 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Reduce":"PMPI_Ireduce",
410 new simgrid::instr::CollTIData(request==MPI_REQUEST_IGNORED ? "reduce":"ireduce", root, 0,
411 datatype->is_replayable() ? count : count * datatype->size(), -1,
412 simgrid::smpi::Datatype::encode(datatype), ""));
413 if(request == MPI_REQUEST_IGNORED)
414 simgrid::smpi::Colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
416 simgrid::smpi::Colls::ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, request);
419 retval = MPI_SUCCESS;
420 TRACE_smpi_comm_out(rank);
427 int PMPI_Reduce_local(void *inbuf, void *inoutbuf, int count, MPI_Datatype datatype, MPI_Op op){
431 if (not datatype->is_valid() || op == MPI_OP_NULL) {
432 retval = MPI_ERR_ARG;
434 op->apply(inbuf, inoutbuf, &count, datatype);
435 retval = MPI_SUCCESS;
441 int PMPI_Allreduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
443 return PMPI_Iallreduce(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
446 int PMPI_Iallreduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
452 if (comm == MPI_COMM_NULL) {
453 retval = MPI_ERR_COMM;
454 } else if (not datatype->is_valid()) {
455 retval = MPI_ERR_TYPE;
456 } else if (op == MPI_OP_NULL) {
459 char* sendtmpbuf = static_cast<char*>(sendbuf);
460 if( sendbuf == MPI_IN_PLACE ) {
461 sendtmpbuf = static_cast<char*>(xbt_malloc(count*datatype->get_extent()));
462 simgrid::smpi::Datatype::copy(recvbuf, count, datatype,sendtmpbuf, count, datatype);
464 int rank = simgrid::s4u::this_actor::get_pid();
466 TRACE_smpi_comm_in(rank, __func__,
467 new simgrid::instr::CollTIData(request==MPI_REQUEST_IGNORED ? "allreduce":"iallreduce", -1, 0,
468 datatype->is_replayable() ? count : count * datatype->size(), -1,
469 simgrid::smpi::Datatype::encode(datatype), ""));
471 if(request == MPI_REQUEST_IGNORED)
472 simgrid::smpi::Colls::allreduce(sendtmpbuf, recvbuf, count, datatype, op, comm);
474 simgrid::smpi::Colls::iallreduce(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
476 if( sendbuf == MPI_IN_PLACE )
477 xbt_free(sendtmpbuf);
479 retval = MPI_SUCCESS;
480 TRACE_smpi_comm_out(rank);
487 int PMPI_Scan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
493 if (comm == MPI_COMM_NULL) {
494 retval = MPI_ERR_COMM;
495 } else if (not datatype->is_valid()) {
496 retval = MPI_ERR_TYPE;
497 } else if (op == MPI_OP_NULL) {
500 int rank = simgrid::s4u::this_actor::get_pid();
501 void* sendtmpbuf = sendbuf;
502 if (sendbuf == MPI_IN_PLACE) {
503 sendtmpbuf = static_cast<void*>(xbt_malloc(count * datatype->size()));
504 memcpy(sendtmpbuf, recvbuf, count * datatype->size());
506 TRACE_smpi_comm_in(rank, __func__, new simgrid::instr::Pt2PtTIData(
507 "scan", -1, datatype->is_replayable() ? count : count * datatype->size(),
508 simgrid::smpi::Datatype::encode(datatype)));
510 retval = simgrid::smpi::Colls::scan(sendtmpbuf, recvbuf, count, datatype, op, comm);
512 TRACE_smpi_comm_out(rank);
513 if (sendbuf == MPI_IN_PLACE)
514 xbt_free(sendtmpbuf);
521 int PMPI_Exscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm){
526 if (comm == MPI_COMM_NULL) {
527 retval = MPI_ERR_COMM;
528 } else if (not datatype->is_valid()) {
529 retval = MPI_ERR_TYPE;
530 } else if (op == MPI_OP_NULL) {
533 int rank = simgrid::s4u::this_actor::get_pid();
534 void* sendtmpbuf = sendbuf;
535 if (sendbuf == MPI_IN_PLACE) {
536 sendtmpbuf = static_cast<void*>(xbt_malloc(count * datatype->size()));
537 memcpy(sendtmpbuf, recvbuf, count * datatype->size());
540 TRACE_smpi_comm_in(rank, __func__, new simgrid::instr::Pt2PtTIData(
541 "exscan", -1, datatype->is_replayable() ? count : count * datatype->size(),
542 simgrid::smpi::Datatype::encode(datatype)));
544 retval = simgrid::smpi::Colls::exscan(sendtmpbuf, recvbuf, count, datatype, op, comm);
546 TRACE_smpi_comm_out(rank);
547 if (sendbuf == MPI_IN_PLACE)
548 xbt_free(sendtmpbuf);
555 int PMPI_Reduce_scatter(void *sendbuf, void *recvbuf, int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
560 if (comm == MPI_COMM_NULL) {
561 retval = MPI_ERR_COMM;
562 } else if (not datatype->is_valid()) {
563 retval = MPI_ERR_TYPE;
564 } else if (op == MPI_OP_NULL) {
566 } else if (recvcounts == nullptr) {
567 retval = MPI_ERR_ARG;
569 int rank = simgrid::s4u::this_actor::get_pid();
570 std::vector<int>* trace_recvcounts = new std::vector<int>;
571 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
574 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
575 trace_recvcounts->push_back(recvcounts[i] * dt_send_size);
576 totalcount += recvcounts[i];
579 void* sendtmpbuf = sendbuf;
580 if (sendbuf == MPI_IN_PLACE) {
581 sendtmpbuf = static_cast<void*>(xbt_malloc(totalcount * datatype->size()));
582 memcpy(sendtmpbuf, recvbuf, totalcount * datatype->size());
585 TRACE_smpi_comm_in(rank, __func__, new simgrid::instr::VarCollTIData(
586 "reducescatter", -1, dt_send_size, nullptr, -1, trace_recvcounts,
587 simgrid::smpi::Datatype::encode(datatype), ""));
589 simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
590 retval = MPI_SUCCESS;
591 TRACE_smpi_comm_out(rank);
593 if (sendbuf == MPI_IN_PLACE)
594 xbt_free(sendtmpbuf);
601 int PMPI_Reduce_scatter_block(void *sendbuf, void *recvbuf, int recvcount,
602 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
607 if (comm == MPI_COMM_NULL) {
608 retval = MPI_ERR_COMM;
609 } else if (not datatype->is_valid()) {
610 retval = MPI_ERR_TYPE;
611 } else if (op == MPI_OP_NULL) {
613 } else if (recvcount < 0) {
614 retval = MPI_ERR_ARG;
616 int count = comm->size();
618 int rank = simgrid::s4u::this_actor::get_pid();
619 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
620 std::vector<int>* trace_recvcounts = new std::vector<int>(recvcount * dt_send_size); // copy data to avoid bad free
622 void* sendtmpbuf = sendbuf;
623 if (sendbuf == MPI_IN_PLACE) {
624 sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * count * datatype->size()));
625 memcpy(sendtmpbuf, recvbuf, recvcount * count * datatype->size());
628 TRACE_smpi_comm_in(rank, __func__,
629 new simgrid::instr::VarCollTIData("reducescatter", -1, 0, nullptr, -1, trace_recvcounts,
630 simgrid::smpi::Datatype::encode(datatype), ""));
632 int* recvcounts = new int[count];
633 for (int i = 0; i < count; i++)
634 recvcounts[i] = recvcount;
635 simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
637 retval = MPI_SUCCESS;
639 TRACE_smpi_comm_out(rank);
641 if (sendbuf == MPI_IN_PLACE)
642 xbt_free(sendtmpbuf);
648 int PMPI_Alltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
649 MPI_Datatype recvtype, MPI_Comm comm){
650 return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
653 int PMPI_Ialltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
654 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request *request)
659 if (comm == MPI_COMM_NULL) {
660 retval = MPI_ERR_COMM;
661 } else if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL) {
662 retval = MPI_ERR_TYPE;
663 } else if (request == nullptr){
664 retval = MPI_ERR_ARG;
666 int rank = simgrid::s4u::this_actor::get_pid();
667 void* sendtmpbuf = static_cast<char*>(sendbuf);
668 int sendtmpcount = sendcount;
669 MPI_Datatype sendtmptype = sendtype;
670 if (sendbuf == MPI_IN_PLACE) {
671 sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * comm->size() * recvtype->size()));
672 memcpy(sendtmpbuf, recvbuf, recvcount * comm->size() * recvtype->size());
673 sendtmpcount = recvcount;
674 sendtmptype = recvtype;
677 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Alltoall":"PMPI_Ialltoall",
678 new simgrid::instr::CollTIData(
679 request==MPI_REQUEST_IGNORED ? "alltoall" : "ialltoall", -1, -1.0,
680 sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(),
681 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
682 simgrid::smpi::Datatype::encode(sendtmptype), simgrid::smpi::Datatype::encode(recvtype)));
683 if(request == MPI_REQUEST_IGNORED)
684 retval = simgrid::smpi::Colls::alltoall(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, comm);
686 retval = simgrid::smpi::Colls::ialltoall(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, comm, request);
688 TRACE_smpi_comm_out(rank);
690 if (sendbuf == MPI_IN_PLACE)
691 xbt_free(sendtmpbuf);
698 int PMPI_Alltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf,
699 int* recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm)
701 return PMPI_Ialltoallv(sendbuf, sendcounts, senddisps, sendtype, recvbuf, recvcounts, recvdisps, recvtype, comm, MPI_REQUEST_IGNORED);
704 int PMPI_Ialltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf,
705 int* recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request *request)
711 if (comm == MPI_COMM_NULL) {
712 retval = MPI_ERR_COMM;
713 } else if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL) {
714 retval = MPI_ERR_TYPE;
715 } else if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
716 recvdisps == nullptr) {
717 retval = MPI_ERR_ARG;
718 } else if (request == nullptr){
719 retval = MPI_ERR_ARG;
721 int rank = simgrid::s4u::this_actor::get_pid();
722 int size = comm->size();
725 std::vector<int>* trace_sendcounts = new std::vector<int>;
726 std::vector<int>* trace_recvcounts = new std::vector<int>;
727 int dt_size_recv = recvtype->size();
729 void* sendtmpbuf = static_cast<char*>(sendbuf);
730 int* sendtmpcounts = sendcounts;
731 int* sendtmpdisps = senddisps;
732 MPI_Datatype sendtmptype = sendtype;
734 for (int i = 0; i < size; i++) { // copy data to avoid bad free
735 recv_size += recvcounts[i] * dt_size_recv;
736 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
737 if (((recvdisps[i] + recvcounts[i]) * dt_size_recv) > maxsize)
738 maxsize = (recvdisps[i] + recvcounts[i]) * dt_size_recv;
741 if (sendbuf == MPI_IN_PLACE) {
742 sendtmpbuf = static_cast<void*>(xbt_malloc(maxsize));
743 memcpy(sendtmpbuf, recvbuf, maxsize);
744 sendtmpcounts = static_cast<int*>(xbt_malloc(size * sizeof(int)));
745 memcpy(sendtmpcounts, recvcounts, size * sizeof(int));
746 sendtmpdisps = static_cast<int*>(xbt_malloc(size * sizeof(int)));
747 memcpy(sendtmpdisps, recvdisps, size * sizeof(int));
748 sendtmptype = recvtype;
751 int dt_size_send = sendtmptype->size();
753 for (int i = 0; i < size; i++) { // copy data to avoid bad free
754 send_size += sendtmpcounts[i] * dt_size_send;
755 trace_sendcounts->push_back(sendtmpcounts[i] * dt_size_send);
758 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Alltoallv":"PMPI_Ialltoallv",
759 new simgrid::instr::VarCollTIData(request==MPI_REQUEST_IGNORED ? "alltoallv":"ialltoallv", -1, send_size, trace_sendcounts, recv_size,
760 trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
761 simgrid::smpi::Datatype::encode(recvtype)));
763 if(request == MPI_REQUEST_IGNORED)
764 retval = simgrid::smpi::Colls::alltoallv(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptype, recvbuf, recvcounts,
765 recvdisps, recvtype, comm);
767 retval = simgrid::smpi::Colls::ialltoallv(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptype, recvbuf, recvcounts,
768 recvdisps, recvtype, comm, request);
769 TRACE_smpi_comm_out(rank);
771 if (sendbuf == MPI_IN_PLACE) {
772 xbt_free(sendtmpbuf);
773 xbt_free(sendtmpcounts);
774 xbt_free(sendtmpdisps);
782 int PMPI_Alltoallw(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype* sendtypes, void* recvbuf,
783 int* recvcounts, int* recvdisps, MPI_Datatype* recvtypes, MPI_Comm comm)
785 return PMPI_Ialltoallw(sendbuf, sendcounts, senddisps, sendtypes, recvbuf, recvcounts, recvdisps, recvtypes, comm, MPI_REQUEST_IGNORED);
788 int PMPI_Ialltoallw(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype* sendtypes, void* recvbuf,
789 int* recvcounts, int* recvdisps, MPI_Datatype* recvtypes, MPI_Comm comm, MPI_Request *request)
795 if (comm == MPI_COMM_NULL) {
796 retval = MPI_ERR_COMM;
797 } else if ((sendbuf != MPI_IN_PLACE && sendtypes == nullptr) || recvtypes == nullptr) {
798 retval = MPI_ERR_TYPE;
799 } else if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
800 recvdisps == nullptr) {
801 retval = MPI_ERR_ARG;
802 } else if (request == nullptr){
803 retval = MPI_ERR_ARG;
805 int rank = simgrid::s4u::this_actor::get_pid();
806 int size = comm->size();
809 std::vector<int>* trace_sendcounts = new std::vector<int>;
810 std::vector<int>* trace_recvcounts = new std::vector<int>;
812 void* sendtmpbuf = static_cast<char*>(sendbuf);
813 int* sendtmpcounts = sendcounts;
814 int* sendtmpdisps = senddisps;
815 MPI_Datatype* sendtmptypes = sendtypes;
816 unsigned long maxsize = 0;
817 for (int i = 0; i < size; i++) { // copy data to avoid bad free
818 if(recvtypes[i]==MPI_DATATYPE_NULL){
819 delete trace_recvcounts;
820 delete trace_sendcounts;
823 recv_size += recvcounts[i] * recvtypes[i]->size();
824 trace_recvcounts->push_back(recvcounts[i] * recvtypes[i]->size());
825 if ((recvdisps[i] + (recvcounts[i] * recvtypes[i]->size())) > maxsize)
826 maxsize = recvdisps[i] + (recvcounts[i] * recvtypes[i]->size());
829 if (sendbuf == MPI_IN_PLACE) {
830 sendtmpbuf = static_cast<void*>(xbt_malloc(maxsize));
831 memcpy(sendtmpbuf, recvbuf, maxsize);
832 sendtmpcounts = static_cast<int*>(xbt_malloc(size * sizeof(int)));
833 memcpy(sendtmpcounts, recvcounts, size * sizeof(int));
834 sendtmpdisps = static_cast<int*>(xbt_malloc(size * sizeof(int)));
835 memcpy(sendtmpdisps, recvdisps, size * sizeof(int));
836 sendtmptypes = static_cast<MPI_Datatype*>(xbt_malloc(size * sizeof(MPI_Datatype)));
837 memcpy(sendtmptypes, recvtypes, size * sizeof(MPI_Datatype));
840 for (int i = 0; i < size; i++) { // copy data to avoid bad free
841 send_size += sendtmpcounts[i] * sendtmptypes[i]->size();
842 trace_sendcounts->push_back(sendtmpcounts[i] * sendtmptypes[i]->size());
845 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Alltoallw":"PMPI_Ialltoallw",
846 new simgrid::instr::VarCollTIData(request==MPI_REQUEST_IGNORED ? "alltoallv":"ialltoallv", -1, send_size, trace_sendcounts, recv_size,
847 trace_recvcounts, simgrid::smpi::Datatype::encode(sendtmptypes[0]),
848 simgrid::smpi::Datatype::encode(recvtypes[0])));
850 if(request == MPI_REQUEST_IGNORED)
851 retval = simgrid::smpi::Colls::alltoallw(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptypes, recvbuf, recvcounts,
852 recvdisps, recvtypes, comm);
854 retval = simgrid::smpi::Colls::ialltoallw(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptypes, recvbuf, recvcounts,
855 recvdisps, recvtypes, comm, request);
856 TRACE_smpi_comm_out(rank);
858 if (sendbuf == MPI_IN_PLACE) {
859 xbt_free(sendtmpbuf);
860 xbt_free(sendtmpcounts);
861 xbt_free(sendtmpdisps);
862 xbt_free(sendtmptypes);