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) {
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) && (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) || ((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) || (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) || (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)) || (recvbuf == nullptr)){
320 return MPI_ERR_BUFFER;
321 } else if (root < 0 || root >= comm->size()){
323 } else if (request == nullptr){
327 if (recvbuf == MPI_IN_PLACE) {
329 recvcount = sendcount;
331 int rank = simgrid::s4u::this_actor::get_pid();
334 rank, request==MPI_REQUEST_IGNORED?"PMPI_Scatter":"PMPI_Iscatter",
335 new simgrid::instr::CollTIData(
336 request==MPI_REQUEST_IGNORED ? "scatter" : "iscatter", root, -1.0,
337 (comm->rank() != root || sendtype->is_replayable()) ? sendcount : sendcount * sendtype->size(),
338 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
339 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
340 if(request == MPI_REQUEST_IGNORED)
341 simgrid::smpi::Colls::scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
343 simgrid::smpi::Colls::iscatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
344 TRACE_smpi_comm_out(rank);
350 int PMPI_Scatterv(void *sendbuf, int *sendcounts, int *displs,
351 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm){
352 return PMPI_Iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, MPI_REQUEST_IGNORED);
355 int PMPI_Iscatterv(void *sendbuf, int *sendcounts, int *displs,
356 MPI_Datatype sendtype, void *recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm, MPI_Request *request)
358 if (comm == MPI_COMM_NULL) {
360 } else if (sendcounts == nullptr || displs == nullptr) {
362 } else if (((comm->rank() == root) && (sendtype == MPI_DATATYPE_NULL)) ||
363 ((recvbuf != MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL))) {
365 } else if (request == nullptr){
367 } else if (recvbuf != MPI_IN_PLACE && recvcount < 0){
368 return MPI_ERR_COUNT;
369 } else if (root < 0 || root >= comm->size()){
372 if (comm->rank() == root){
373 if(recvbuf == MPI_IN_PLACE) {
375 recvcount = sendcounts[comm->rank()];
377 for (int i = 0; i < comm->size(); i++){
379 return MPI_ERR_COUNT;
384 int rank = simgrid::s4u::this_actor::get_pid();
385 int dt_size_send = sendtype->is_replayable() ? 1 : sendtype->size();
387 std::vector<int>* trace_sendcounts = new std::vector<int>;
388 if (comm->rank() == root) {
389 for (int i = 0; i < comm->size(); i++){ // copy data to avoid bad free
390 trace_sendcounts->push_back(sendcounts[i] * dt_size_send);
394 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Scatterv":"PMPI_Iscatterv",
395 new simgrid::instr::VarCollTIData(
396 request==MPI_REQUEST_IGNORED ? "scatterv":"iscatterv", root, dt_size_send, trace_sendcounts,
397 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(), nullptr,
398 simgrid::smpi::Datatype::encode(sendtype), simgrid::smpi::Datatype::encode(recvtype)));
399 if(request == MPI_REQUEST_IGNORED)
400 simgrid::smpi::Colls::scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
402 simgrid::smpi::Colls::iscatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm, request);
403 TRACE_smpi_comm_out(rank);
410 int PMPI_Reduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
412 return PMPI_Ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, MPI_REQUEST_IGNORED);
415 int PMPI_Ireduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm, MPI_Request* request)
417 if (comm == MPI_COMM_NULL) {
419 } if ((sendbuf == nullptr) || ((comm->rank() == root) && recvbuf == nullptr)) {
420 return MPI_ERR_BUFFER;
421 } else if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()){
423 } else if (op == MPI_OP_NULL) {
425 } else if (request == nullptr){
427 } else if (root < 0 || root >= comm->size()){
429 } else if (count < 0){
430 return MPI_ERR_COUNT;
433 int rank = simgrid::s4u::this_actor::get_pid();
435 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Reduce":"PMPI_Ireduce",
436 new simgrid::instr::CollTIData(request==MPI_REQUEST_IGNORED ? "reduce":"ireduce", root, 0,
437 datatype->is_replayable() ? count : count * datatype->size(), -1,
438 simgrid::smpi::Datatype::encode(datatype), ""));
439 if(request == MPI_REQUEST_IGNORED)
440 simgrid::smpi::Colls::reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
442 simgrid::smpi::Colls::ireduce(sendbuf, recvbuf, count, datatype, op, root, comm, request);
444 TRACE_smpi_comm_out(rank);
450 int PMPI_Reduce_local(void *inbuf, void *inoutbuf, int count, MPI_Datatype datatype, MPI_Op op){
454 if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()){
455 retval = MPI_ERR_TYPE;
456 } else if (op == MPI_OP_NULL) {
458 } else if (count < 0){
459 retval = MPI_ERR_COUNT;
461 op->apply(inbuf, inoutbuf, &count, datatype);
462 retval = MPI_SUCCESS;
468 int PMPI_Allreduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
470 return PMPI_Iallreduce(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
473 int PMPI_Iallreduce(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
475 if (comm == MPI_COMM_NULL) {
477 } if ((sendbuf == nullptr) || (recvbuf == nullptr)) {
478 return MPI_ERR_BUFFER;
479 } else if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()) {
481 } else if (count < 0){
482 return MPI_ERR_COUNT;
483 } else if (op == MPI_OP_NULL) {
485 } else if (request == nullptr){
489 char* sendtmpbuf = static_cast<char*>(sendbuf);
490 if( sendbuf == MPI_IN_PLACE ) {
491 sendtmpbuf = static_cast<char*>(xbt_malloc(count*datatype->get_extent()));
492 simgrid::smpi::Datatype::copy(recvbuf, count, datatype,sendtmpbuf, count, datatype);
494 int rank = simgrid::s4u::this_actor::get_pid();
496 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ?"PMPI_Allreduce":"PMPI_Iallreduce",
497 new simgrid::instr::CollTIData(request==MPI_REQUEST_IGNORED ? "allreduce":"iallreduce", -1, 0,
498 datatype->is_replayable() ? count : count * datatype->size(), -1,
499 simgrid::smpi::Datatype::encode(datatype), ""));
501 if(request == MPI_REQUEST_IGNORED)
502 simgrid::smpi::Colls::allreduce(sendtmpbuf, recvbuf, count, datatype, op, comm);
504 simgrid::smpi::Colls::iallreduce(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
506 if( sendbuf == MPI_IN_PLACE )
507 xbt_free(sendtmpbuf);
509 TRACE_smpi_comm_out(rank);
515 int PMPI_Scan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
517 return PMPI_Iscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
520 int PMPI_Iscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request)
526 if (comm == MPI_COMM_NULL) {
527 retval = MPI_ERR_COMM;
528 } else if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()){
529 retval = MPI_ERR_TYPE;
530 } else if (op == MPI_OP_NULL) {
532 } else if (request == nullptr){
533 retval = MPI_ERR_ARG;
534 } else if (count < 0){
535 retval = MPI_ERR_COUNT;
537 int rank = simgrid::s4u::this_actor::get_pid();
538 void* sendtmpbuf = sendbuf;
539 if (sendbuf == MPI_IN_PLACE) {
540 sendtmpbuf = static_cast<void*>(xbt_malloc(count * datatype->size()));
541 memcpy(sendtmpbuf, recvbuf, count * datatype->size());
543 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Scan" : "PMPI_Iscan", new simgrid::instr::Pt2PtTIData(
544 request==MPI_REQUEST_IGNORED ? "scan":"iscan", -1,
545 datatype->is_replayable() ? count : count * datatype->size(),
546 simgrid::smpi::Datatype::encode(datatype)));
548 if(request == MPI_REQUEST_IGNORED)
549 retval = simgrid::smpi::Colls::scan(sendtmpbuf, recvbuf, count, datatype, op, comm);
551 retval = simgrid::smpi::Colls::iscan(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
553 TRACE_smpi_comm_out(rank);
554 if (sendbuf == MPI_IN_PLACE)
555 xbt_free(sendtmpbuf);
562 int PMPI_Exscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
564 return PMPI_Iexscan(sendbuf, recvbuf, count, datatype, op, comm, MPI_REQUEST_IGNORED);
567 int PMPI_Iexscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request* request){
572 if (comm == MPI_COMM_NULL) {
573 retval = MPI_ERR_COMM;
574 } else if (not datatype->is_valid()) {
575 retval = MPI_ERR_TYPE;
576 } else if (op == MPI_OP_NULL) {
578 } else if (request == nullptr){
579 retval = MPI_ERR_ARG;
580 } else if (count < 0){
581 retval = MPI_ERR_COUNT;
583 int rank = simgrid::s4u::this_actor::get_pid();
584 void* sendtmpbuf = sendbuf;
585 if (sendbuf == MPI_IN_PLACE) {
586 sendtmpbuf = static_cast<void*>(xbt_malloc(count * datatype->size()));
587 memcpy(sendtmpbuf, recvbuf, count * datatype->size());
590 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Exscan" : "PMPI_Iexscan", new simgrid::instr::Pt2PtTIData(
591 request==MPI_REQUEST_IGNORED ? "exscan":"iexscan", -1, datatype->is_replayable() ? count : count * datatype->size(),
592 simgrid::smpi::Datatype::encode(datatype)));
594 if(request == MPI_REQUEST_IGNORED)
595 retval = simgrid::smpi::Colls::exscan(sendtmpbuf, recvbuf, count, datatype, op, comm);
597 retval = simgrid::smpi::Colls::iexscan(sendtmpbuf, recvbuf, count, datatype, op, comm, request);
599 TRACE_smpi_comm_out(rank);
600 if (sendbuf == MPI_IN_PLACE)
601 xbt_free(sendtmpbuf);
608 int PMPI_Reduce_scatter(void *sendbuf, void *recvbuf, int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
610 return PMPI_Ireduce_scatter(sendbuf, recvbuf, recvcounts, datatype, op, comm, MPI_REQUEST_IGNORED);
613 int PMPI_Ireduce_scatter(void *sendbuf, void *recvbuf, int *recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
617 if (comm == MPI_COMM_NULL) {
618 retval = MPI_ERR_COMM;
619 } else if ((sendbuf == nullptr) || (recvbuf == nullptr)) {
620 retval = MPI_ERR_BUFFER;
621 } else if (datatype == MPI_DATATYPE_NULL || not datatype->is_valid()){
622 retval = MPI_ERR_TYPE;
623 } else if (op == MPI_OP_NULL) {
625 } else if (recvcounts == nullptr) {
626 retval = MPI_ERR_ARG;
627 } else if (request == nullptr){
628 retval = MPI_ERR_ARG;
630 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
632 return MPI_ERR_COUNT;
636 int rank = simgrid::s4u::this_actor::get_pid();
637 std::vector<int>* trace_recvcounts = new std::vector<int>;
638 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
641 for (int i = 0; i < comm->size(); i++) { // copy data to avoid bad free
642 trace_recvcounts->push_back(recvcounts[i] * dt_send_size);
643 totalcount += recvcounts[i];
646 void* sendtmpbuf = sendbuf;
647 if (sendbuf == MPI_IN_PLACE) {
648 sendtmpbuf = static_cast<void*>(xbt_malloc(totalcount * datatype->size()));
649 memcpy(sendtmpbuf, recvbuf, totalcount * datatype->size());
652 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter": "PMPI_Ireduce_scatter", new simgrid::instr::VarCollTIData(
653 request==MPI_REQUEST_IGNORED ? "reducescatter":"ireducescatter", -1, dt_send_size, nullptr, -1, trace_recvcounts,
654 simgrid::smpi::Datatype::encode(datatype), ""));
656 if(request == MPI_REQUEST_IGNORED)
657 simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
659 simgrid::smpi::Colls::ireduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm, request);
661 retval = MPI_SUCCESS;
662 TRACE_smpi_comm_out(rank);
664 if (sendbuf == MPI_IN_PLACE)
665 xbt_free(sendtmpbuf);
672 int PMPI_Reduce_scatter_block(void *sendbuf, void *recvbuf, int recvcount,
673 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
675 return PMPI_Ireduce_scatter_block(sendbuf, recvbuf, recvcount, datatype, op, comm, MPI_REQUEST_IGNORED);
678 int PMPI_Ireduce_scatter_block(void *sendbuf, void *recvbuf, int recvcount,
679 MPI_Datatype datatype, MPI_Op op, MPI_Comm comm, MPI_Request *request)
684 if (comm == MPI_COMM_NULL) {
685 retval = MPI_ERR_COMM;
686 } else if (not datatype->is_valid()) {
687 retval = MPI_ERR_TYPE;
688 } else if (op == MPI_OP_NULL) {
690 } else if (recvcount < 0) {
691 retval = MPI_ERR_ARG;
692 } else if (request == nullptr){
693 retval = MPI_ERR_ARG;
695 int count = comm->size();
697 int rank = simgrid::s4u::this_actor::get_pid();
698 int dt_send_size = datatype->is_replayable() ? 1 : datatype->size();
699 std::vector<int>* trace_recvcounts = new std::vector<int>(recvcount * dt_send_size); // copy data to avoid bad free
701 void* sendtmpbuf = sendbuf;
702 if (sendbuf == MPI_IN_PLACE) {
703 sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * count * datatype->size()));
704 memcpy(sendtmpbuf, recvbuf, recvcount * count * datatype->size());
707 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED ? "PMPI_Reduce_scatter_block":"PMPI_Ireduce_scatter_block",
708 new simgrid::instr::VarCollTIData(request==MPI_REQUEST_IGNORED ? "reducescatter":"ireducescatter", -1, 0, nullptr, -1, trace_recvcounts,
709 simgrid::smpi::Datatype::encode(datatype), ""));
711 int* recvcounts = new int[count];
712 for (int i = 0; i < count; i++)
713 recvcounts[i] = recvcount;
714 if(request == MPI_REQUEST_IGNORED)
715 simgrid::smpi::Colls::reduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm);
717 simgrid::smpi::Colls::ireduce_scatter(sendtmpbuf, recvbuf, recvcounts, datatype, op, comm, request);
719 retval = MPI_SUCCESS;
721 TRACE_smpi_comm_out(rank);
723 if (sendbuf == MPI_IN_PLACE)
724 xbt_free(sendtmpbuf);
730 int PMPI_Alltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
731 MPI_Datatype recvtype, MPI_Comm comm){
732 return PMPI_Ialltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm, MPI_REQUEST_IGNORED);
735 int PMPI_Ialltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount,
736 MPI_Datatype recvtype, MPI_Comm comm, MPI_Request *request)
741 if (comm == MPI_COMM_NULL) {
742 retval = MPI_ERR_COMM;
743 } else if (sendbuf == nullptr || recvbuf == nullptr) {
744 retval = MPI_ERR_BUFFER;
745 } else if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL) {
746 retval = MPI_ERR_TYPE;
747 } else if ((sendbuf != MPI_IN_PLACE && sendcount < 0) || recvcount < 0){
748 retval = MPI_ERR_COUNT;
749 } else if (request == nullptr){
750 retval = MPI_ERR_ARG;
752 int rank = simgrid::s4u::this_actor::get_pid();
753 void* sendtmpbuf = static_cast<char*>(sendbuf);
754 int sendtmpcount = sendcount;
755 MPI_Datatype sendtmptype = sendtype;
756 if (sendbuf == MPI_IN_PLACE) {
757 sendtmpbuf = static_cast<void*>(xbt_malloc(recvcount * comm->size() * recvtype->size()));
758 memcpy(sendtmpbuf, recvbuf, recvcount * comm->size() * recvtype->size());
759 sendtmpcount = recvcount;
760 sendtmptype = recvtype;
763 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Alltoall":"PMPI_Ialltoall",
764 new simgrid::instr::CollTIData(
765 request==MPI_REQUEST_IGNORED ? "alltoall" : "ialltoall", -1, -1.0,
766 sendtmptype->is_replayable() ? sendtmpcount : sendtmpcount * sendtmptype->size(),
767 recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(),
768 simgrid::smpi::Datatype::encode(sendtmptype), simgrid::smpi::Datatype::encode(recvtype)));
769 if(request == MPI_REQUEST_IGNORED)
770 retval = simgrid::smpi::Colls::alltoall(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, comm);
772 retval = simgrid::smpi::Colls::ialltoall(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount, recvtype, comm, request);
774 TRACE_smpi_comm_out(rank);
776 if (sendbuf == MPI_IN_PLACE)
777 xbt_free(sendtmpbuf);
784 int PMPI_Alltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf,
785 int* recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm)
787 return PMPI_Ialltoallv(sendbuf, sendcounts, senddisps, sendtype, recvbuf, recvcounts, recvdisps, recvtype, comm, MPI_REQUEST_IGNORED);
790 int PMPI_Ialltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf,
791 int* recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm, MPI_Request *request)
794 if (comm == MPI_COMM_NULL) {
795 retval = MPI_ERR_COMM;
796 } else if (sendbuf == nullptr || recvbuf == nullptr) {
797 retval = MPI_ERR_BUFFER;
798 } else if ((sendbuf != MPI_IN_PLACE && sendtype == MPI_DATATYPE_NULL) || recvtype == MPI_DATATYPE_NULL) {
799 retval = MPI_ERR_TYPE;
800 } else if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
801 recvdisps == nullptr) {
802 retval = MPI_ERR_ARG;
803 } else if (request == nullptr){
804 retval = MPI_ERR_ARG;
806 int rank = simgrid::s4u::this_actor::get_pid();
807 int size = comm->size();
808 for (int i = 0; i < size; i++) {
809 if (recvcounts[i] <0 || (sendbuf != MPI_IN_PLACE && sendcounts[i]<0))
810 return MPI_ERR_COUNT;
815 std::vector<int>* trace_sendcounts = new std::vector<int>;
816 std::vector<int>* trace_recvcounts = new std::vector<int>;
817 int dt_size_recv = recvtype->size();
819 void* sendtmpbuf = static_cast<char*>(sendbuf);
820 int* sendtmpcounts = sendcounts;
821 int* sendtmpdisps = senddisps;
822 MPI_Datatype sendtmptype = sendtype;
824 for (int i = 0; i < size; i++) { // copy data to avoid bad free
825 recv_size += recvcounts[i] * dt_size_recv;
826 trace_recvcounts->push_back(recvcounts[i] * dt_size_recv);
827 if (((recvdisps[i] + recvcounts[i]) * dt_size_recv) > maxsize)
828 maxsize = (recvdisps[i] + recvcounts[i]) * dt_size_recv;
831 if (sendbuf == MPI_IN_PLACE) {
832 sendtmpbuf = static_cast<void*>(xbt_malloc(maxsize));
833 memcpy(sendtmpbuf, recvbuf, maxsize);
834 sendtmpcounts = static_cast<int*>(xbt_malloc(size * sizeof(int)));
835 memcpy(sendtmpcounts, recvcounts, size * sizeof(int));
836 sendtmpdisps = static_cast<int*>(xbt_malloc(size * sizeof(int)));
837 memcpy(sendtmpdisps, recvdisps, size * sizeof(int));
838 sendtmptype = recvtype;
841 int dt_size_send = sendtmptype->size();
843 for (int i = 0; i < size; i++) { // copy data to avoid bad free
844 send_size += sendtmpcounts[i] * dt_size_send;
845 trace_sendcounts->push_back(sendtmpcounts[i] * dt_size_send);
848 TRACE_smpi_comm_in(rank, request==MPI_REQUEST_IGNORED?"PMPI_Alltoallv":"PMPI_Ialltoallv",
849 new simgrid::instr::VarCollTIData(request==MPI_REQUEST_IGNORED ? "alltoallv":"ialltoallv", -1, send_size, trace_sendcounts, recv_size,
850 trace_recvcounts, simgrid::smpi::Datatype::encode(sendtmptype),
851 simgrid::smpi::Datatype::encode(recvtype)));
853 if(request == MPI_REQUEST_IGNORED)
854 retval = simgrid::smpi::Colls::alltoallv(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptype, recvbuf, recvcounts,
855 recvdisps, recvtype, comm);
857 retval = simgrid::smpi::Colls::ialltoallv(sendtmpbuf, sendtmpcounts, sendtmpdisps, sendtmptype, recvbuf, recvcounts,
858 recvdisps, recvtype, comm, request);
859 TRACE_smpi_comm_out(rank);
861 if (sendbuf == MPI_IN_PLACE) {
862 xbt_free(sendtmpbuf);
863 xbt_free(sendtmpcounts);
864 xbt_free(sendtmpdisps);
871 int PMPI_Alltoallw(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype* sendtypes, void* recvbuf,
872 int* recvcounts, int* recvdisps, MPI_Datatype* recvtypes, MPI_Comm comm)
874 return PMPI_Ialltoallw(sendbuf, sendcounts, senddisps, sendtypes, recvbuf, recvcounts, recvdisps, recvtypes, comm, MPI_REQUEST_IGNORED);
877 int PMPI_Ialltoallw(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype* sendtypes, void* recvbuf,
878 int* recvcounts, int* recvdisps, MPI_Datatype* recvtypes, MPI_Comm comm, MPI_Request *request)
881 if (comm == MPI_COMM_NULL) {
882 retval = MPI_ERR_COMM;
883 } else if (sendbuf == nullptr || recvbuf == nullptr) {
884 retval = MPI_ERR_BUFFER;
885 } else if ((sendbuf != MPI_IN_PLACE && sendtypes == nullptr) || recvtypes == nullptr) {
886 retval = MPI_ERR_TYPE;
887 } else if ((sendbuf != MPI_IN_PLACE && (sendcounts == nullptr || senddisps == nullptr)) || recvcounts == nullptr ||
888 recvdisps == nullptr) {
889 retval = MPI_ERR_ARG;
890 } else if (request == nullptr){
891 retval = MPI_ERR_ARG;
894 int rank = simgrid::s4u::this_actor::get_pid();
895 int size = comm->size();
896 for (int i = 0; i < size; i++) { // copy data to avoid bad free
897 if (recvcounts[i] <0 || (sendbuf != MPI_IN_PLACE && sendcounts[i]<0))
898 return MPI_ERR_COUNT;
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);