1 /* Copyright (c) 2007-2018. 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 static int getPid(MPI_Comm, int);
15 static int getPid(MPI_Comm comm, int id)
17 simgrid::s4u::ActorPtr actor = comm->group()->actor(id);
18 return (actor == nullptr) ? MPI_UNDEFINED : actor->getPid();
21 /* PMPI User level calls */
22 extern "C" { // Obviously, the C MPI interface should use the C linkage
24 int PMPI_Send_init(void *buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm, MPI_Request * request)
29 if (request == nullptr) {
31 } else if (comm == MPI_COMM_NULL) {
32 retval = MPI_ERR_COMM;
33 } else if (not datatype->is_valid()) {
34 retval = MPI_ERR_TYPE;
35 } else if (dst == MPI_PROC_NULL) {
38 *request = simgrid::smpi::Request::send_init(buf, count, datatype, dst, tag, comm);
42 if (retval != MPI_SUCCESS && request != nullptr)
43 *request = MPI_REQUEST_NULL;
47 int PMPI_Recv_init(void *buf, int count, MPI_Datatype datatype, int src, int tag, MPI_Comm comm, MPI_Request * request)
52 if (request == nullptr) {
54 } else if (comm == MPI_COMM_NULL) {
55 retval = MPI_ERR_COMM;
56 } else if (not datatype->is_valid()) {
57 retval = MPI_ERR_TYPE;
58 } else if (src == MPI_PROC_NULL) {
61 *request = simgrid::smpi::Request::recv_init(buf, count, datatype, src, tag, comm);
65 if (retval != MPI_SUCCESS && request != nullptr)
66 *request = MPI_REQUEST_NULL;
70 int PMPI_Ssend_init(void* buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm, MPI_Request* request)
75 if (request == nullptr) {
77 } else if (comm == MPI_COMM_NULL) {
78 retval = MPI_ERR_COMM;
79 } else if (not datatype->is_valid()) {
80 retval = MPI_ERR_TYPE;
81 } else if (dst == MPI_PROC_NULL) {
84 *request = simgrid::smpi::Request::ssend_init(buf, count, datatype, dst, tag, comm);
88 if (retval != MPI_SUCCESS && request != nullptr)
89 *request = MPI_REQUEST_NULL;
93 int PMPI_Start(MPI_Request * request)
98 if (request == nullptr || *request == MPI_REQUEST_NULL) {
99 retval = MPI_ERR_REQUEST;
102 retval = MPI_SUCCESS;
108 int PMPI_Startall(int count, MPI_Request * requests)
112 if (requests == nullptr) {
113 retval = MPI_ERR_ARG;
115 retval = MPI_SUCCESS;
116 for (int i = 0; i < count; i++) {
117 if(requests[i] == MPI_REQUEST_NULL) {
118 retval = MPI_ERR_REQUEST;
121 if(retval != MPI_ERR_REQUEST) {
122 simgrid::smpi::Request::startall(count, requests);
129 int PMPI_Request_free(MPI_Request * request)
134 if (*request == MPI_REQUEST_NULL) {
135 retval = MPI_ERR_ARG;
137 simgrid::smpi::Request::unref(request);
138 retval = MPI_SUCCESS;
144 int PMPI_Irecv(void *buf, int count, MPI_Datatype datatype, int src, int tag, MPI_Comm comm, MPI_Request * request)
150 if (request == nullptr) {
151 retval = MPI_ERR_ARG;
152 } else if (comm == MPI_COMM_NULL) {
153 retval = MPI_ERR_COMM;
154 } else if (src == MPI_PROC_NULL) {
155 *request = MPI_REQUEST_NULL;
156 retval = MPI_SUCCESS;
157 } else if (src!=MPI_ANY_SOURCE && (src >= comm->group()->size() || src <0)){
158 retval = MPI_ERR_RANK;
159 } else if ((count < 0) || (buf==nullptr && count > 0)) {
160 retval = MPI_ERR_COUNT;
161 } else if (not datatype->is_valid()) {
162 retval = MPI_ERR_TYPE;
163 } else if(tag<0 && tag != MPI_ANY_TAG){
164 retval = MPI_ERR_TAG;
167 int my_proc_id = simgrid::s4u::Actor::self()->getPid();
169 TRACE_smpi_comm_in(my_proc_id, __FUNCTION__,
170 new simgrid::instr::Pt2PtTIData("Irecv", src,
171 datatype->is_replayable() ? count : count * datatype->size(),
172 encode_datatype(datatype)));
174 *request = simgrid::smpi::Request::irecv(buf, count, datatype, src, tag, comm);
175 retval = MPI_SUCCESS;
177 TRACE_smpi_comm_out(my_proc_id);
181 if (retval != MPI_SUCCESS && request != nullptr)
182 *request = MPI_REQUEST_NULL;
187 int PMPI_Isend(void *buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm, MPI_Request * request)
192 if (request == nullptr) {
193 retval = MPI_ERR_ARG;
194 } else if (comm == MPI_COMM_NULL) {
195 retval = MPI_ERR_COMM;
196 } else if (dst == MPI_PROC_NULL) {
197 *request = MPI_REQUEST_NULL;
198 retval = MPI_SUCCESS;
199 } else if (dst >= comm->group()->size() || dst <0){
200 retval = MPI_ERR_RANK;
201 } else if ((count < 0) || (buf==nullptr && count > 0)) {
202 retval = MPI_ERR_COUNT;
203 } else if (not datatype->is_valid()) {
204 retval = MPI_ERR_TYPE;
205 } else if(tag<0 && tag != MPI_ANY_TAG){
206 retval = MPI_ERR_TAG;
208 int my_proc_id = simgrid::s4u::Actor::self()->getPid();
209 int trace_dst = getPid(comm, dst);
210 TRACE_smpi_comm_in(my_proc_id, __FUNCTION__,
211 new simgrid::instr::Pt2PtTIData("Isend", dst,
212 datatype->is_replayable() ? count : count * datatype->size(),
213 encode_datatype(datatype)));
215 TRACE_smpi_send(my_proc_id, my_proc_id, trace_dst, tag, count * datatype->size());
217 *request = simgrid::smpi::Request::isend(buf, count, datatype, dst, tag, comm);
218 retval = MPI_SUCCESS;
220 TRACE_smpi_comm_out(my_proc_id);
224 if (retval != MPI_SUCCESS && request!=nullptr)
225 *request = MPI_REQUEST_NULL;
229 int PMPI_Issend(void* buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm, MPI_Request* request)
234 if (request == nullptr) {
235 retval = MPI_ERR_ARG;
236 } else if (comm == MPI_COMM_NULL) {
237 retval = MPI_ERR_COMM;
238 } else if (dst == MPI_PROC_NULL) {
239 *request = MPI_REQUEST_NULL;
240 retval = MPI_SUCCESS;
241 } else if (dst >= comm->group()->size() || dst <0){
242 retval = MPI_ERR_RANK;
243 } else if ((count < 0)|| (buf==nullptr && count > 0)) {
244 retval = MPI_ERR_COUNT;
245 } else if (not datatype->is_valid()) {
246 retval = MPI_ERR_TYPE;
247 } else if(tag<0 && tag != MPI_ANY_TAG){
248 retval = MPI_ERR_TAG;
250 int my_proc_id = simgrid::s4u::Actor::self()->getPid();
251 int trace_dst = getPid(comm, dst);
252 TRACE_smpi_comm_in(my_proc_id, __FUNCTION__,
253 new simgrid::instr::Pt2PtTIData("ISsend", dst,
254 datatype->is_replayable() ? count : count * datatype->size(),
255 encode_datatype(datatype)));
256 TRACE_smpi_send(my_proc_id, my_proc_id, trace_dst, tag, count * datatype->size());
258 *request = simgrid::smpi::Request::issend(buf, count, datatype, dst, tag, comm);
259 retval = MPI_SUCCESS;
261 TRACE_smpi_comm_out(my_proc_id);
265 if (retval != MPI_SUCCESS && request!=nullptr)
266 *request = MPI_REQUEST_NULL;
270 int PMPI_Recv(void *buf, int count, MPI_Datatype datatype, int src, int tag, MPI_Comm comm, MPI_Status * status)
275 if (comm == MPI_COMM_NULL) {
276 retval = MPI_ERR_COMM;
277 } else if (src == MPI_PROC_NULL) {
278 simgrid::smpi::Status::empty(status);
279 status->MPI_SOURCE = MPI_PROC_NULL;
280 retval = MPI_SUCCESS;
281 } else if (src!=MPI_ANY_SOURCE && (src >= comm->group()->size() || src <0)){
282 retval = MPI_ERR_RANK;
283 } else if ((count < 0) || (buf==nullptr && count > 0)) {
284 retval = MPI_ERR_COUNT;
285 } else if (not datatype->is_valid()) {
286 retval = MPI_ERR_TYPE;
287 } else if(tag<0 && tag != MPI_ANY_TAG){
288 retval = MPI_ERR_TAG;
290 int my_proc_id = simgrid::s4u::Actor::self()->getPid();
291 TRACE_smpi_comm_in(my_proc_id, __FUNCTION__,
292 new simgrid::instr::Pt2PtTIData("recv", src,
293 datatype->is_replayable() ? count : count * datatype->size(),
294 encode_datatype(datatype)));
296 simgrid::smpi::Request::recv(buf, count, datatype, src, tag, comm, status);
297 retval = MPI_SUCCESS;
299 // the src may not have been known at the beginning of the recv (MPI_ANY_SOURCE)
300 if (status != MPI_STATUS_IGNORE) {
301 int src_traced = getPid(comm, status->MPI_SOURCE);
302 if (not TRACE_smpi_view_internals()) {
303 TRACE_smpi_recv(src_traced, my_proc_id, tag);
306 TRACE_smpi_comm_out(my_proc_id);
313 int PMPI_Send(void *buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm)
319 if (comm == MPI_COMM_NULL) {
320 retval = MPI_ERR_COMM;
321 } else if (dst == MPI_PROC_NULL) {
322 retval = MPI_SUCCESS;
323 } else if (dst >= comm->group()->size() || dst <0){
324 retval = MPI_ERR_RANK;
325 } else if ((count < 0) || (buf == nullptr && count > 0)) {
326 retval = MPI_ERR_COUNT;
327 } else if (not datatype->is_valid()) {
328 retval = MPI_ERR_TYPE;
329 } else if(tag < 0 && tag != MPI_ANY_TAG){
330 retval = MPI_ERR_TAG;
332 int my_proc_id = simgrid::s4u::Actor::self()->getPid();
333 int dst_traced = getPid(comm, dst);
334 TRACE_smpi_comm_in(my_proc_id, __FUNCTION__,
335 new simgrid::instr::Pt2PtTIData("send", dst,
336 datatype->is_replayable() ? count : count * datatype->size(),
337 encode_datatype(datatype)));
338 if (not TRACE_smpi_view_internals()) {
339 TRACE_smpi_send(my_proc_id, my_proc_id, dst_traced, tag, count * datatype->size());
342 simgrid::smpi::Request::send(buf, count, datatype, dst, tag, comm);
343 retval = MPI_SUCCESS;
345 TRACE_smpi_comm_out(my_proc_id);
352 int PMPI_Ssend(void* buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm) {
357 if (comm == MPI_COMM_NULL) {
358 retval = MPI_ERR_COMM;
359 } else if (dst == MPI_PROC_NULL) {
360 retval = MPI_SUCCESS;
361 } else if (dst >= comm->group()->size() || dst <0){
362 retval = MPI_ERR_RANK;
363 } else if ((count < 0) || (buf==nullptr && count > 0)) {
364 retval = MPI_ERR_COUNT;
365 } else if (not datatype->is_valid()) {
366 retval = MPI_ERR_TYPE;
367 } else if(tag<0 && tag != MPI_ANY_TAG){
368 retval = MPI_ERR_TAG;
370 int my_proc_id = simgrid::s4u::Actor::self()->getPid();
371 int dst_traced = getPid(comm, dst);
372 TRACE_smpi_comm_in(my_proc_id, __FUNCTION__,
373 new simgrid::instr::Pt2PtTIData("Ssend", dst,
374 datatype->is_replayable() ? count : count * datatype->size(),
375 encode_datatype(datatype)));
376 TRACE_smpi_send(my_proc_id, my_proc_id, dst_traced, tag, count * datatype->size());
378 simgrid::smpi::Request::ssend(buf, count, datatype, dst, tag, comm);
379 retval = MPI_SUCCESS;
381 TRACE_smpi_comm_out(my_proc_id);
388 int PMPI_Sendrecv(void* sendbuf, int sendcount, MPI_Datatype sendtype, int dst, int sendtag, void* recvbuf,
389 int recvcount, MPI_Datatype recvtype, int src, int recvtag, MPI_Comm comm, MPI_Status* status)
395 if (comm == MPI_COMM_NULL) {
396 retval = MPI_ERR_COMM;
397 } else if (not sendtype->is_valid() || not recvtype->is_valid()) {
398 retval = MPI_ERR_TYPE;
399 } else if (src == MPI_PROC_NULL || dst == MPI_PROC_NULL) {
400 simgrid::smpi::Status::empty(status);
401 status->MPI_SOURCE = MPI_PROC_NULL;
402 retval = MPI_SUCCESS;
403 }else if (dst >= comm->group()->size() || dst <0 ||
404 (src!=MPI_ANY_SOURCE && (src >= comm->group()->size() || src <0))){
405 retval = MPI_ERR_RANK;
406 } else if ((sendcount < 0 || recvcount<0) ||
407 (sendbuf==nullptr && sendcount > 0) || (recvbuf==nullptr && recvcount>0)) {
408 retval = MPI_ERR_COUNT;
409 } else if((sendtag<0 && sendtag != MPI_ANY_TAG)||(recvtag<0 && recvtag != MPI_ANY_TAG)){
410 retval = MPI_ERR_TAG;
412 int my_proc_id = simgrid::s4u::Actor::self()->getPid();
413 int dst_traced = getPid(comm, dst);
414 int src_traced = getPid(comm, src);
416 // FIXME: Hack the way to trace this one
417 std::vector<int>* dst_hack = new std::vector<int>;
418 std::vector<int>* src_hack = new std::vector<int>;
419 dst_hack->push_back(dst_traced);
420 src_hack->push_back(src_traced);
421 TRACE_smpi_comm_in(my_proc_id, __FUNCTION__,
422 new simgrid::instr::VarCollTIData(
423 "sendRecv", -1, sendtype->is_replayable() ? sendcount : sendcount * sendtype->size(),
424 dst_hack, recvtype->is_replayable() ? recvcount : recvcount * recvtype->size(), src_hack,
425 encode_datatype(sendtype), encode_datatype(recvtype)));
427 TRACE_smpi_send(my_proc_id, my_proc_id, dst_traced, sendtag, sendcount * sendtype->size());
429 simgrid::smpi::Request::sendrecv(sendbuf, sendcount, sendtype, dst, sendtag, recvbuf, recvcount, recvtype, src,
430 recvtag, comm, status);
431 retval = MPI_SUCCESS;
433 TRACE_smpi_recv(src_traced, my_proc_id, recvtag);
434 TRACE_smpi_comm_out(my_proc_id);
441 int PMPI_Sendrecv_replace(void* buf, int count, MPI_Datatype datatype, int dst, int sendtag, int src, int recvtag,
442 MPI_Comm comm, MPI_Status* status)
445 if (not datatype->is_valid()) {
447 } else if (count < 0) {
448 return MPI_ERR_COUNT;
450 int size = datatype->get_extent() * count;
451 void* recvbuf = xbt_new0(char, size);
452 retval = MPI_Sendrecv(buf, count, datatype, dst, sendtag, recvbuf, count, datatype, src, recvtag, comm, status);
453 if(retval==MPI_SUCCESS){
454 simgrid::smpi::Datatype::copy(recvbuf, count, datatype, buf, count, datatype);
462 int PMPI_Test(MPI_Request * request, int *flag, MPI_Status * status)
466 if (request == nullptr || flag == nullptr) {
467 retval = MPI_ERR_ARG;
468 } else if (*request == MPI_REQUEST_NULL) {
470 simgrid::smpi::Status::empty(status);
471 retval = MPI_SUCCESS;
473 int my_proc_id = ((*request)->comm() != MPI_COMM_NULL) ? simgrid::s4u::Actor::self()->getPid() : -1;
475 TRACE_smpi_testing_in(my_proc_id);
477 *flag = simgrid::smpi::Request::test(request,status);
479 TRACE_smpi_testing_out(my_proc_id);
480 retval = MPI_SUCCESS;
486 int PMPI_Testany(int count, MPI_Request requests[], int *index, int *flag, MPI_Status * status)
491 if (index == nullptr || flag == nullptr) {
492 retval = MPI_ERR_ARG;
494 *flag = simgrid::smpi::Request::testany(count, requests, index, status);
495 retval = MPI_SUCCESS;
501 int PMPI_Testall(int count, MPI_Request* requests, int* flag, MPI_Status* statuses)
506 if (flag == nullptr) {
507 retval = MPI_ERR_ARG;
509 *flag = simgrid::smpi::Request::testall(count, requests, statuses);
510 retval = MPI_SUCCESS;
516 int PMPI_Probe(int source, int tag, MPI_Comm comm, MPI_Status* status) {
520 if (status == nullptr) {
521 retval = MPI_ERR_ARG;
522 } else if (comm == MPI_COMM_NULL) {
523 retval = MPI_ERR_COMM;
524 } else if (source == MPI_PROC_NULL) {
525 simgrid::smpi::Status::empty(status);
526 status->MPI_SOURCE = MPI_PROC_NULL;
527 retval = MPI_SUCCESS;
529 simgrid::smpi::Request::probe(source, tag, comm, status);
530 retval = MPI_SUCCESS;
536 int PMPI_Iprobe(int source, int tag, MPI_Comm comm, int* flag, MPI_Status* status) {
540 if (flag == nullptr) {
541 retval = MPI_ERR_ARG;
542 } else if (comm == MPI_COMM_NULL) {
543 retval = MPI_ERR_COMM;
544 } else if (source == MPI_PROC_NULL) {
546 simgrid::smpi::Status::empty(status);
547 status->MPI_SOURCE = MPI_PROC_NULL;
548 retval = MPI_SUCCESS;
550 simgrid::smpi::Request::iprobe(source, tag, comm, flag, status);
551 retval = MPI_SUCCESS;
557 // TODO: cheinrich: Move declaration to other file? Rename this function - it's used for PMPI_Wait*?
558 static void trace_smpi_recv_helper(MPI_Request* request, MPI_Status* status);
559 static void trace_smpi_recv_helper(MPI_Request* request, MPI_Status* status)
561 MPI_Request req = *request;
562 if (req != MPI_REQUEST_NULL) { // Received requests become null
563 int src_traced = req->src();
564 // the src may not have been known at the beginning of the recv (MPI_ANY_SOURCE)
565 int dst_traced = req->dst();
566 if (req->flags() & RECV) { // Is this request a wait for RECV?
567 if (src_traced == MPI_ANY_SOURCE)
568 src_traced = (status != MPI_STATUSES_IGNORE) ? req->comm()->group()->rank(status->MPI_SOURCE) : req->src();
569 TRACE_smpi_recv(src_traced, dst_traced, req->tag());
574 int PMPI_Wait(MPI_Request * request, MPI_Status * status)
580 simgrid::smpi::Status::empty(status);
582 if (request == nullptr) {
583 retval = MPI_ERR_ARG;
584 } else if (*request == MPI_REQUEST_NULL) {
585 retval = MPI_SUCCESS;
587 int my_proc_id = (*request)->comm() != MPI_COMM_NULL
588 ? simgrid::s4u::Actor::self()->getPid()
589 : -1; // TODO: cheinrich: Check if this correct or if it should be MPI_UNDEFINED
591 TRACE_smpi_comm_in(my_proc_id, __FUNCTION__, new simgrid::instr::NoOpTIData("wait"));
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(my_proc_id);
598 trace_smpi_recv_helper(request, status);
605 int PMPI_Waitany(int count, MPI_Request requests[], int *index, MPI_Status * status)
607 if (index == nullptr)
615 int rank_traced = simgrid::s4u::Actor::self()->getPid(); // FIXME: In PMPI_Wait, we check if the comm is null?
616 TRACE_smpi_comm_in(rank_traced, __FUNCTION__, new simgrid::instr::CpuTIData("waitAny", static_cast<double>(count)));
618 *index = simgrid::smpi::Request::waitany(count, requests, status);
620 if(*index!=MPI_UNDEFINED){
621 trace_smpi_recv_helper(&requests[*index], status);
622 TRACE_smpi_comm_out(rank_traced);
629 int PMPI_Waitall(int count, MPI_Request requests[], MPI_Status status[])
633 int rank_traced = simgrid::s4u::Actor::self()->getPid(); // FIXME: In PMPI_Wait, we check if the comm is null?
634 TRACE_smpi_comm_in(rank_traced, __FUNCTION__, new simgrid::instr::CpuTIData("waitAll", static_cast<double>(count)));
636 int retval = simgrid::smpi::Request::waitall(count, requests, status);
638 for (int i = 0; i < count; i++) {
639 trace_smpi_recv_helper(&requests[i], &status[i]);
641 TRACE_smpi_comm_out(rank_traced);
647 int PMPI_Waitsome(int incount, MPI_Request requests[], int *outcount, int *indices, MPI_Status status[])
652 if (outcount == nullptr) {
653 retval = MPI_ERR_ARG;
655 *outcount = simgrid::smpi::Request::waitsome(incount, requests, indices, status);
656 retval = MPI_SUCCESS;
662 int PMPI_Testsome(int incount, MPI_Request requests[], int* outcount, int* indices, MPI_Status status[])
667 if (outcount == nullptr) {
668 retval = MPI_ERR_ARG;
670 *outcount = simgrid::smpi::Request::testsome(incount, requests, indices, status);
671 retval = MPI_SUCCESS;
677 MPI_Request PMPI_Request_f2c(MPI_Fint request){
678 return static_cast<MPI_Request>(simgrid::smpi::Request::f2c(request));
681 MPI_Fint PMPI_Request_c2f(MPI_Request request) {
682 return request->c2f();