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 TRACE_smpi_comm_in(rank, __FUNCTION__,
163 new simgrid::instr::Pt2PtTIData("Irecv", comm->group()->index(src),
164 datatype->is_basic() ? count : count * datatype->size(),
165 encode_datatype(datatype)));
167 *request = simgrid::smpi::Request::irecv(buf, count, datatype, src, tag, comm);
168 retval = MPI_SUCCESS;
170 TRACE_smpi_comm_out(rank);
174 if (retval != MPI_SUCCESS && request != nullptr)
175 *request = MPI_REQUEST_NULL;
180 int PMPI_Isend(void *buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm, MPI_Request * request)
185 if (request == nullptr) {
186 retval = MPI_ERR_ARG;
187 } else if (comm == MPI_COMM_NULL) {
188 retval = MPI_ERR_COMM;
189 } else if (dst == MPI_PROC_NULL) {
190 *request = MPI_REQUEST_NULL;
191 retval = MPI_SUCCESS;
192 } else if (dst >= comm->group()->size() || dst <0){
193 retval = MPI_ERR_RANK;
194 } else if ((count < 0) || (buf==nullptr && count > 0)) {
195 retval = MPI_ERR_COUNT;
196 } else if (not datatype->is_valid()) {
197 retval = MPI_ERR_TYPE;
198 } else if(tag<0 && tag != MPI_ANY_TAG){
199 retval = MPI_ERR_TAG;
201 int rank = smpi_process()->index();
202 int trace_dst = comm->group()->index(dst);
203 TRACE_smpi_comm_in(rank, __FUNCTION__,
204 new simgrid::instr::Pt2PtTIData("Isend", trace_dst,
205 datatype->is_basic() ? count : count * datatype->size(),
206 encode_datatype(datatype)));
208 TRACE_smpi_send(rank, rank, trace_dst, tag, count * datatype->size());
210 *request = simgrid::smpi::Request::isend(buf, count, datatype, dst, tag, comm);
211 retval = MPI_SUCCESS;
213 TRACE_smpi_comm_out(rank);
217 if (retval != MPI_SUCCESS && request!=nullptr)
218 *request = MPI_REQUEST_NULL;
222 int PMPI_Issend(void* buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm, MPI_Request* request)
227 if (request == nullptr) {
228 retval = MPI_ERR_ARG;
229 } else if (comm == MPI_COMM_NULL) {
230 retval = MPI_ERR_COMM;
231 } else if (dst == MPI_PROC_NULL) {
232 *request = MPI_REQUEST_NULL;
233 retval = MPI_SUCCESS;
234 } else if (dst >= comm->group()->size() || dst <0){
235 retval = MPI_ERR_RANK;
236 } else if ((count < 0)|| (buf==nullptr && count > 0)) {
237 retval = MPI_ERR_COUNT;
238 } else if (not datatype->is_valid()) {
239 retval = MPI_ERR_TYPE;
240 } else if(tag<0 && tag != MPI_ANY_TAG){
241 retval = MPI_ERR_TAG;
243 int rank = smpi_process()->index();
244 int trace_dst = comm->group()->index(dst);
245 TRACE_smpi_comm_in(rank, __FUNCTION__,
246 new simgrid::instr::Pt2PtTIData("ISsend", trace_dst,
247 datatype->is_basic() ? count : count * datatype->size(),
248 encode_datatype(datatype)));
249 TRACE_smpi_send(rank, rank, trace_dst, tag, count * datatype->size());
251 *request = simgrid::smpi::Request::issend(buf, count, datatype, dst, tag, comm);
252 retval = MPI_SUCCESS;
254 TRACE_smpi_comm_out(rank);
258 if (retval != MPI_SUCCESS && request!=nullptr)
259 *request = MPI_REQUEST_NULL;
263 int PMPI_Recv(void *buf, int count, MPI_Datatype datatype, int src, int tag, MPI_Comm comm, MPI_Status * status)
268 if (comm == MPI_COMM_NULL) {
269 retval = MPI_ERR_COMM;
270 } else if (src == MPI_PROC_NULL) {
271 simgrid::smpi::Status::empty(status);
272 status->MPI_SOURCE = MPI_PROC_NULL;
273 retval = MPI_SUCCESS;
274 } else if (src!=MPI_ANY_SOURCE && (src >= comm->group()->size() || src <0)){
275 retval = MPI_ERR_RANK;
276 } else if ((count < 0) || (buf==nullptr && count > 0)) {
277 retval = MPI_ERR_COUNT;
278 } else if (not datatype->is_valid()) {
279 retval = MPI_ERR_TYPE;
280 } else if(tag<0 && tag != MPI_ANY_TAG){
281 retval = MPI_ERR_TAG;
283 int rank = smpi_process()->index();
284 int src_traced = comm->group()->index(src);
285 TRACE_smpi_comm_in(rank, __FUNCTION__,
286 new simgrid::instr::Pt2PtTIData("recv", src_traced,
287 datatype->is_basic() ? count : count * datatype->size(),
288 encode_datatype(datatype)));
290 simgrid::smpi::Request::recv(buf, count, datatype, src, tag, comm, status);
291 retval = MPI_SUCCESS;
293 // the src may not have been known at the beginning of the recv (MPI_ANY_SOURCE)
294 if (status != MPI_STATUS_IGNORE) {
295 src_traced = comm->group()->index(status->MPI_SOURCE);
296 if (not TRACE_smpi_view_internals()) {
297 TRACE_smpi_recv(src_traced, rank, tag);
300 TRACE_smpi_comm_out(rank);
307 int PMPI_Send(void *buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm)
313 if (comm == MPI_COMM_NULL) {
314 retval = MPI_ERR_COMM;
315 } else if (dst == MPI_PROC_NULL) {
316 retval = MPI_SUCCESS;
317 } else if (dst >= comm->group()->size() || dst <0){
318 retval = MPI_ERR_RANK;
319 } else if ((count < 0) || (buf == nullptr && count > 0)) {
320 retval = MPI_ERR_COUNT;
321 } else if (not datatype->is_valid()) {
322 retval = MPI_ERR_TYPE;
323 } else if(tag < 0 && tag != MPI_ANY_TAG){
324 retval = MPI_ERR_TAG;
326 int rank = smpi_process()->index();
327 int dst_traced = comm->group()->index(dst);
328 TRACE_smpi_comm_in(rank, __FUNCTION__,
329 new simgrid::instr::Pt2PtTIData("send", dst_traced,
330 datatype->is_basic() ? count : count * datatype->size(),
331 encode_datatype(datatype)));
332 if (not TRACE_smpi_view_internals()) {
333 TRACE_smpi_send(rank, rank, dst_traced, tag,count*datatype->size());
336 simgrid::smpi::Request::send(buf, count, datatype, dst, tag, comm);
337 retval = MPI_SUCCESS;
339 TRACE_smpi_comm_out(rank);
346 int PMPI_Ssend(void* buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm) {
351 if (comm == MPI_COMM_NULL) {
352 retval = MPI_ERR_COMM;
353 } else if (dst == MPI_PROC_NULL) {
354 retval = MPI_SUCCESS;
355 } else if (dst >= comm->group()->size() || dst <0){
356 retval = MPI_ERR_RANK;
357 } else if ((count < 0) || (buf==nullptr && count > 0)) {
358 retval = MPI_ERR_COUNT;
359 } else if (not datatype->is_valid()) {
360 retval = MPI_ERR_TYPE;
361 } else if(tag<0 && tag != MPI_ANY_TAG){
362 retval = MPI_ERR_TAG;
364 int rank = smpi_process()->index();
365 int dst_traced = comm->group()->index(dst);
366 TRACE_smpi_comm_in(rank, __FUNCTION__,
367 new simgrid::instr::Pt2PtTIData("Ssend", dst_traced,
368 datatype->is_basic() ? count : count * datatype->size(),
369 encode_datatype(datatype)));
370 TRACE_smpi_send(rank, rank, dst_traced, tag, count * datatype->size());
372 simgrid::smpi::Request::ssend(buf, count, datatype, dst, tag, comm);
373 retval = MPI_SUCCESS;
375 TRACE_smpi_comm_out(rank);
382 int PMPI_Sendrecv(void* sendbuf, int sendcount, MPI_Datatype sendtype, int dst, int sendtag, void* recvbuf,
383 int recvcount, MPI_Datatype recvtype, int src, int recvtag, MPI_Comm comm, MPI_Status* status)
389 if (comm == MPI_COMM_NULL) {
390 retval = MPI_ERR_COMM;
391 } else if (not sendtype->is_valid() || not recvtype->is_valid()) {
392 retval = MPI_ERR_TYPE;
393 } else if (src == MPI_PROC_NULL || dst == MPI_PROC_NULL) {
394 simgrid::smpi::Status::empty(status);
395 status->MPI_SOURCE = MPI_PROC_NULL;
396 retval = MPI_SUCCESS;
397 }else if (dst >= comm->group()->size() || dst <0 ||
398 (src!=MPI_ANY_SOURCE && (src >= comm->group()->size() || src <0))){
399 retval = MPI_ERR_RANK;
400 } else if ((sendcount < 0 || recvcount<0) ||
401 (sendbuf==nullptr && sendcount > 0) || (recvbuf==nullptr && recvcount>0)) {
402 retval = MPI_ERR_COUNT;
403 } else if((sendtag<0 && sendtag != MPI_ANY_TAG)||(recvtag<0 && recvtag != MPI_ANY_TAG)){
404 retval = MPI_ERR_TAG;
406 int rank = smpi_process()->index();
407 int dst_traced = comm->group()->index(dst);
408 int src_traced = comm->group()->index(src);
410 // FIXME: Hack the way to trace this one
411 std::vector<int>* dst_hack = new std::vector<int>;
412 std::vector<int>* src_hack = new std::vector<int>;
413 dst_hack->push_back(dst_traced);
414 src_hack->push_back(src_traced);
415 TRACE_smpi_comm_in(rank, __FUNCTION__,
416 new simgrid::instr::VarCollTIData(
417 "sendRecv", -1, sendtype->is_basic() ? sendcount : sendcount * sendtype->size(), dst_hack,
418 recvtype->is_basic() ? recvcount : recvcount * recvtype->size(), src_hack,
419 encode_datatype(sendtype), encode_datatype(recvtype)));
421 TRACE_smpi_send(rank, rank, dst_traced, sendtag, sendcount * sendtype->size());
423 simgrid::smpi::Request::sendrecv(sendbuf, sendcount, sendtype, dst, sendtag, recvbuf, recvcount, recvtype, src,
424 recvtag, comm, status);
425 retval = MPI_SUCCESS;
427 TRACE_smpi_recv(src_traced, rank, recvtag);
428 TRACE_smpi_comm_out(rank);
435 int PMPI_Sendrecv_replace(void* buf, int count, MPI_Datatype datatype, int dst, int sendtag, int src, int recvtag,
436 MPI_Comm comm, MPI_Status* status)
439 if (not datatype->is_valid()) {
441 } else if (count < 0) {
442 return MPI_ERR_COUNT;
444 int size = datatype->get_extent() * count;
445 void* recvbuf = xbt_new0(char, size);
446 retval = MPI_Sendrecv(buf, count, datatype, dst, sendtag, recvbuf, count, datatype, src, recvtag, comm, status);
447 if(retval==MPI_SUCCESS){
448 simgrid::smpi::Datatype::copy(recvbuf, count, datatype, buf, count, datatype);
456 int PMPI_Test(MPI_Request * request, int *flag, MPI_Status * status)
460 if (request == nullptr || flag == nullptr) {
461 retval = MPI_ERR_ARG;
462 } else if (*request == MPI_REQUEST_NULL) {
464 simgrid::smpi::Status::empty(status);
465 retval = MPI_SUCCESS;
467 int rank = ((*request)->comm() != MPI_COMM_NULL) ? smpi_process()->index() : -1;
469 TRACE_smpi_testing_in(rank);
471 *flag = simgrid::smpi::Request::test(request,status);
473 TRACE_smpi_testing_out(rank);
474 retval = MPI_SUCCESS;
480 int PMPI_Testany(int count, MPI_Request requests[], int *index, int *flag, MPI_Status * status)
485 if (index == nullptr || flag == nullptr) {
486 retval = MPI_ERR_ARG;
488 *flag = simgrid::smpi::Request::testany(count, requests, index, status);
489 retval = MPI_SUCCESS;
495 int PMPI_Testall(int count, MPI_Request* requests, int* flag, MPI_Status* statuses)
500 if (flag == nullptr) {
501 retval = MPI_ERR_ARG;
503 *flag = simgrid::smpi::Request::testall(count, requests, statuses);
504 retval = MPI_SUCCESS;
510 int PMPI_Probe(int source, int tag, MPI_Comm comm, MPI_Status* status) {
514 if (status == nullptr) {
515 retval = MPI_ERR_ARG;
516 } else if (comm == MPI_COMM_NULL) {
517 retval = MPI_ERR_COMM;
518 } else if (source == MPI_PROC_NULL) {
519 simgrid::smpi::Status::empty(status);
520 status->MPI_SOURCE = MPI_PROC_NULL;
521 retval = MPI_SUCCESS;
523 simgrid::smpi::Request::probe(source, tag, comm, status);
524 retval = MPI_SUCCESS;
530 int PMPI_Iprobe(int source, int tag, MPI_Comm comm, int* flag, MPI_Status* status) {
534 if (flag == nullptr) {
535 retval = MPI_ERR_ARG;
536 } else if (comm == MPI_COMM_NULL) {
537 retval = MPI_ERR_COMM;
538 } else if (source == MPI_PROC_NULL) {
540 simgrid::smpi::Status::empty(status);
541 status->MPI_SOURCE = MPI_PROC_NULL;
542 retval = MPI_SUCCESS;
544 simgrid::smpi::Request::iprobe(source, tag, comm, flag, status);
545 retval = MPI_SUCCESS;
551 int PMPI_Wait(MPI_Request * request, MPI_Status * status)
557 simgrid::smpi::Status::empty(status);
559 if (request == nullptr) {
560 retval = MPI_ERR_ARG;
561 } else if (*request == MPI_REQUEST_NULL) {
562 retval = MPI_SUCCESS;
564 int rank = (*request)->comm() != MPI_COMM_NULL ? smpi_process()->index() : -1;
566 int src_traced = (*request)->src();
567 int dst_traced = (*request)->dst();
568 int tag_traced= (*request)->tag();
569 MPI_Comm comm = (*request)->comm();
570 int is_wait_for_receive = ((*request)->flags() & RECV);
571 TRACE_smpi_comm_in(rank, __FUNCTION__, new simgrid::instr::NoOpTIData("wait"));
573 simgrid::smpi::Request::wait(request, status);
574 retval = MPI_SUCCESS;
576 //the src may not have been known at the beginning of the recv (MPI_ANY_SOURCE)
577 TRACE_smpi_comm_out(rank);
578 if (is_wait_for_receive) {
579 if(src_traced==MPI_ANY_SOURCE)
580 src_traced = (status != MPI_STATUS_IGNORE) ? comm->group()->rank(status->MPI_SOURCE) : src_traced;
581 TRACE_smpi_recv(src_traced, dst_traced, tag_traced);
589 int PMPI_Waitany(int count, MPI_Request requests[], int *index, MPI_Status * status)
591 if (index == nullptr)
598 //save requests information for tracing
599 struct savedvalstype {
606 savedvalstype* savedvals = xbt_new0(savedvalstype, count);
608 for (int i = 0; i < count; i++) {
609 MPI_Request req = requests[i]; //already received requests are no longer valid
611 savedvals[i]=(savedvalstype){req->src(), req->dst(), (req->flags() & RECV), req->tag(), req->comm()};
614 int rank_traced = smpi_process()->index();
615 TRACE_smpi_comm_in(rank_traced, __FUNCTION__, new simgrid::instr::CpuTIData("waitAny", static_cast<double>(count)));
617 *index = simgrid::smpi::Request::waitany(count, requests, status);
619 if(*index!=MPI_UNDEFINED){
620 int src_traced = savedvals[*index].src;
621 //the src may not have been known at the beginning of the recv (MPI_ANY_SOURCE)
622 int dst_traced = savedvals[*index].dst;
623 int is_wait_for_receive = savedvals[*index].recv;
624 if (is_wait_for_receive) {
625 if(savedvals[*index].src==MPI_ANY_SOURCE)
626 src_traced = (status != MPI_STATUSES_IGNORE) ? savedvals[*index].comm->group()->rank(status->MPI_SOURCE)
627 : savedvals[*index].src;
628 TRACE_smpi_recv(src_traced, dst_traced, savedvals[*index].tag);
630 TRACE_smpi_comm_out(rank_traced);
638 int PMPI_Waitall(int count, MPI_Request requests[], MPI_Status status[])
641 //save information from requests
642 struct savedvalstype {
650 savedvalstype* savedvals=xbt_new0(savedvalstype, count);
652 for (int i = 0; i < count; i++) {
653 MPI_Request req = requests[i];
654 if(req!=MPI_REQUEST_NULL){
655 savedvals[i]=(savedvalstype){req->src(), req->dst(), (req->flags() & RECV), req->tag(), 1, req->comm()};
657 savedvals[i].valid=0;
660 int rank_traced = smpi_process()->index();
661 TRACE_smpi_comm_in(rank_traced, __FUNCTION__, new simgrid::instr::CpuTIData("waitAll", static_cast<double>(count)));
663 int retval = simgrid::smpi::Request::waitall(count, requests, status);
665 for (int i = 0; i < count; i++) {
666 if(savedvals[i].valid){
667 // the src may not have been known at the beginning of the recv (MPI_ANY_SOURCE)
668 int src_traced = savedvals[i].src;
669 int dst_traced = savedvals[i].dst;
670 int is_wait_for_receive = savedvals[i].recv;
671 if (is_wait_for_receive) {
672 if(src_traced==MPI_ANY_SOURCE)
673 src_traced = (status != MPI_STATUSES_IGNORE) ? savedvals[i].comm->group()->rank(status[i].MPI_SOURCE)
675 TRACE_smpi_recv(src_traced, dst_traced,savedvals[i].tag);
679 TRACE_smpi_comm_out(rank_traced);
686 int PMPI_Waitsome(int incount, MPI_Request requests[], int *outcount, int *indices, MPI_Status status[])
691 if (outcount == nullptr) {
692 retval = MPI_ERR_ARG;
694 *outcount = simgrid::smpi::Request::waitsome(incount, requests, indices, status);
695 retval = MPI_SUCCESS;
701 int PMPI_Testsome(int incount, MPI_Request requests[], int* outcount, int* indices, MPI_Status status[])
706 if (outcount == nullptr) {
707 retval = MPI_ERR_ARG;
709 *outcount = simgrid::smpi::Request::testsome(incount, requests, indices, status);
710 retval = MPI_SUCCESS;
716 MPI_Request PMPI_Request_f2c(MPI_Fint request){
717 return static_cast<MPI_Request>(simgrid::smpi::Request::f2c(request));
720 MPI_Fint PMPI_Request_c2f(MPI_Request request) {
721 return request->c2f();