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 (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()) {
59 retval = MPI_ERR_TYPE;
60 } else if (count < 0){
61 retval = MPI_ERR_COUNT;
62 } else if (root < 0 || root >= comm->size()){
63 retval = MPI_ERR_ROOT;
64 } else if (request == nullptr){
67 int rank = simgrid::s4u::this_actor::get_pid();
68 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Bcast":"PMPI_Ibcast",
69 new simgrid::instr::CollTIData(request==MPI_REQUEST_IGNORED?"bcast":"ibcast", root, -1.0,
70 datatype->is_replayable() ? count : count * datatype->size(), -1,
71 simgrid::smpi::Datatype::encode(datatype), ""));
72 if (comm->size() > 1){
73 if(request==MPI_REQUEST_IGNORED)
74 simgrid::smpi::Colls::bcast(buf, count, datatype, root, comm);
76 simgrid::smpi::Colls::ibcast(buf, count, datatype, root, comm, request);
78 if(request!=MPI_REQUEST_IGNORED)
79 *request = MPI_REQUEST_NULL;
83 TRACE_smpi_comm_out(rank);
89 int PMPI_Gather(void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
90 int root, MPI_Comm comm){
91 return PMPI_Igather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
94 int PMPI_Igather(void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
95 int root, MPI_Comm comm, MPI_Request *request)
101 if (comm == MPI_COMM_NULL) {
102 retval = MPI_ERR_COMM;
103 } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
104 ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL))){
105 retval = MPI_ERR_TYPE;
106 } else if ((( sendbuf != MPI_IN_PLACE) && (sendcount <0)) || ((comm->rank() == root) && (recvcount <0))){
107 retval = MPI_ERR_COUNT;
108 } else if (root < 0 || root >= comm->size()){
109 retval = MPI_ERR_ROOT;
110 } else if (request == nullptr){
111 retval = MPI_ERR_ARG;
114 char* sendtmpbuf = static_cast<char*>(sendbuf);
115 int sendtmpcount = sendcount;
116 MPI_Datatype sendtmptype = sendtype;
117 if( (comm->rank() == root) && (sendbuf == MPI_IN_PLACE )) {
119 sendtmptype=recvtype;
121 int rank = simgrid::s4u::this_actor::get_pid();
124 rank, request==MPI_REQUEST_IGNORED?"PMPI_Gather":"PMPI_Igather",
125 new simgrid::instr::CollTIData(
126 request==MPI_REQUEST_IGNORED ? "gather":"igather", root, -1.0, sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(),
127 (comm->rank() != root || recvtype->is_replayable()) ? recvcount : recvcount * recvtype->size(),
128 simgrid::smpi::Datatype::encode(sendtmptype), simgrid::smpi::Datatype::encode(recvtype)));
129 if(request == MPI_REQUEST_IGNORED)
130 simgrid::smpi::Colls::gather(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, root, comm);
132 simgrid::smpi::Colls::igather(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, root, comm, request);
134 retval = MPI_SUCCESS;
135 TRACE_smpi_comm_out(rank);
142 int PMPI_Gatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, int *recvcounts, int *displs,
143 MPI_Datatype recvtype, int root, MPI_Comm comm){
144 return PMPI_Igatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, root, comm, MPI_REQUEST_IGNORED);
147 int PMPI_Igatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, int *recvcounts, int *displs,
148 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request *request)
154 if (comm == MPI_COMM_NULL) {
155 retval = MPI_ERR_COMM;
156 } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
157 ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL))){
158 retval = MPI_ERR_TYPE;
159 } else if (( sendbuf != MPI_IN_PLACE) && (sendcount <0)){
160 retval = MPI_ERR_COUNT;
161 } else if ((comm->rank() == root) && (recvcounts == nullptr || displs == nullptr)) {
162 retval = MPI_ERR_ARG;
163 } else if (root < 0 || root >= comm->size()){
164 retval = MPI_ERR_ROOT;
165 } else if (request == nullptr){
166 retval = MPI_ERR_ARG;
168 char* sendtmpbuf = static_cast<char*>(sendbuf);
169 int sendtmpcount = sendcount;
170 MPI_Datatype sendtmptype = sendtype;
171 if( (comm->rank() == root) && (sendbuf == MPI_IN_PLACE )) {
173 sendtmptype=recvtype;
176 int rank = simgrid::s4u::this_actor::get_pid();
177 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
179 std::vector<int>* trace_recvcounts = new std::vector<int>;
180 if (comm->rank() == root) {
181 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
182 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
185 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Gatherv":"PMPI_Igatherv",
186 new simgrid::instr::VarCollTIData(
187 request==MPI_REQUEST_IGNORED ? "gatherv":"igatherv", root,
188 sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(), nullptr,
189 dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtmptype),
190 simgrid::smpi::Datatype::encode(recvtype)));
191 if(request == MPI_REQUEST_IGNORED)
192 retval = simgrid::smpi::Colls::gatherv(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcounts, displs, recvtype, root, comm);
194 retval = simgrid::smpi::Colls::igatherv(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcounts, displs, recvtype, root, comm, request);
195 TRACE_smpi_comm_out(rank);
202 int PMPI_Allgather(void *sendbuf, int sendcount, MPI_Datatype sendtype,
203 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm){
204 return PMPI_Iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
207 int PMPI_Iallgather(void *sendbuf, int sendcount, MPI_Datatype sendtype,
208 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
210 int retval = MPI_SUCCESS;
214 if (comm == MPI_COMM_NULL) {
215 retval = MPI_ERR_COMM;
216 } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
217 (recvtype == MPI_DATATYPE_NULL)){
218 retval = MPI_ERR_TYPE;
219 } else if ((( sendbuf != MPI_IN_PLACE) && (sendcount <0)) ||
221 retval = MPI_ERR_COUNT;
222 } else if (request == nullptr){
223 retval = MPI_ERR_ARG;
225 if(sendbuf == MPI_IN_PLACE) {
226 sendbuf=static_cast<char*>(recvbuf)+recvtype->get_extent()*recvcount*comm->rank();
230 int rank = simgrid::s4u::this_actor::get_pid();
232 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Allgather":"PMPI_Iallggather",
233 new simgrid::instr::CollTIData(
234 request==MPI_REQUEST_IGNORED ? "allgather" : "iallgather", -1, -1.0, sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
235 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
236 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
237 if(request == MPI_REQUEST_IGNORED)
238 simgrid::smpi::Colls::allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
240 simgrid::smpi::Colls::iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, request);
241 TRACE_smpi_comm_out(rank);
247 int PMPI_Allgatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype,
248 void *recvbuf, int *recvcounts, int *displs, MPI_Datatype recvtype, MPI_Comm comm){
249 return PMPI_Iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, MPI_REQUEST_IGNORED);
252 int PMPI_Iallgatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype,
253 void *recvbuf, int *recvcounts, int *displs, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
259 if (comm == MPI_COMM_NULL) {
260 retval = MPI_ERR_COMM;
261 } else if (((sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) || (recvtype == MPI_DATATYPE_NULL)) {
262 retval = MPI_ERR_TYPE;
263 } else if (( sendbuf != MPI_IN_PLACE) && (sendcount <0)){
264 retval = MPI_ERR_COUNT;
265 } else if (recvcounts == nullptr || displs == nullptr) {
266 retval = MPI_ERR_ARG;
267 } else if (request == nullptr){
268 retval = MPI_ERR_ARG;
271 if(sendbuf == MPI_IN_PLACE) {
272 sendbuf=static_cast<char*>(recvbuf)+recvtype->get_extent()*displs[comm->rank()];
273 sendcount=recvcounts[comm->rank()];
276 int rank = simgrid::s4u::this_actor::get_pid();
277 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
279 std::vector<int>* trace_recvcounts = new std::vector<int>;
280 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
281 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
283 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Allgatherv":"PMPI_Iallgatherv",
284 new simgrid::instr::VarCollTIData(
285 request==MPI_REQUEST_IGNORED ? "allgatherv" : "iallgatherv", -1, sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
286 nullptr, dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
287 simgrid::smpi::Datatype::encode(recvtype)));
288 if(request == MPI_REQUEST_IGNORED)
289 simgrid::smpi::Colls::allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
291 simgrid::smpi::Colls::iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, request);
292 retval = MPI_SUCCESS;
293 TRACE_smpi_comm_out(rank);
300 int PMPI_Scatter(void *sendbuf, int sendcount, MPI_Datatype sendtype,
301 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
302 return PMPI_Iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
305 int PMPI_Iscatter(void *sendbuf, int sendcount, MPI_Datatype sendtype,
306 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
312 if (comm == MPI_COMM_NULL) {
313 retval = MPI_ERR_COMM;
314 } else if (((comm->rank() == root) && (sendtype == MPI_DATATYPE_NULL || not sendtype->is_valid())) ||
315 ((recvbuf != MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL || not recvtype->is_valid()))) {
316 retval = MPI_ERR_TYPE;
317 } else if (((comm->rank() == root) && (sendcount < 0)) ||
318 ((recvbuf != MPI_IN_PLACE) && (recvcount < 0))) {
319 retval = MPI_ERR_COUNT;
320 } else if ((sendbuf == recvbuf) ||
321 ((comm->rank()==root) && sendcount>0 && (sendbuf == nullptr))){
322 retval = MPI_ERR_BUFFER;
323 } else if (root < 0 || root >= comm->size()){
324 retval = MPI_ERR_ROOT;
325 } else if (request == nullptr){
326 retval = MPI_ERR_ARG;
329 if (recvbuf == MPI_IN_PLACE) {
331 recvcount = sendcount;
333 int rank = simgrid::s4u::this_actor::get_pid();
336 rank, request==MPI_REQUEST_IGNORED?"PMPI_Scatter":"PMPI_Iscatter",
337 new simgrid::instr::CollTIData(
338 request==MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
339 (comm->rank() != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
340 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
341 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
342 if(request == MPI_REQUEST_IGNORED)
343 simgrid::smpi::Colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
345 simgrid::smpi::Colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
346 retval = MPI_SUCCESS;
347 TRACE_smpi_comm_out(rank);
354 int PMPI_Scatterv(void *sendbuf, int *sendcounts, int *displs,
355 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
356 return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
359 int PMPI_Iscatterv(void *sendbuf, int *sendcounts, int *displs,
360 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request *request)
366 if (comm == MPI_COMM_NULL) {
367 retval = MPI_ERR_COMM;
368 } else if (sendcounts == nullptr || displs == nullptr) {
369 retval = MPI_ERR_ARG;
370 } else if (((comm->rank() == root) && (sendtype == MPI_DATATYPE_NULL)) ||
371 ((recvbuf != MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL))) {
372 retval = MPI_ERR_TYPE;
373 } else if (request == nullptr){
374 retval = MPI_ERR_ARG;
375 } else if (recvbuf != MPI_IN_PLACE && recvcount < 0){
376 retval = MPI_ERR_COUNT;
377 } else if (root < 0 || root >= comm->size()){
378 retval = MPI_ERR_ROOT;
380 if (recvbuf == MPI_IN_PLACE) {
382 if(sendcounts[comm->rank()]<0)
383 return MPI_ERR_COUNT;
384 recvcount = sendcounts[comm->rank()];
386 int rank = simgrid::s4u::this_actor::get_pid();
387 int dt_size_send = sendtype->is_replayable() ? 1 : sendtype->size();
389 std::vector<int>* trace_sendcounts = new std::vector<int>;
390 if (comm->rank() == root) {
391 for (int i = 0; i < comm->size(); i++){ // copy data to avoid bad free
392 trace_sendcounts->push_back(sendcounts[i] * dt_size_send);
394 return MPI_ERR_COUNT;
398 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Scatterv":"PMPI_Iscatterv",
399 new simgrid::instr::VarCollTIData(
400 request==MPI_REQUEST_IGNORED ? "scatterv":"iscatterv", root, dt_size_send, trace_sendcounts,
401 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(), nullptr,
402 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
403 if(request == MPI_REQUEST_IGNORED)
404 retval = simgrid::smpi::Colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
406 retval = simgrid::smpi::Colls::iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
408 TRACE_smpi_comm_out(rank);
415 int PMPI_Reduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
417 return PMPI_Ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, MPI_REQUEST_IGNORED);
420 int PMPI_Ireduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm, MPI_Request* request)
426 if (comm == MPI_COMM_NULL) {
427 retval = MPI_ERR_COMM;
428 } else if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()){
429 retval = MPI_ERR_TYPE;
430 } else if (op == MPI_OP_NULL) {
432 } else if (request == nullptr){
433 retval = MPI_ERR_ARG;
434 } else if (root < 0 || root >= comm->size()){
435 retval = MPI_ERR_ROOT;
436 } else if (count < 0){
437 retval = MPI_ERR_COUNT;
439 int rank = simgrid::s4u::this_actor::get_pid();
441 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Reduce":"PMPI_Ireduce",
442 new simgrid::instr::CollTIData(request==MPI_REQUEST_IGNORED ? "reduce":"ireduce", root, 0,
443 datatype->is_replayable() ? count : count * datatype->size(), -1,
444 simgrid::smpi::Datatype::encode(datatype), ""));
445 if(request == MPI_REQUEST_IGNORED)
446 simgrid::smpi::Colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
448 simgrid::smpi::Colls::ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, request);
451 retval = MPI_SUCCESS;
452 TRACE_smpi_comm_out(rank);
459 int PMPI_Reduce_local(void *inbuf, void *inoutbuf, int count, MPI_Datatype datatype, MPI_Op op){
463 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()){
464 retval = MPI_ERR_TYPE;
465 } else if (op == MPI_OP_NULL) {
467 } else if (count < 0){
468 retval = MPI_ERR_COUNT;
470 op->apply(inbuf, inoutbuf, &count, datatype);
471 retval = MPI_SUCCESS;
477 int PMPI_Allreduce(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(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
488 if (comm == MPI_COMM_NULL) {
489 retval = MPI_ERR_COMM;
490 } else if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()) {
491 retval = MPI_ERR_TYPE;
492 } else if (count < 0){
493 retval = MPI_ERR_COUNT;
494 } else if (op == MPI_OP_NULL) {
496 } else if (request == nullptr){
497 retval = MPI_ERR_ARG;
499 char* sendtmpbuf = static_cast<char*>(sendbuf);
500 if( sendbuf == MPI_IN_PLACE ) {
501 sendtmpbuf = static_cast<char*>(xbt_malloc(count*datatype->get_extent()));
502 simgrid::smpi::Datatype::copy(recvbuf, count, datatype,sendtmpbuf, count, datatype);
504 int rank = simgrid::s4u::this_actor::get_pid();
506 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ?"PMPI_Allreduce":"PMPI_Iallreduce",
507 new simgrid::instr::CollTIData(request==MPI_REQUEST_IGNORED ? "allreduce":"iallreduce", -1, 0,
508 datatype->is_replayable() ? count : count * datatype->size(), -1,
509 simgrid::smpi::Datatype::encode(datatype), ""));
511 if(request == MPI_REQUEST_IGNORED)
512 simgrid::smpi::Colls::allreduce(sendtmpbuf, recvbuf, count, datatype, op, comm);
514 simgrid::smpi::Colls::iallreduce(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
516 if( sendbuf == MPI_IN_PLACE )
517 xbt_free(sendtmpbuf);
519 retval = MPI_SUCCESS;
520 TRACE_smpi_comm_out(rank);
527 int PMPI_Scan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
529 return PMPI_Iscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
532 int PMPI_Iscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request)
538 if (comm == MPI_COMM_NULL) {
539 retval = MPI_ERR_COMM;
540 } else if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()){
541 retval = MPI_ERR_TYPE;
542 } else if (op == MPI_OP_NULL) {
544 } else if (request == nullptr){
545 retval = MPI_ERR_ARG;
546 } else if (count < 0){
547 retval = MPI_ERR_COUNT;
549 int rank = simgrid::s4u::this_actor::get_pid();
550 void* sendtmpbuf = sendbuf;
551 if (sendbuf == MPI_IN_PLACE) {
552 sendtmpbuf = static_cast<void*>(xbt_malloc(count * datatype->size()));
553 memcpy(sendtmpbuf, recvbuf, count * datatype->size());
555 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Scan" : "PMPI_Iscan", new simgrid::instr::Pt2PtTIData(
556 request==MPI_REQUEST_IGNORED ? "scan":"iscan", -1,
557 datatype->is_replayable() ? count : count * datatype->size(),
558 simgrid::smpi::Datatype::encode(datatype)));
560 if(request == MPI_REQUEST_IGNORED)
561 retval = simgrid::smpi::Colls::scan(sendtmpbuf, recvbuf, count, datatype, op, comm);
563 retval = simgrid::smpi::Colls::iscan(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
565 TRACE_smpi_comm_out(rank);
566 if (sendbuf == MPI_IN_PLACE)
567 xbt_free(sendtmpbuf);
574 int PMPI_Exscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
576 return PMPI_Iexscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
579 int PMPI_Iexscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request){
584 if (comm == MPI_COMM_NULL) {
585 retval = MPI_ERR_COMM;
586 } else if (not datatype->is_valid()) {
587 retval = MPI_ERR_TYPE;
588 } else if (op == MPI_OP_NULL) {
590 } else if (request == nullptr){
591 retval = MPI_ERR_ARG;
592 } else if (count < 0){
593 retval = MPI_ERR_COUNT;
595 int rank = simgrid::s4u::this_actor::get_pid();
596 void* sendtmpbuf = sendbuf;
597 if (sendbuf == MPI_IN_PLACE) {
598 sendtmpbuf = static_cast<void*>(xbt_malloc(count * datatype->size()));
599 memcpy(sendtmpbuf, recvbuf, count * datatype->size());
602 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Exscan" : "PMPI_Iexscan", new simgrid::instr::Pt2PtTIData(
603 request==MPI_REQUEST_IGNORED ? "exscan":"iexscan", -1, datatype->is_replayable() ? count : count * datatype->size(),
604 simgrid::smpi::Datatype::encode(datatype)));
606 if(request == MPI_REQUEST_IGNORED)
607 retval = simgrid::smpi::Colls::exscan(sendtmpbuf, recvbuf, count, datatype, op, comm);
609 retval = simgrid::smpi::Colls::iexscan(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
611 TRACE_smpi_comm_out(rank);
612 if (sendbuf == MPI_IN_PLACE)
613 xbt_free(sendtmpbuf);
620 int PMPI_Reduce_scatter(void *sendbuf, void *recvbuf, int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
622 return PMPI_Ireduce_scatter(sendbuf, recvbuf, recvcounts, datatype, op, comm, MPI_REQUEST_IGNORED);
625 int PMPI_Ireduce_scatter(void *sendbuf, void *recvbuf, int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
630 if (comm == MPI_COMM_NULL) {
631 retval = MPI_ERR_COMM;
632 } else if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()){
633 retval = MPI_ERR_TYPE;
634 } else if (op == MPI_OP_NULL) {
636 } else if (recvcounts == nullptr) {
637 retval = MPI_ERR_ARG;
638 } else if (request == nullptr){
639 retval = MPI_ERR_ARG;
641 int rank = simgrid::s4u::this_actor::get_pid();
642 std::vector<int>* trace_recvcounts = new std::vector<int>;
643 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
646 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
648 return MPI_ERR_COUNT;
649 trace_recvcounts->push_back(recvcounts[i] * dt_send_size);
650 totalcount += recvcounts[i];
653 void* sendtmpbuf = sendbuf;
654 if (sendbuf == MPI_IN_PLACE) {
655 sendtmpbuf = static_cast<void*>(xbt_malloc(totalcount * datatype->size()));
656 memcpy(sendtmpbuf, recvbuf, totalcount * datatype->size());
659 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter": "PMPI_Ireduce_scatter", new simgrid::instr::VarCollTIData(
660 request==MPI_REQUEST_IGNORED ? "reducescatter":"ireducescatter", -1, dt_send_size, nullptr, -1, trace_recvcounts,
661 simgrid::smpi::Datatype::encode(datatype), ""));
663 if(request == MPI_REQUEST_IGNORED)
664 simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
666 simgrid::smpi::Colls::ireduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm, request);
668 retval = MPI_SUCCESS;
669 TRACE_smpi_comm_out(rank);
671 if (sendbuf == MPI_IN_PLACE)
672 xbt_free(sendtmpbuf);
679 int PMPI_Reduce_scatter_block(void *sendbuf, void *recvbuf, int recvcount,
680 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
682 return PMPI_Ireduce_scatter_block(sendbuf, recvbuf, recvcount, datatype, op, comm, MPI_REQUEST_IGNORED);
685 int PMPI_Ireduce_scatter_block(void *sendbuf, void *recvbuf, int recvcount,
686 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
691 if (comm == MPI_COMM_NULL) {
692 retval = MPI_ERR_COMM;
693 } else if (not datatype->is_valid()) {
694 retval = MPI_ERR_TYPE;
695 } else if (op == MPI_OP_NULL) {
697 } else if (recvcount < 0) {
698 retval = MPI_ERR_ARG;
699 } else if (request == nullptr){
700 retval = MPI_ERR_ARG;
702 int count = comm->size();
704 int rank = simgrid::s4u::this_actor::get_pid();
705 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
706 std::vector<int>* trace_recvcounts = new std::vector<int>(recvcount * dt_send_size); // copy data to avoid bad free
708 void* sendtmpbuf = sendbuf;
709 if (sendbuf == MPI_IN_PLACE) {
710 sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * count * datatype->size()));
711 memcpy(sendtmpbuf, recvbuf, recvcount * count * datatype->size());
714 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter_block":"PMPI_Ireduce_scatter_block",
715 new simgrid::instr::VarCollTIData(request==MPI_REQUEST_IGNORED ? "reducescatter":"ireducescatter", -1, 0, nullptr, -1, trace_recvcounts,
716 simgrid::smpi::Datatype::encode(datatype), ""));
718 int* recvcounts = new int[count];
719 for (int i = 0; i < count; i++)
720 recvcounts[i] = recvcount;
721 if(request == MPI_REQUEST_IGNORED)
722 simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
724 simgrid::smpi::Colls::ireduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm, request);
726 retval = MPI_SUCCESS;
728 TRACE_smpi_comm_out(rank);
730 if (sendbuf == MPI_IN_PLACE)
731 xbt_free(sendtmpbuf);
737 int PMPI_Alltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
738 MPI_Datatype recvtype, MPI_Comm comm){
739 return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
742 int PMPI_Ialltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
743 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request *request)
748 if (comm == MPI_COMM_NULL) {
749 retval = MPI_ERR_COMM;
750 } else if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL) {
751 retval = MPI_ERR_TYPE;
752 } else if ((sendbuf != MPI_IN_PLACE && sendcount < 0) || recvcount < 0){
753 retval = MPI_ERR_COUNT;
754 } else if (request == nullptr){
755 retval = MPI_ERR_ARG;
757 int rank = simgrid::s4u::this_actor::get_pid();
758 void* sendtmpbuf = static_cast<char*>(sendbuf);
759 int sendtmpcount = sendcount;
760 MPI_Datatype sendtmptype = sendtype;
761 if (sendbuf == MPI_IN_PLACE) {
762 sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * comm->size() * recvtype->size()));
763 memcpy(sendtmpbuf, recvbuf, recvcount * comm->size() * recvtype->size());
764 sendtmpcount = recvcount;
765 sendtmptype = recvtype;
768 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Alltoall":"PMPI_Ialltoall",
769 new simgrid::instr::CollTIData(
770 request==MPI_REQUEST_IGNORED ? "alltoall" : "ialltoall", -1, -1.0,
771 sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(),
772 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
773 simgrid::smpi::Datatype::encode(sendtmptype), simgrid::smpi::Datatype::encode(recvtype)));
774 if(request == MPI_REQUEST_IGNORED)
775 retval = simgrid::smpi::Colls::alltoall(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, comm);
777 retval = simgrid::smpi::Colls::ialltoall(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, comm, request);
779 TRACE_smpi_comm_out(rank);
781 if (sendbuf == MPI_IN_PLACE)
782 xbt_free(sendtmpbuf);
789 int PMPI_Alltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf,
790 int* recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm)
792 return PMPI_Ialltoallv(sendbuf, sendcounts, senddisps, sendtype, recvbuf, recvcounts, recvdisps, recvtype, comm, MPI_REQUEST_IGNORED);
795 int PMPI_Ialltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf,
796 int* recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request *request)
802 if (comm == MPI_COMM_NULL) {
803 retval = MPI_ERR_COMM;
804 } else if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL) {
805 retval = MPI_ERR_TYPE;
806 } else if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
807 recvdisps == nullptr) {
808 retval = MPI_ERR_ARG;
809 } else if (request == nullptr){
810 retval = MPI_ERR_ARG;
812 int rank = simgrid::s4u::this_actor::get_pid();
813 int size = comm->size();
816 std::vector<int>* trace_sendcounts = new std::vector<int>;
817 std::vector<int>* trace_recvcounts = new std::vector<int>;
818 int dt_size_recv = recvtype->size();
820 void* sendtmpbuf = static_cast<char*>(sendbuf);
821 int* sendtmpcounts = sendcounts;
822 int* sendtmpdisps = senddisps;
823 MPI_Datatype sendtmptype = sendtype;
825 for (int i = 0; i < size; i++) { // copy data to avoid bad free
826 if (recvcounts[i] <0 || (sendbuf != MPI_IN_PLACE && sendcounts[i]<0))
827 return MPI_ERR_COUNT;
828 recv_size += recvcounts[i] * dt_size_recv;
829 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
830 if (((recvdisps[i] + recvcounts[i]) * dt_size_recv) > maxsize)
831 maxsize = (recvdisps[i] + recvcounts[i]) * dt_size_recv;
834 if (sendbuf == MPI_IN_PLACE) {
835 sendtmpbuf = static_cast<void*>(xbt_malloc(maxsize));
836 memcpy(sendtmpbuf, recvbuf, maxsize);
837 sendtmpcounts = static_cast<int*>(xbt_malloc(size * sizeof(int)));
838 memcpy(sendtmpcounts, recvcounts, size * sizeof(int));
839 sendtmpdisps = static_cast<int*>(xbt_malloc(size * sizeof(int)));
840 memcpy(sendtmpdisps, recvdisps, size * sizeof(int));
841 sendtmptype = recvtype;
844 int dt_size_send = sendtmptype->size();
846 for (int i = 0; i < size; i++) { // copy data to avoid bad free
847 send_size += sendtmpcounts[i] * dt_size_send;
848 trace_sendcounts->push_back(sendtmpcounts[i] * dt_size_send);
851 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Alltoallv":"PMPI_Ialltoallv",
852 new simgrid::instr::VarCollTIData(request==MPI_REQUEST_IGNORED ? "alltoallv":"ialltoallv", -1, send_size, trace_sendcounts, recv_size,
853 trace_recvcounts, simgrid::smpi::Datatype::encode(sendtmptype),
854 simgrid::smpi::Datatype::encode(recvtype)));
856 if(request == MPI_REQUEST_IGNORED)
857 retval = simgrid::smpi::Colls::alltoallv(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptype, recvbuf, recvcounts,
858 recvdisps, recvtype, comm);
860 retval = simgrid::smpi::Colls::ialltoallv(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptype, recvbuf, recvcounts,
861 recvdisps, recvtype, comm, request);
862 TRACE_smpi_comm_out(rank);
864 if (sendbuf == MPI_IN_PLACE) {
865 xbt_free(sendtmpbuf);
866 xbt_free(sendtmpcounts);
867 xbt_free(sendtmpdisps);
875 int PMPI_Alltoallw(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype* sendtypes, void* recvbuf,
876 int* recvcounts, int* recvdisps, MPI_Datatype* recvtypes, MPI_Comm comm)
878 return PMPI_Ialltoallw(sendbuf, sendcounts, senddisps, sendtypes, recvbuf, recvcounts, recvdisps, recvtypes, comm, MPI_REQUEST_IGNORED);
881 int PMPI_Ialltoallw(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype* sendtypes, void* recvbuf,
882 int* recvcounts, int* recvdisps, MPI_Datatype* recvtypes, MPI_Comm comm, MPI_Request *request)
888 if (comm == MPI_COMM_NULL) {
889 retval = MPI_ERR_COMM;
890 } else if ((sendbuf != MPI_IN_PLACE && sendtypes == nullptr) || recvtypes == nullptr) {
891 retval = MPI_ERR_TYPE;
892 } else if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
893 recvdisps == nullptr) {
894 retval = MPI_ERR_ARG;
895 } else if (request == nullptr){
896 retval = MPI_ERR_ARG;
898 int rank = simgrid::s4u::this_actor::get_pid();
899 int size = comm->size();
902 std::vector<int>* trace_sendcounts = new std::vector<int>;
903 std::vector<int>* trace_recvcounts = new std::vector<int>;
905 void* sendtmpbuf = static_cast<char*>(sendbuf);
906 int* sendtmpcounts = sendcounts;
907 int* sendtmpdisps = senddisps;
908 MPI_Datatype* sendtmptypes = sendtypes;
909 unsigned long maxsize = 0;
910 for (int i = 0; i < size; i++) { // copy data to avoid bad free
911 if(recvtypes[i]==MPI_DATATYPE_NULL){
912 delete trace_recvcounts;
913 delete trace_sendcounts;
916 recv_size += recvcounts[i] * recvtypes[i]->size();
917 trace_recvcounts->push_back(recvcounts[i] * recvtypes[i]->size());
918 if ((recvdisps[i] + (recvcounts[i] * recvtypes[i]->size())) > maxsize)
919 maxsize = recvdisps[i] + (recvcounts[i] * recvtypes[i]->size());
922 if (sendbuf == MPI_IN_PLACE) {
923 sendtmpbuf = static_cast<void*>(xbt_malloc(maxsize));
924 memcpy(sendtmpbuf, recvbuf, maxsize);
925 sendtmpcounts = static_cast<int*>(xbt_malloc(size * sizeof(int)));
926 memcpy(sendtmpcounts, recvcounts, size * sizeof(int));
927 sendtmpdisps = static_cast<int*>(xbt_malloc(size * sizeof(int)));
928 memcpy(sendtmpdisps, recvdisps, size * sizeof(int));
929 sendtmptypes = static_cast<MPI_Datatype*>(xbt_malloc(size * sizeof(MPI_Datatype)));
930 memcpy(sendtmptypes, recvtypes, size * sizeof(MPI_Datatype));
933 for (int i = 0; i < size; i++) { // copy data to avoid bad free
934 send_size += sendtmpcounts[i] * sendtmptypes[i]->size();
935 trace_sendcounts->push_back(sendtmpcounts[i] * sendtmptypes[i]->size());
938 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Alltoallw":"PMPI_Ialltoallw",
939 new simgrid::instr::VarCollTIData(request==MPI_REQUEST_IGNORED ? "alltoallv":"ialltoallv", -1, send_size, trace_sendcounts, recv_size,
940 trace_recvcounts, simgrid::smpi::Datatype::encode(sendtmptypes[0]),
941 simgrid::smpi::Datatype::encode(recvtypes[0])));
943 if(request == MPI_REQUEST_IGNORED)
944 retval = simgrid::smpi::Colls::alltoallw(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptypes, recvbuf, recvcounts,
945 recvdisps, recvtypes, comm);
947 retval = simgrid::smpi::Colls::ialltoallw(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptypes, recvbuf, recvcounts,
948 recvdisps, recvtypes, comm, request);
949 TRACE_smpi_comm_out(rank);
951 if (sendbuf == MPI_IN_PLACE) {
952 xbt_free(sendtmpbuf);
953 xbt_free(sendtmpcounts);
954 xbt_free(sendtmpdisps);
955 xbt_free(sendtmptypes);