1 /* Copyright (c) 2007-2017. 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_comm.hpp"
8 #include "smpi_datatype.hpp"
9 #include "smpi_process.hpp"
10 #include "smpi_request.hpp"
12 XBT_LOG_EXTERNAL_DEFAULT_CATEGORY(smpi_pmpi);
14 /* PMPI User level calls */
15 extern "C" { // Obviously, the C MPI interface should use the C linkage
17 int PMPI_Send_init(void *buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm, MPI_Request * request)
22 if (request == nullptr) {
24 } else if (comm == MPI_COMM_NULL) {
25 retval = MPI_ERR_COMM;
26 } else if (not datatype->is_valid()) {
27 retval = MPI_ERR_TYPE;
28 } else if (dst == MPI_PROC_NULL) {
31 *request = simgrid::smpi::Request::send_init(buf, count, datatype, dst, tag, comm);
35 if (retval != MPI_SUCCESS && request != nullptr)
36 *request = MPI_REQUEST_NULL;
40 int PMPI_Recv_init(void *buf, int count, MPI_Datatype datatype, int src, int tag, MPI_Comm comm, MPI_Request * request)
45 if (request == nullptr) {
47 } else if (comm == MPI_COMM_NULL) {
48 retval = MPI_ERR_COMM;
49 } else if (not datatype->is_valid()) {
50 retval = MPI_ERR_TYPE;
51 } else if (src == MPI_PROC_NULL) {
54 *request = simgrid::smpi::Request::recv_init(buf, count, datatype, src, tag, comm);
58 if (retval != MPI_SUCCESS && request != nullptr)
59 *request = MPI_REQUEST_NULL;
63 int PMPI_Ssend_init(void* buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm, MPI_Request* request)
68 if (request == nullptr) {
70 } else if (comm == MPI_COMM_NULL) {
71 retval = MPI_ERR_COMM;
72 } else if (not datatype->is_valid()) {
73 retval = MPI_ERR_TYPE;
74 } else if (dst == MPI_PROC_NULL) {
77 *request = simgrid::smpi::Request::ssend_init(buf, count, datatype, dst, tag, comm);
81 if (retval != MPI_SUCCESS && request != nullptr)
82 *request = MPI_REQUEST_NULL;
86 int PMPI_Start(MPI_Request * request)
91 if (request == nullptr || *request == MPI_REQUEST_NULL) {
92 retval = MPI_ERR_REQUEST;
101 int PMPI_Startall(int count, MPI_Request * requests)
105 if (requests == nullptr) {
106 retval = MPI_ERR_ARG;
108 retval = MPI_SUCCESS;
109 for (int i = 0; i < count; i++) {
110 if(requests[i] == MPI_REQUEST_NULL) {
111 retval = MPI_ERR_REQUEST;
114 if(retval != MPI_ERR_REQUEST) {
115 simgrid::smpi::Request::startall(count, requests);
122 int PMPI_Request_free(MPI_Request * request)
127 if (*request == MPI_REQUEST_NULL) {
128 retval = MPI_ERR_ARG;
130 simgrid::smpi::Request::unref(request);
131 retval = MPI_SUCCESS;
137 int PMPI_Irecv(void *buf, int count, MPI_Datatype datatype, int src, int tag, MPI_Comm comm, MPI_Request * request)
143 if (request == nullptr) {
144 retval = MPI_ERR_ARG;
145 } else if (comm == MPI_COMM_NULL) {
146 retval = MPI_ERR_COMM;
147 } else if (src == MPI_PROC_NULL) {
148 *request = MPI_REQUEST_NULL;
149 retval = MPI_SUCCESS;
150 } else if (src!=MPI_ANY_SOURCE && (src >= comm->group()->size() || src <0)){
151 retval = MPI_ERR_RANK;
152 } else if ((count < 0) || (buf==nullptr && count > 0)) {
153 retval = MPI_ERR_COUNT;
154 } else if (not datatype->is_valid()) {
155 retval = MPI_ERR_TYPE;
156 } else if(tag<0 && tag != MPI_ANY_TAG){
157 retval = MPI_ERR_TAG;
160 int rank = smpi_process()->index();
162 instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
163 extra->type = TRACING_IRECV;
164 extra->src = comm->group()->index(src);
166 extra->datatype1 = encode_datatype(datatype);
167 extra->send_size = datatype->is_basic() ? count : count * datatype->size();
168 TRACE_smpi_comm_in(rank, __FUNCTION__, extra);
170 *request = simgrid::smpi::Request::irecv(buf, count, datatype, src, tag, comm);
171 retval = MPI_SUCCESS;
173 TRACE_smpi_comm_out(rank);
177 if (retval != MPI_SUCCESS && request != nullptr)
178 *request = MPI_REQUEST_NULL;
183 int PMPI_Isend(void *buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm, MPI_Request * request)
188 if (request == nullptr) {
189 retval = MPI_ERR_ARG;
190 } else if (comm == MPI_COMM_NULL) {
191 retval = MPI_ERR_COMM;
192 } else if (dst == MPI_PROC_NULL) {
193 *request = MPI_REQUEST_NULL;
194 retval = MPI_SUCCESS;
195 } else if (dst >= comm->group()->size() || dst <0){
196 retval = MPI_ERR_RANK;
197 } else if ((count < 0) || (buf==nullptr && count > 0)) {
198 retval = MPI_ERR_COUNT;
199 } else if (not datatype->is_valid()) {
200 retval = MPI_ERR_TYPE;
201 } else if(tag<0 && tag != MPI_ANY_TAG){
202 retval = MPI_ERR_TAG;
204 int rank = smpi_process()->index();
205 int dst_traced = comm->group()->index(dst);
206 instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
207 extra->type = TRACING_ISEND;
209 extra->dst = dst_traced;
210 extra->datatype1 = encode_datatype(datatype);
211 extra->send_size = datatype->is_basic() ? count : count * datatype->size();
212 TRACE_smpi_comm_in(rank, __FUNCTION__, extra);
213 TRACE_smpi_send(rank, rank, dst_traced, tag, count*datatype->size());
215 *request = simgrid::smpi::Request::isend(buf, count, datatype, dst, tag, comm);
216 retval = MPI_SUCCESS;
218 TRACE_smpi_comm_out(rank);
222 if (retval != MPI_SUCCESS && request!=nullptr)
223 *request = MPI_REQUEST_NULL;
227 int PMPI_Issend(void* buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm, MPI_Request* request)
232 if (request == nullptr) {
233 retval = MPI_ERR_ARG;
234 } else if (comm == MPI_COMM_NULL) {
235 retval = MPI_ERR_COMM;
236 } else if (dst == MPI_PROC_NULL) {
237 *request = MPI_REQUEST_NULL;
238 retval = MPI_SUCCESS;
239 } else if (dst >= comm->group()->size() || dst <0){
240 retval = MPI_ERR_RANK;
241 } else if ((count < 0)|| (buf==nullptr && count > 0)) {
242 retval = MPI_ERR_COUNT;
243 } else if (not datatype->is_valid()) {
244 retval = MPI_ERR_TYPE;
245 } else if(tag<0 && tag != MPI_ANY_TAG){
246 retval = MPI_ERR_TAG;
248 int rank = smpi_process()->index();
249 int dst_traced = comm->group()->index(dst);
250 instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
251 extra->type = TRACING_ISSEND;
253 extra->dst = dst_traced;
254 extra->datatype1 = encode_datatype(datatype);
255 extra->send_size = datatype->is_basic() ? count : count * datatype->size();
256 TRACE_smpi_comm_in(rank, __FUNCTION__, extra);
257 TRACE_smpi_send(rank, rank, dst_traced, tag, count*datatype->size());
259 *request = simgrid::smpi::Request::issend(buf, count, datatype, dst, tag, comm);
260 retval = MPI_SUCCESS;
262 TRACE_smpi_comm_out(rank);
266 if (retval != MPI_SUCCESS && request!=nullptr)
267 *request = MPI_REQUEST_NULL;
271 int PMPI_Recv(void *buf, int count, MPI_Datatype datatype, int src, int tag, MPI_Comm comm, MPI_Status * status)
276 if (comm == MPI_COMM_NULL) {
277 retval = MPI_ERR_COMM;
278 } else if (src == MPI_PROC_NULL) {
279 simgrid::smpi::Status::empty(status);
280 status->MPI_SOURCE = MPI_PROC_NULL;
281 retval = MPI_SUCCESS;
282 } else if (src!=MPI_ANY_SOURCE && (src >= comm->group()->size() || src <0)){
283 retval = MPI_ERR_RANK;
284 } else if ((count < 0) || (buf==nullptr && count > 0)) {
285 retval = MPI_ERR_COUNT;
286 } else if (not datatype->is_valid()) {
287 retval = MPI_ERR_TYPE;
288 } else if(tag<0 && tag != MPI_ANY_TAG){
289 retval = MPI_ERR_TAG;
291 int rank = smpi_process()->index();
292 int src_traced = comm->group()->index(src);
293 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
294 extra->type = TRACING_RECV;
295 extra->src = src_traced;
297 extra->datatype1 = encode_datatype(datatype);
298 extra->send_size = datatype->is_basic() ? count : count * datatype->size();
299 TRACE_smpi_comm_in(rank, __FUNCTION__, extra);
301 simgrid::smpi::Request::recv(buf, count, datatype, src, tag, comm, status);
302 retval = MPI_SUCCESS;
304 // the src may not have been known at the beginning of the recv (MPI_ANY_SOURCE)
305 if (status != MPI_STATUS_IGNORE) {
306 src_traced = comm->group()->index(status->MPI_SOURCE);
307 if (not TRACE_smpi_view_internals()) {
308 TRACE_smpi_recv(src_traced, rank, tag);
311 TRACE_smpi_comm_out(rank);
318 int PMPI_Send(void *buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm)
324 if (comm == MPI_COMM_NULL) {
325 retval = MPI_ERR_COMM;
326 } else if (dst == MPI_PROC_NULL) {
327 retval = MPI_SUCCESS;
328 } else if (dst >= comm->group()->size() || dst <0){
329 retval = MPI_ERR_RANK;
330 } else if ((count < 0) || (buf == nullptr && count > 0)) {
331 retval = MPI_ERR_COUNT;
332 } else if (not datatype->is_valid()) {
333 retval = MPI_ERR_TYPE;
334 } else if(tag < 0 && tag != MPI_ANY_TAG){
335 retval = MPI_ERR_TAG;
337 int rank = smpi_process()->index();
338 int dst_traced = comm->group()->index(dst);
339 instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
340 extra->type = TRACING_SEND;
342 extra->dst = dst_traced;
343 extra->datatype1 = encode_datatype(datatype);
344 extra->send_size = datatype->is_basic() ? count : count * datatype->size();
345 TRACE_smpi_comm_in(rank, __FUNCTION__, extra);
346 if (not TRACE_smpi_view_internals()) {
347 TRACE_smpi_send(rank, rank, dst_traced, tag,count*datatype->size());
350 simgrid::smpi::Request::send(buf, count, datatype, dst, tag, comm);
351 retval = MPI_SUCCESS;
353 TRACE_smpi_comm_out(rank);
360 int PMPI_Ssend(void* buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm) {
365 if (comm == MPI_COMM_NULL) {
366 retval = MPI_ERR_COMM;
367 } else if (dst == MPI_PROC_NULL) {
368 retval = MPI_SUCCESS;
369 } else if (dst >= comm->group()->size() || dst <0){
370 retval = MPI_ERR_RANK;
371 } else if ((count < 0) || (buf==nullptr && count > 0)) {
372 retval = MPI_ERR_COUNT;
373 } else if (not datatype->is_valid()) {
374 retval = MPI_ERR_TYPE;
375 } else if(tag<0 && tag != MPI_ANY_TAG){
376 retval = MPI_ERR_TAG;
378 int rank = smpi_process()->index();
379 int dst_traced = comm->group()->index(dst);
380 instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
381 extra->type = TRACING_SSEND;
383 extra->dst = dst_traced;
384 extra->datatype1 = encode_datatype(datatype);
385 extra->send_size = datatype->is_basic() ? count : count * datatype->size();
387 TRACE_smpi_comm_in(rank, __FUNCTION__, extra);
388 TRACE_smpi_send(rank, rank, dst_traced, tag, count * datatype->size());
390 simgrid::smpi::Request::ssend(buf, count, datatype, dst, tag, comm);
391 retval = MPI_SUCCESS;
393 TRACE_smpi_comm_out(rank);
400 int PMPI_Sendrecv(void* sendbuf, int sendcount, MPI_Datatype sendtype, int dst, int sendtag, void* recvbuf,
401 int recvcount, MPI_Datatype recvtype, int src, int recvtag, MPI_Comm comm, MPI_Status* status)
407 if (comm == MPI_COMM_NULL) {
408 retval = MPI_ERR_COMM;
409 } else if (not sendtype->is_valid() || not recvtype->is_valid()) {
410 retval = MPI_ERR_TYPE;
411 } else if (src == MPI_PROC_NULL || dst == MPI_PROC_NULL) {
412 simgrid::smpi::Status::empty(status);
413 status->MPI_SOURCE = MPI_PROC_NULL;
414 retval = MPI_SUCCESS;
415 }else if (dst >= comm->group()->size() || dst <0 ||
416 (src!=MPI_ANY_SOURCE && (src >= comm->group()->size() || src <0))){
417 retval = MPI_ERR_RANK;
418 } else if ((sendcount < 0 || recvcount<0) ||
419 (sendbuf==nullptr && sendcount > 0) || (recvbuf==nullptr && recvcount>0)) {
420 retval = MPI_ERR_COUNT;
421 } else if((sendtag<0 && sendtag != MPI_ANY_TAG)||(recvtag<0 && recvtag != MPI_ANY_TAG)){
422 retval = MPI_ERR_TAG;
424 int rank = smpi_process()->index();
425 int dst_traced = comm->group()->index(dst);
426 int src_traced = comm->group()->index(src);
427 instr_extra_data extra = xbt_new0(s_instr_extra_data_t, 1);
428 extra->type = TRACING_SENDRECV;
429 extra->src = src_traced;
430 extra->dst = dst_traced;
431 extra->datatype1 = encode_datatype(sendtype);
432 extra->send_size = sendtype->is_basic() ? sendcount : sendcount * sendtype->size();
433 extra->datatype2 = encode_datatype(recvtype);
434 extra->recv_size = recvtype->is_basic() ? recvcount : recvcount * recvtype->size();
436 TRACE_smpi_comm_in(rank, __FUNCTION__, extra);
437 TRACE_smpi_send(rank, rank, dst_traced, sendtag, sendcount * sendtype->size());
439 simgrid::smpi::Request::sendrecv(sendbuf, sendcount, sendtype, dst, sendtag, recvbuf, recvcount, recvtype, src,
440 recvtag, comm, status);
441 retval = MPI_SUCCESS;
443 TRACE_smpi_comm_out(rank);
444 TRACE_smpi_recv(src_traced, rank, recvtag);
451 int PMPI_Sendrecv_replace(void* buf, int count, MPI_Datatype datatype, int dst, int sendtag, int src, int recvtag,
452 MPI_Comm comm, MPI_Status* status)
455 if (not datatype->is_valid()) {
457 } else if (count < 0) {
458 return MPI_ERR_COUNT;
460 int size = datatype->get_extent() * count;
461 void* recvbuf = xbt_new0(char, size);
462 retval = MPI_Sendrecv(buf, count, datatype, dst, sendtag, recvbuf, count, datatype, src, recvtag, comm, status);
463 if(retval==MPI_SUCCESS){
464 simgrid::smpi::Datatype::copy(recvbuf, count, datatype, buf, count, datatype);
472 int PMPI_Test(MPI_Request * request, int *flag, MPI_Status * status)
476 if (request == nullptr || flag == nullptr) {
477 retval = MPI_ERR_ARG;
478 } else if (*request == MPI_REQUEST_NULL) {
480 simgrid::smpi::Status::empty(status);
481 retval = MPI_SUCCESS;
483 int rank = ((*request)->comm() != MPI_COMM_NULL) ? smpi_process()->index() : -1;
485 instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
486 extra->type = TRACING_TEST;
487 TRACE_smpi_testing_in(rank, extra);
489 *flag = simgrid::smpi::Request::test(request,status);
491 TRACE_smpi_testing_out(rank);
492 retval = MPI_SUCCESS;
498 int PMPI_Testany(int count, MPI_Request requests[], int *index, int *flag, MPI_Status * status)
503 if (index == nullptr || flag == nullptr) {
504 retval = MPI_ERR_ARG;
506 *flag = simgrid::smpi::Request::testany(count, requests, index, status);
507 retval = MPI_SUCCESS;
513 int PMPI_Testall(int count, MPI_Request* requests, int* flag, MPI_Status* statuses)
518 if (flag == nullptr) {
519 retval = MPI_ERR_ARG;
521 *flag = simgrid::smpi::Request::testall(count, requests, statuses);
522 retval = MPI_SUCCESS;
528 int PMPI_Probe(int source, int tag, MPI_Comm comm, MPI_Status* status) {
532 if (status == nullptr) {
533 retval = MPI_ERR_ARG;
534 } else if (comm == MPI_COMM_NULL) {
535 retval = MPI_ERR_COMM;
536 } else if (source == MPI_PROC_NULL) {
537 simgrid::smpi::Status::empty(status);
538 status->MPI_SOURCE = MPI_PROC_NULL;
539 retval = MPI_SUCCESS;
541 simgrid::smpi::Request::probe(source, tag, comm, status);
542 retval = MPI_SUCCESS;
548 int PMPI_Iprobe(int source, int tag, MPI_Comm comm, int* flag, MPI_Status* status) {
552 if (flag == nullptr) {
553 retval = MPI_ERR_ARG;
554 } else if (comm == MPI_COMM_NULL) {
555 retval = MPI_ERR_COMM;
556 } else if (source == MPI_PROC_NULL) {
558 simgrid::smpi::Status::empty(status);
559 status->MPI_SOURCE = MPI_PROC_NULL;
560 retval = MPI_SUCCESS;
562 simgrid::smpi::Request::iprobe(source, tag, comm, flag, status);
563 retval = MPI_SUCCESS;
569 int PMPI_Wait(MPI_Request * request, MPI_Status * status)
575 simgrid::smpi::Status::empty(status);
577 if (request == nullptr) {
578 retval = MPI_ERR_ARG;
579 } else if (*request == MPI_REQUEST_NULL) {
580 retval = MPI_SUCCESS;
582 int rank = (*request)->comm() != MPI_COMM_NULL ? smpi_process()->index() : -1;
584 int src_traced = (*request)->src();
585 int dst_traced = (*request)->dst();
586 int tag_traced= (*request)->tag();
587 MPI_Comm comm = (*request)->comm();
588 int is_wait_for_receive = ((*request)->flags() & RECV);
589 instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
590 extra->type = TRACING_WAIT;
591 TRACE_smpi_comm_in(rank, __FUNCTION__, extra);
593 simgrid::smpi::Request::wait(request, status);
594 retval = MPI_SUCCESS;
596 //the src may not have been known at the beginning of the recv (MPI_ANY_SOURCE)
597 TRACE_smpi_comm_out(rank);
598 if (is_wait_for_receive) {
599 if(src_traced==MPI_ANY_SOURCE)
600 src_traced = (status != MPI_STATUS_IGNORE) ? comm->group()->rank(status->MPI_SOURCE) : src_traced;
601 TRACE_smpi_recv(src_traced, dst_traced, tag_traced);
609 int PMPI_Waitany(int count, MPI_Request requests[], int *index, MPI_Status * status)
611 if (index == nullptr)
618 //save requests information for tracing
619 struct savedvalstype {
626 savedvalstype* savedvals = xbt_new0(savedvalstype, count);
628 for (int i = 0; i < count; i++) {
629 MPI_Request req = requests[i]; //already received requests are no longer valid
631 savedvals[i]=(savedvalstype){req->src(), req->dst(), (req->flags() & RECV), req->tag(), req->comm()};
634 int rank_traced = smpi_process()->index();
635 instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
636 extra->type = TRACING_WAITANY;
637 extra->send_size=count;
638 TRACE_smpi_comm_in(rank_traced, __FUNCTION__, extra);
640 *index = simgrid::smpi::Request::waitany(count, requests, status);
642 if(*index!=MPI_UNDEFINED){
643 int src_traced = savedvals[*index].src;
644 //the src may not have been known at the beginning of the recv (MPI_ANY_SOURCE)
645 int dst_traced = savedvals[*index].dst;
646 int is_wait_for_receive = savedvals[*index].recv;
647 if (is_wait_for_receive) {
648 if(savedvals[*index].src==MPI_ANY_SOURCE)
649 src_traced = (status != MPI_STATUSES_IGNORE) ? savedvals[*index].comm->group()->rank(status->MPI_SOURCE)
650 : savedvals[*index].src;
651 TRACE_smpi_recv(src_traced, dst_traced, savedvals[*index].tag);
653 TRACE_smpi_comm_out(rank_traced);
661 int PMPI_Waitall(int count, MPI_Request requests[], MPI_Status status[])
664 //save information from requests
665 struct savedvalstype {
673 savedvalstype* savedvals=xbt_new0(savedvalstype, count);
675 for (int i = 0; i < count; i++) {
676 MPI_Request req = requests[i];
677 if(req!=MPI_REQUEST_NULL){
678 savedvals[i]=(savedvalstype){req->src(), req->dst(), (req->flags() & RECV), req->tag(), 1, req->comm()};
680 savedvals[i].valid=0;
683 int rank_traced = smpi_process()->index();
684 instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
685 extra->type = TRACING_WAITALL;
686 extra->send_size=count;
687 TRACE_smpi_comm_in(rank_traced, __FUNCTION__, extra);
689 int retval = simgrid::smpi::Request::waitall(count, requests, status);
691 for (int i = 0; i < count; i++) {
692 if(savedvals[i].valid){
693 // the src may not have been known at the beginning of the recv (MPI_ANY_SOURCE)
694 int src_traced = savedvals[i].src;
695 int dst_traced = savedvals[i].dst;
696 int is_wait_for_receive = savedvals[i].recv;
697 if (is_wait_for_receive) {
698 if(src_traced==MPI_ANY_SOURCE)
699 src_traced = (status != MPI_STATUSES_IGNORE) ? savedvals[i].comm->group()->rank(status[i].MPI_SOURCE)
701 TRACE_smpi_recv(src_traced, dst_traced,savedvals[i].tag);
705 TRACE_smpi_comm_out(rank_traced);
712 int PMPI_Waitsome(int incount, MPI_Request requests[], int *outcount, int *indices, MPI_Status status[])
717 if (outcount == nullptr) {
718 retval = MPI_ERR_ARG;
720 *outcount = simgrid::smpi::Request::waitsome(incount, requests, indices, status);
721 retval = MPI_SUCCESS;
727 int PMPI_Testsome(int incount, MPI_Request requests[], int* outcount, int* indices, MPI_Status status[])
732 if (outcount == nullptr) {
733 retval = MPI_ERR_ARG;
735 *outcount = simgrid::smpi::Request::testsome(incount, requests, indices, status);
736 retval = MPI_SUCCESS;
742 MPI_Request PMPI_Request_f2c(MPI_Fint request){
743 return static_cast<MPI_Request>(simgrid::smpi::Request::f2c(request));
746 MPI_Fint PMPI_Request_c2f(MPI_Request request) {
747 return request->c2f();