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_Barrier(MPI_Comm comm)
20 return PMPI_Ibarrier(comm, MPI_REQUEST_IGNORED);
23 int PMPI_Ibarrier(MPI_Comm comm, MPI_Request *request)
27 if (comm == MPI_COMM_NULL) {
28 retval = MPI_ERR_COMM;
29 } else if(request == nullptr){
32 int rank = simgrid::s4u::this_actor::get_pid();
33 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED? "PMPI_Barrier" : "PMPI_Ibarrier", new simgrid::instr::NoOpTIData(request==MPI_REQUEST_IGNORED? "barrier" : "ibarrier"));
34 if(request==MPI_REQUEST_IGNORED){
35 simgrid::smpi::Colls::barrier(comm);
36 //Barrier can be used to synchronize RMA calls. Finish all requests from comm before.
37 comm->finish_rma_calls();
39 simgrid::smpi::Colls::ibarrier(comm, request);
40 TRACE_smpi_comm_out(rank);
46 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
48 return PMPI_Ibcast(buf, count, datatype, root, comm, MPI_REQUEST_IGNORED);
51 int PMPI_Ibcast(void *buf, int count, MPI_Datatype datatype,
52 int root, MPI_Comm comm, MPI_Request* request)
54 if (comm == MPI_COMM_NULL) {
56 } if (buf == nullptr && count > 0) {
57 return MPI_ERR_BUFFER;
58 } else if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()) {
60 } else if (count < 0){
62 } else if (root < 0 || root >= comm->size()){
64 } else if (request == nullptr){
68 int rank = simgrid::s4u::this_actor::get_pid();
69 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Bcast":"PMPI_Ibcast",
70 new simgrid::instr::CollTIData(request==MPI_REQUEST_IGNORED?"bcast":"ibcast", root, -1.0,
71 datatype->is_replayable() ? count : count * datatype->size(), -1,
72 simgrid::smpi::Datatype::encode(datatype), ""));
73 if (comm->size() > 1){
74 if(request==MPI_REQUEST_IGNORED)
75 simgrid::smpi::Colls::bcast(buf, count, datatype, root, comm);
77 simgrid::smpi::Colls::ibcast(buf, count, datatype, root, comm, request);
79 if(request!=MPI_REQUEST_IGNORED)
80 *request = MPI_REQUEST_NULL;
82 TRACE_smpi_comm_out(rank);
88 int PMPI_Gather(void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
89 int root, MPI_Comm comm){
90 return PMPI_Igather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
93 int PMPI_Igather(void *sendbuf, int sendcount, MPI_Datatype sendtype,void *recvbuf, int recvcount, MPI_Datatype recvtype,
94 int root, MPI_Comm comm, MPI_Request *request)
96 if (comm == MPI_COMM_NULL) {
98 } else if ((sendbuf == nullptr) || ((comm->rank() == root) && recvbuf == nullptr)) {
99 return MPI_ERR_BUFFER;
100 } else if (((sendbuf != MPI_IN_PLACE && sendcount > 0) && (sendtype == MPI_DATATYPE_NULL)) ||
101 ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL))){
103 } else if ((( sendbuf != MPI_IN_PLACE) && (sendcount <0)) || ((comm->rank() == root) && (recvcount <0))){
104 return MPI_ERR_COUNT;
105 } else if (root < 0 || root >= comm->size()){
107 } else if (request == nullptr){
111 char* sendtmpbuf = static_cast<char*>(sendbuf);
112 int sendtmpcount = sendcount;
113 MPI_Datatype sendtmptype = sendtype;
114 if( (comm->rank() == root) && (sendbuf == MPI_IN_PLACE )) {
116 sendtmptype=recvtype;
118 int rank = simgrid::s4u::this_actor::get_pid();
121 rank, request==MPI_REQUEST_IGNORED?"PMPI_Gather":"PMPI_Igather",
122 new simgrid::instr::CollTIData(
123 request==MPI_REQUEST_IGNORED ? "gather":"igather", root, -1.0, sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(),
124 (comm->rank() != root || recvtype->is_replayable()) ? recvcount : recvcount * recvtype->size(),
125 simgrid::smpi::Datatype::encode(sendtmptype), simgrid::smpi::Datatype::encode(recvtype)));
126 if(request == MPI_REQUEST_IGNORED)
127 simgrid::smpi::Colls::gather(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, root, comm);
129 simgrid::smpi::Colls::igather(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, root, comm, request);
131 TRACE_smpi_comm_out(rank);
137 int PMPI_Gatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, int *recvcounts, int *displs,
138 MPI_Datatype recvtype, int root, MPI_Comm comm){
139 return PMPI_Igatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, root, comm, MPI_REQUEST_IGNORED);
142 int PMPI_Igatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, int *recvcounts, int *displs,
143 MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request *request)
145 if (comm == MPI_COMM_NULL) {
147 } else if ((sendbuf == nullptr && sendcount > 0) || ((comm->rank() == root) && recvbuf == nullptr)) {
148 return MPI_ERR_BUFFER;
149 } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
150 ((comm->rank() == root) && (recvtype == MPI_DATATYPE_NULL))){
152 } else if (( sendbuf != MPI_IN_PLACE) && (sendcount <0)){
153 return MPI_ERR_COUNT;
154 } else if ((comm->rank() == root) && (recvcounts == nullptr || displs == nullptr)) {
156 } else if (root < 0 || root >= comm->size()){
158 } else if (request == nullptr){
161 for (int i = 0; i < comm->size(); i++){
162 if((comm->rank() == root) && (recvcounts[i]<0))
163 return MPI_ERR_COUNT;
167 char* sendtmpbuf = static_cast<char*>(sendbuf);
168 int sendtmpcount = sendcount;
169 MPI_Datatype sendtmptype = sendtype;
170 if( (comm->rank() == root) && (sendbuf == MPI_IN_PLACE )) {
172 sendtmptype=recvtype;
175 int rank = simgrid::s4u::this_actor::get_pid();
176 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
178 std::vector<int>* trace_recvcounts = new std::vector<int>;
179 if (comm->rank() == root) {
180 for (int i = 0; i < comm->size(); i++) // copy data to avoid bad free
181 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
184 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Gatherv":"PMPI_Igatherv",
185 new simgrid::instr::VarCollTIData(
186 request==MPI_REQUEST_IGNORED ? "gatherv":"igatherv", root,
187 sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(), nullptr,
188 dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtmptype),
189 simgrid::smpi::Datatype::encode(recvtype)));
190 if(request == MPI_REQUEST_IGNORED)
191 simgrid::smpi::Colls::gatherv(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcounts, displs, recvtype, root, comm);
193 simgrid::smpi::Colls::igatherv(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcounts, displs, recvtype, root, comm, request);
194 TRACE_smpi_comm_out(rank);
200 int PMPI_Allgather(void *sendbuf, int sendcount, MPI_Datatype sendtype,
201 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm){
202 return PMPI_Iallgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
205 int PMPI_Iallgather(void *sendbuf, int sendcount, MPI_Datatype sendtype,
206 void *recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request* request)
208 int retval = MPI_SUCCESS;
212 if (comm == MPI_COMM_NULL) {
213 retval = MPI_ERR_COMM;
214 } else if ((sendbuf == nullptr && sendcount > 0) || (recvbuf == nullptr)){
215 retval = MPI_ERR_BUFFER;
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)
255 if (comm == MPI_COMM_NULL) {
257 } else if ((sendbuf == nullptr && sendcount > 0) || (recvbuf == nullptr)){
258 return MPI_ERR_BUFFER;
259 } else if (((sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) || (recvtype == MPI_DATATYPE_NULL)) {
261 } else if (( sendbuf != MPI_IN_PLACE) && (sendcount <0)){
262 return MPI_ERR_COUNT;
263 } else if (recvcounts == nullptr || displs == nullptr) {
265 } else if (request == nullptr){
268 for (int i = 0; i < comm->size(); i++){ // copy data to avoid bad free
269 if (recvcounts[i] < 0)
270 return MPI_ERR_COUNT;
274 if(sendbuf == MPI_IN_PLACE) {
275 sendbuf=static_cast<char*>(recvbuf)+recvtype->get_extent()*displs[comm->rank()];
276 sendcount=recvcounts[comm->rank()];
279 int rank = simgrid::s4u::this_actor::get_pid();
280 int dt_size_recv = recvtype->is_replayable() ? 1 : recvtype->size();
282 std::vector<int>* trace_recvcounts = new std::vector<int>;
283 for (int i = 0; i < comm->size(); i++){ // copy data to avoid bad free
284 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
287 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Allgatherv":"PMPI_Iallgatherv",
288 new simgrid::instr::VarCollTIData(
289 request==MPI_REQUEST_IGNORED ? "allgatherv" : "iallgatherv", -1, sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
290 nullptr, dt_size_recv, trace_recvcounts, simgrid::smpi::Datatype::encode(sendtype),
291 simgrid::smpi::Datatype::encode(recvtype)));
292 if(request == MPI_REQUEST_IGNORED)
293 simgrid::smpi::Colls::allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
295 simgrid::smpi::Colls::iallgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm, request);
296 TRACE_smpi_comm_out(rank);
302 int PMPI_Scatter(void *sendbuf, int sendcount, MPI_Datatype sendtype,
303 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
304 return PMPI_Iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
307 int PMPI_Iscatter(void *sendbuf, int sendcount, MPI_Datatype sendtype,
308 void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request* request)
310 if (comm == MPI_COMM_NULL) {
312 } else if (((comm->rank() == root) && (sendtype == MPI_DATATYPE_NULL || not sendtype->is_valid())) ||
313 ((recvbuf != MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL || not recvtype->is_valid()))) {
315 } else if (((comm->rank() == root) && (sendcount < 0)) ||
316 ((recvbuf != MPI_IN_PLACE) && (recvcount < 0))) {
317 return MPI_ERR_COUNT;
318 } else if ((sendbuf == recvbuf) ||
319 ((comm->rank()==root) && sendcount>0 && (sendbuf == nullptr)) ||
320 (recvcount > 0 && recvbuf == nullptr)){
321 return MPI_ERR_BUFFER;
322 } else if (root < 0 || root >= comm->size()){
324 } else if (request == nullptr){
328 if (recvbuf == MPI_IN_PLACE) {
330 recvcount = sendcount;
332 int rank = simgrid::s4u::this_actor::get_pid();
335 rank, request==MPI_REQUEST_IGNORED?"PMPI_Scatter":"PMPI_Iscatter",
336 new simgrid::instr::CollTIData(
337 request==MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
338 (comm->rank() != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
339 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
340 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
341 if(request == MPI_REQUEST_IGNORED)
342 simgrid::smpi::Colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
344 simgrid::smpi::Colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
345 TRACE_smpi_comm_out(rank);
351 int PMPI_Scatterv(void *sendbuf, int *sendcounts, int *displs,
352 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
353 return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
356 int PMPI_Iscatterv(void *sendbuf, int *sendcounts, int *displs,
357 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request *request)
359 if (comm == MPI_COMM_NULL) {
361 } else if (sendcounts == nullptr || displs == nullptr) {
363 } else if (((comm->rank() == root) && (sendtype == MPI_DATATYPE_NULL)) ||
364 ((recvbuf != MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL))) {
366 } else if (request == nullptr){
368 } else if (recvbuf != MPI_IN_PLACE && recvcount < 0){
369 return MPI_ERR_COUNT;
370 } else if (root < 0 || root >= comm->size()){
373 if (comm->rank() == root){
374 if(recvbuf == MPI_IN_PLACE) {
376 recvcount = sendcounts[comm->rank()];
378 for (int i = 0; i < comm->size(); i++){
380 return MPI_ERR_COUNT;
385 int rank = simgrid::s4u::this_actor::get_pid();
386 int dt_size_send = sendtype->is_replayable() ? 1 : sendtype->size();
388 std::vector<int>* trace_sendcounts = new std::vector<int>;
389 if (comm->rank() == root) {
390 for (int i = 0; i < comm->size(); i++){ // copy data to avoid bad free
391 trace_sendcounts->push_back(sendcounts[i] * dt_size_send);
395 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Scatterv":"PMPI_Iscatterv",
396 new simgrid::instr::VarCollTIData(
397 request==MPI_REQUEST_IGNORED ? "scatterv":"iscatterv", root, dt_size_send, trace_sendcounts,
398 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(), nullptr,
399 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
400 if(request == MPI_REQUEST_IGNORED)
401 simgrid::smpi::Colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
403 simgrid::smpi::Colls::iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
404 TRACE_smpi_comm_out(rank);
411 int PMPI_Reduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
413 return PMPI_Ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, MPI_REQUEST_IGNORED);
416 int PMPI_Ireduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm, MPI_Request* request)
418 if (comm == MPI_COMM_NULL) {
420 } else if ((sendbuf == nullptr && count > 0) || ((comm->rank() == root) && recvbuf == nullptr)) {
421 return MPI_ERR_BUFFER;
422 } else if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()){
424 } else if (op == MPI_OP_NULL) {
426 } else if (request == nullptr){
428 } else if (root < 0 || root >= comm->size()){
430 } else if (count < 0){
431 return MPI_ERR_COUNT;
434 int rank = simgrid::s4u::this_actor::get_pid();
436 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Reduce":"PMPI_Ireduce",
437 new simgrid::instr::CollTIData(request==MPI_REQUEST_IGNORED ? "reduce":"ireduce", root, 0,
438 datatype->is_replayable() ? count : count * datatype->size(), -1,
439 simgrid::smpi::Datatype::encode(datatype), ""));
440 if(request == MPI_REQUEST_IGNORED)
441 simgrid::smpi::Colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
443 simgrid::smpi::Colls::ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, request);
445 TRACE_smpi_comm_out(rank);
451 int PMPI_Reduce_local(void *inbuf, void *inoutbuf, int count, MPI_Datatype datatype, MPI_Op op){
455 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()){
456 retval = MPI_ERR_TYPE;
457 } else if (op == MPI_OP_NULL) {
459 } else if (count < 0){
460 retval = MPI_ERR_COUNT;
462 op->apply(inbuf, inoutbuf, &count, datatype);
463 retval = MPI_SUCCESS;
469 int PMPI_Allreduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
471 return PMPI_Iallreduce(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
474 int PMPI_Iallreduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
476 if (comm == MPI_COMM_NULL) {
478 } else if ((sendbuf == nullptr && count > 0) || (recvbuf == nullptr)) {
479 return MPI_ERR_BUFFER;
480 } else if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()) {
482 } else if (count < 0){
483 return MPI_ERR_COUNT;
484 } else if (op == MPI_OP_NULL) {
486 } else if (request == nullptr){
490 char* sendtmpbuf = static_cast<char*>(sendbuf);
491 if( sendbuf == MPI_IN_PLACE ) {
492 sendtmpbuf = static_cast<char*>(xbt_malloc(count*datatype->get_extent()));
493 simgrid::smpi::Datatype::copy(recvbuf, count, datatype,sendtmpbuf, count, datatype);
495 int rank = simgrid::s4u::this_actor::get_pid();
497 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ?"PMPI_Allreduce":"PMPI_Iallreduce",
498 new simgrid::instr::CollTIData(request==MPI_REQUEST_IGNORED ? "allreduce":"iallreduce", -1, 0,
499 datatype->is_replayable() ? count : count * datatype->size(), -1,
500 simgrid::smpi::Datatype::encode(datatype), ""));
502 if(request == MPI_REQUEST_IGNORED)
503 simgrid::smpi::Colls::allreduce(sendtmpbuf, recvbuf, count, datatype, op, comm);
505 simgrid::smpi::Colls::iallreduce(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
507 if( sendbuf == MPI_IN_PLACE )
508 xbt_free(sendtmpbuf);
510 TRACE_smpi_comm_out(rank);
516 int PMPI_Scan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
518 return PMPI_Iscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
521 int PMPI_Iscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request)
527 if (comm == MPI_COMM_NULL) {
528 retval = MPI_ERR_COMM;
529 } else if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()){
530 retval = MPI_ERR_TYPE;
531 } else if (op == MPI_OP_NULL) {
533 } else if (request == nullptr){
534 retval = MPI_ERR_ARG;
535 } else if (count < 0){
536 retval = MPI_ERR_COUNT;
538 int rank = simgrid::s4u::this_actor::get_pid();
539 void* sendtmpbuf = sendbuf;
540 if (sendbuf == MPI_IN_PLACE) {
541 sendtmpbuf = static_cast<void*>(xbt_malloc(count * datatype->size()));
542 memcpy(sendtmpbuf, recvbuf, count * datatype->size());
544 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Scan" : "PMPI_Iscan", new simgrid::instr::Pt2PtTIData(
545 request==MPI_REQUEST_IGNORED ? "scan":"iscan", -1,
546 datatype->is_replayable() ? count : count * datatype->size(),
547 simgrid::smpi::Datatype::encode(datatype)));
549 if(request == MPI_REQUEST_IGNORED)
550 retval = simgrid::smpi::Colls::scan(sendtmpbuf, recvbuf, count, datatype, op, comm);
552 retval = simgrid::smpi::Colls::iscan(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
554 TRACE_smpi_comm_out(rank);
555 if (sendbuf == MPI_IN_PLACE)
556 xbt_free(sendtmpbuf);
563 int PMPI_Exscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
565 return PMPI_Iexscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
568 int PMPI_Iexscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request){
573 if (comm == MPI_COMM_NULL) {
574 retval = MPI_ERR_COMM;
575 } else if (not datatype->is_valid()) {
576 retval = MPI_ERR_TYPE;
577 } else if (op == MPI_OP_NULL) {
579 } else if (request == nullptr){
580 retval = MPI_ERR_ARG;
581 } else if (count < 0){
582 retval = MPI_ERR_COUNT;
584 int rank = simgrid::s4u::this_actor::get_pid();
585 void* sendtmpbuf = sendbuf;
586 if (sendbuf == MPI_IN_PLACE) {
587 sendtmpbuf = static_cast<void*>(xbt_malloc(count * datatype->size()));
588 memcpy(sendtmpbuf, recvbuf, count * datatype->size());
591 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Exscan" : "PMPI_Iexscan", new simgrid::instr::Pt2PtTIData(
592 request==MPI_REQUEST_IGNORED ? "exscan":"iexscan", -1, datatype->is_replayable() ? count : count * datatype->size(),
593 simgrid::smpi::Datatype::encode(datatype)));
595 if(request == MPI_REQUEST_IGNORED)
596 retval = simgrid::smpi::Colls::exscan(sendtmpbuf, recvbuf, count, datatype, op, comm);
598 retval = simgrid::smpi::Colls::iexscan(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
600 TRACE_smpi_comm_out(rank);
601 if (sendbuf == MPI_IN_PLACE)
602 xbt_free(sendtmpbuf);
609 int PMPI_Reduce_scatter(void *sendbuf, void *recvbuf, int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
611 return PMPI_Ireduce_scatter(sendbuf, recvbuf, recvcounts, datatype, op, comm, MPI_REQUEST_IGNORED);
614 int PMPI_Ireduce_scatter(void *sendbuf, void *recvbuf, int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
618 if (comm == MPI_COMM_NULL) {
619 retval = MPI_ERR_COMM;
620 } else if ((sendbuf == nullptr) || (recvbuf == nullptr)) {
621 retval = MPI_ERR_BUFFER;
622 } else if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()){
623 retval = MPI_ERR_TYPE;
624 } else if (op == MPI_OP_NULL) {
626 } else if (recvcounts == nullptr) {
627 retval = MPI_ERR_ARG;
628 } else if (request == nullptr){
629 retval = MPI_ERR_ARG;
631 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
633 return MPI_ERR_COUNT;
637 int rank = simgrid::s4u::this_actor::get_pid();
638 std::vector<int>* trace_recvcounts = new std::vector<int>;
639 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
642 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
643 trace_recvcounts->push_back(recvcounts[i] * dt_send_size);
644 totalcount += recvcounts[i];
647 void* sendtmpbuf = sendbuf;
648 if (sendbuf == MPI_IN_PLACE) {
649 sendtmpbuf = static_cast<void*>(xbt_malloc(totalcount * datatype->size()));
650 memcpy(sendtmpbuf, recvbuf, totalcount * datatype->size());
653 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter": "PMPI_Ireduce_scatter", new simgrid::instr::VarCollTIData(
654 request==MPI_REQUEST_IGNORED ? "reducescatter":"ireducescatter", -1, dt_send_size, nullptr, -1, trace_recvcounts,
655 simgrid::smpi::Datatype::encode(datatype), ""));
657 if(request == MPI_REQUEST_IGNORED)
658 simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
660 simgrid::smpi::Colls::ireduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm, request);
662 retval = MPI_SUCCESS;
663 TRACE_smpi_comm_out(rank);
665 if (sendbuf == MPI_IN_PLACE)
666 xbt_free(sendtmpbuf);
673 int PMPI_Reduce_scatter_block(void *sendbuf, void *recvbuf, int recvcount,
674 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
676 return PMPI_Ireduce_scatter_block(sendbuf, recvbuf, recvcount, datatype, op, comm, MPI_REQUEST_IGNORED);
679 int PMPI_Ireduce_scatter_block(void *sendbuf, void *recvbuf, int recvcount,
680 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
685 if (comm == MPI_COMM_NULL) {
686 retval = MPI_ERR_COMM;
687 } else if (not datatype->is_valid()) {
688 retval = MPI_ERR_TYPE;
689 } else if (op == MPI_OP_NULL) {
691 } else if (recvcount < 0) {
692 retval = MPI_ERR_ARG;
693 } else if (request == nullptr){
694 retval = MPI_ERR_ARG;
696 int count = comm->size();
698 int rank = simgrid::s4u::this_actor::get_pid();
699 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
700 std::vector<int>* trace_recvcounts = new std::vector<int>(recvcount * dt_send_size); // copy data to avoid bad free
702 void* sendtmpbuf = sendbuf;
703 if (sendbuf == MPI_IN_PLACE) {
704 sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * count * datatype->size()));
705 memcpy(sendtmpbuf, recvbuf, recvcount * count * datatype->size());
708 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter_block":"PMPI_Ireduce_scatter_block",
709 new simgrid::instr::VarCollTIData(request==MPI_REQUEST_IGNORED ? "reducescatter":"ireducescatter", -1, 0, nullptr, -1, trace_recvcounts,
710 simgrid::smpi::Datatype::encode(datatype), ""));
712 int* recvcounts = new int[count];
713 for (int i = 0; i < count; i++)
714 recvcounts[i] = recvcount;
715 if(request == MPI_REQUEST_IGNORED)
716 simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
718 simgrid::smpi::Colls::ireduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm, request);
720 retval = MPI_SUCCESS;
722 TRACE_smpi_comm_out(rank);
724 if (sendbuf == MPI_IN_PLACE)
725 xbt_free(sendtmpbuf);
731 int PMPI_Alltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
732 MPI_Datatype recvtype, MPI_Comm comm){
733 return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
736 int PMPI_Ialltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
737 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request *request)
742 if (comm == MPI_COMM_NULL) {
743 retval = MPI_ERR_COMM;
744 } else if ((sendbuf == nullptr && sendcount > 0) || (recvbuf == nullptr && recvcount > 0)) {
745 retval = MPI_ERR_BUFFER;
746 } else if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL) {
747 retval = MPI_ERR_TYPE;
748 } else if ((sendbuf != MPI_IN_PLACE && sendcount < 0) || recvcount < 0){
749 retval = MPI_ERR_COUNT;
750 } else if (request == nullptr){
751 retval = MPI_ERR_ARG;
753 int rank = simgrid::s4u::this_actor::get_pid();
754 void* sendtmpbuf = static_cast<char*>(sendbuf);
755 int sendtmpcount = sendcount;
756 MPI_Datatype sendtmptype = sendtype;
757 if (sendbuf == MPI_IN_PLACE) {
758 sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * comm->size() * recvtype->size()));
759 memcpy(sendtmpbuf, recvbuf, recvcount * comm->size() * recvtype->size());
760 sendtmpcount = recvcount;
761 sendtmptype = recvtype;
764 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Alltoall":"PMPI_Ialltoall",
765 new simgrid::instr::CollTIData(
766 request==MPI_REQUEST_IGNORED ? "alltoall" : "ialltoall", -1, -1.0,
767 sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(),
768 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
769 simgrid::smpi::Datatype::encode(sendtmptype), simgrid::smpi::Datatype::encode(recvtype)));
770 if(request == MPI_REQUEST_IGNORED)
771 retval = simgrid::smpi::Colls::alltoall(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, comm);
773 retval = simgrid::smpi::Colls::ialltoall(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, comm, request);
775 TRACE_smpi_comm_out(rank);
777 if (sendbuf == MPI_IN_PLACE)
778 xbt_free(sendtmpbuf);
785 int PMPI_Alltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf,
786 int* recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm)
788 return PMPI_Ialltoallv(sendbuf, sendcounts, senddisps, sendtype, recvbuf, recvcounts, recvdisps, recvtype, comm, MPI_REQUEST_IGNORED);
791 int PMPI_Ialltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf,
792 int* recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request *request)
795 if (comm == MPI_COMM_NULL) {
796 retval = MPI_ERR_COMM;
797 } else if (sendbuf == nullptr || recvbuf == nullptr) {
798 retval = MPI_ERR_BUFFER;
799 } else if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL) {
800 retval = MPI_ERR_TYPE;
801 } else if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
802 recvdisps == nullptr) {
803 retval = MPI_ERR_ARG;
804 } else if (request == nullptr){
805 retval = MPI_ERR_ARG;
807 int rank = simgrid::s4u::this_actor::get_pid();
808 int size = comm->size();
809 for (int i = 0; i < size; i++) {
810 if (recvcounts[i] <0 || (sendbuf != MPI_IN_PLACE && sendcounts[i]<0))
811 return MPI_ERR_COUNT;
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 recv_size += recvcounts[i] * dt_size_recv;
827 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
828 if (((recvdisps[i] + recvcounts[i]) * dt_size_recv) > maxsize)
829 maxsize = (recvdisps[i] + recvcounts[i]) * dt_size_recv;
832 if (sendbuf == MPI_IN_PLACE) {
833 sendtmpbuf = static_cast<void*>(xbt_malloc(maxsize));
834 memcpy(sendtmpbuf, recvbuf, maxsize);
835 sendtmpcounts = static_cast<int*>(xbt_malloc(size * sizeof(int)));
836 memcpy(sendtmpcounts, recvcounts, size * sizeof(int));
837 sendtmpdisps = static_cast<int*>(xbt_malloc(size * sizeof(int)));
838 memcpy(sendtmpdisps, recvdisps, size * sizeof(int));
839 sendtmptype = recvtype;
842 int dt_size_send = sendtmptype->size();
844 for (int i = 0; i < size; i++) { // copy data to avoid bad free
845 send_size += sendtmpcounts[i] * dt_size_send;
846 trace_sendcounts->push_back(sendtmpcounts[i] * dt_size_send);
849 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Alltoallv":"PMPI_Ialltoallv",
850 new simgrid::instr::VarCollTIData(request==MPI_REQUEST_IGNORED ? "alltoallv":"ialltoallv", -1, send_size, trace_sendcounts, recv_size,
851 trace_recvcounts, simgrid::smpi::Datatype::encode(sendtmptype),
852 simgrid::smpi::Datatype::encode(recvtype)));
854 if(request == MPI_REQUEST_IGNORED)
855 retval = simgrid::smpi::Colls::alltoallv(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptype, recvbuf, recvcounts,
856 recvdisps, recvtype, comm);
858 retval = simgrid::smpi::Colls::ialltoallv(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptype, recvbuf, recvcounts,
859 recvdisps, recvtype, comm, request);
860 TRACE_smpi_comm_out(rank);
862 if (sendbuf == MPI_IN_PLACE) {
863 xbt_free(sendtmpbuf);
864 xbt_free(sendtmpcounts);
865 xbt_free(sendtmpdisps);
872 int PMPI_Alltoallw(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype* sendtypes, void* recvbuf,
873 int* recvcounts, int* recvdisps, MPI_Datatype* recvtypes, MPI_Comm comm)
875 return PMPI_Ialltoallw(sendbuf, sendcounts, senddisps, sendtypes, recvbuf, recvcounts, recvdisps, recvtypes, comm, MPI_REQUEST_IGNORED);
878 int PMPI_Ialltoallw(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype* sendtypes, void* recvbuf,
879 int* recvcounts, int* recvdisps, MPI_Datatype* recvtypes, MPI_Comm comm, MPI_Request *request)
882 if (comm == MPI_COMM_NULL) {
883 retval = MPI_ERR_COMM;
884 } else if (sendbuf == nullptr || recvbuf == nullptr) {
885 retval = MPI_ERR_BUFFER;
886 } else if ((sendbuf != MPI_IN_PLACE && sendtypes == nullptr) || recvtypes == nullptr) {
887 retval = MPI_ERR_TYPE;
888 } else if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
889 recvdisps == nullptr) {
890 retval = MPI_ERR_ARG;
891 } else if (request == nullptr){
892 retval = MPI_ERR_ARG;
895 int rank = simgrid::s4u::this_actor::get_pid();
896 int size = comm->size();
897 for (int i = 0; i < size; i++) { // copy data to avoid bad free
898 if (recvcounts[i] <0 || (sendbuf != MPI_IN_PLACE && sendcounts[i]<0))
899 return MPI_ERR_COUNT;
903 std::vector<int>* trace_sendcounts = new std::vector<int>;
904 std::vector<int>* trace_recvcounts = new std::vector<int>;
906 void* sendtmpbuf = static_cast<char*>(sendbuf);
907 int* sendtmpcounts = sendcounts;
908 int* sendtmpdisps = senddisps;
909 MPI_Datatype* sendtmptypes = sendtypes;
910 unsigned long maxsize = 0;
911 for (int i = 0; i < size; i++) { // copy data to avoid bad free
912 if(recvtypes[i]==MPI_DATATYPE_NULL){
913 delete trace_recvcounts;
914 delete trace_sendcounts;
917 recv_size += recvcounts[i] * recvtypes[i]->size();
918 trace_recvcounts->push_back(recvcounts[i] * recvtypes[i]->size());
919 if ((recvdisps[i] + (recvcounts[i] * recvtypes[i]->size())) > maxsize)
920 maxsize = recvdisps[i] + (recvcounts[i] * recvtypes[i]->size());
923 if (sendbuf == MPI_IN_PLACE) {
924 sendtmpbuf = static_cast<void*>(xbt_malloc(maxsize));
925 memcpy(sendtmpbuf, recvbuf, maxsize);
926 sendtmpcounts = static_cast<int*>(xbt_malloc(size * sizeof(int)));
927 memcpy(sendtmpcounts, recvcounts, size * sizeof(int));
928 sendtmpdisps = static_cast<int*>(xbt_malloc(size * sizeof(int)));
929 memcpy(sendtmpdisps, recvdisps, size * sizeof(int));
930 sendtmptypes = static_cast<MPI_Datatype*>(xbt_malloc(size * sizeof(MPI_Datatype)));
931 memcpy(sendtmptypes, recvtypes, size * sizeof(MPI_Datatype));
934 for (int i = 0; i < size; i++) { // copy data to avoid bad free
935 send_size += sendtmpcounts[i] * sendtmptypes[i]->size();
936 trace_sendcounts->push_back(sendtmpcounts[i] * sendtmptypes[i]->size());
939 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Alltoallw":"PMPI_Ialltoallw",
940 new simgrid::instr::VarCollTIData(request==MPI_REQUEST_IGNORED ? "alltoallv":"ialltoallv", -1, send_size, trace_sendcounts, recv_size,
941 trace_recvcounts, simgrid::smpi::Datatype::encode(sendtmptypes[0]),
942 simgrid::smpi::Datatype::encode(recvtypes[0])));
944 if(request == MPI_REQUEST_IGNORED)
945 retval = simgrid::smpi::Colls::alltoallw(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptypes, recvbuf, recvcounts,
946 recvdisps, recvtypes, comm);
948 retval = simgrid::smpi::Colls::ialltoallw(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptypes, recvbuf, recvcounts,
949 recvdisps, recvtypes, comm, request);
950 TRACE_smpi_comm_out(rank);
952 if (sendbuf == MPI_IN_PLACE) {
953 xbt_free(sendtmpbuf);
954 xbt_free(sendtmpcounts);
955 xbt_free(sendtmpdisps);
956 xbt_free(sendtmptypes);