12 #include "Teuchos_StandardParameterEntryValidators.hpp"
23 if (sendType == DISTRIBUTOR_ISEND) {
26 else if (sendType == DISTRIBUTOR_SEND) {
29 else if (sendType == DISTRIBUTOR_ALLTOALL) {
32 #if defined(HAVE_TPETRACORE_MPI_ADVANCE)
33 else if (sendType == DISTRIBUTOR_MPIADVANCE_ALLTOALL) {
34 return "MpiAdvanceAlltoall";
36 else if (sendType == DISTRIBUTOR_MPIADVANCE_NBRALLTOALLV) {
37 return "MpiAdvanceNbralltoallv";
41 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid "
42 "EDistributorSendType enum value " << sendType <<
".");
49 if (s ==
"Isend")
return DISTRIBUTOR_ISEND;
50 if (s ==
"Send")
return DISTRIBUTOR_SEND;
51 if (s ==
"Alltoall")
return DISTRIBUTOR_ALLTOALL;
52 #if defined(HAVE_TPETRACORE_MPI_ADVANCE)
53 if (s ==
"MpiAdvanceAlltoall")
return DISTRIBUTOR_MPIADVANCE_ALLTOALL;
54 if (s ==
"MpiAdvanceNbralltoallv")
return DISTRIBUTOR_MPIADVANCE_NBRALLTOALLV;
56 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid string to convert to EDistributorSendType enum value: " << s);
62 if (std::find(valids.begin(), valids.end(), s) == valids.end()) {
63 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Invalid string for EDistributorSendType enum value: " << s);
72 case Details::DISTRIBUTOR_NOT_INITIALIZED:
73 return "Not initialized yet";
74 case Details::DISTRIBUTOR_INITIALIZED_BY_CREATE_FROM_SENDS:
75 return "By createFromSends";
76 case Details::DISTRIBUTOR_INITIALIZED_BY_CREATE_FROM_RECVS:
77 return "By createFromRecvs";
78 case Details::DISTRIBUTOR_INITIALIZED_BY_CREATE_FROM_SENDS_N_RECVS:
79 return "By createFromSendsAndRecvs";
80 case Details::DISTRIBUTOR_INITIALIZED_BY_REVERSE:
81 return "By createReverseDistributor";
82 case Details::DISTRIBUTOR_INITIALIZED_BY_COPY:
83 return "By copy constructor";
89 DistributorPlan::DistributorPlan(Teuchos::RCP<
const Teuchos::Comm<int>> comm)
91 #if defined(HAVE_TPETRACORE_MPI_ADVANCE)
92 mpixComm_(Teuchos::null),
94 howInitialized_(DISTRIBUTOR_NOT_INITIALIZED),
95 reversePlan_(Teuchos::null),
97 sendMessageToSelf_(false),
98 numSendsToOtherProcs_(0),
101 totalReceiveLength_(0)
104 DistributorPlan::DistributorPlan(
const DistributorPlan& otherPlan)
105 : comm_(otherPlan.comm_),
106 #if defined(HAVE_TPETRACORE_MPI_ADVANCE)
107 mpixComm_(otherPlan.mpixComm_),
109 howInitialized_(DISTRIBUTOR_INITIALIZED_BY_COPY),
110 reversePlan_(otherPlan.reversePlan_),
111 sendType_(otherPlan.sendType_),
112 sendMessageToSelf_(otherPlan.sendMessageToSelf_),
113 numSendsToOtherProcs_(otherPlan.numSendsToOtherProcs_),
114 procIdsToSendTo_(otherPlan.procIdsToSendTo_),
115 startsTo_(otherPlan.startsTo_),
116 lengthsTo_(otherPlan.lengthsTo_),
117 maxSendLength_(otherPlan.maxSendLength_),
118 indicesTo_(otherPlan.indicesTo_),
119 numReceives_(otherPlan.numReceives_),
120 totalReceiveLength_(otherPlan.totalReceiveLength_),
121 lengthsFrom_(otherPlan.lengthsFrom_),
122 procsFrom_(otherPlan.procsFrom_),
123 startsFrom_(otherPlan.startsFrom_),
124 indicesFrom_(otherPlan.indicesFrom_)
127 size_t DistributorPlan::createFromSends(
const Teuchos::ArrayView<const int>& exportProcIDs) {
128 using Teuchos::outArg;
129 using Teuchos::REDUCE_MAX;
130 using Teuchos::reduceAll;
132 const char rawPrefix[] =
"Tpetra::DistributorPlan::createFromSends";
134 const size_t numExports = exportProcIDs.size();
135 const int myProcID = comm_->getRank();
136 const int numProcs = comm_->getSize();
180 for (
size_t i = 0; i < numExports; ++i) {
181 const int exportID = exportProcIDs[i];
182 if (exportID >= numProcs || exportID < 0) {
188 reduceAll<int, int> (*comm_, REDUCE_MAX, badID, outArg (gbl_badID));
189 TEUCHOS_TEST_FOR_EXCEPTION
190 (gbl_badID >= 0, std::runtime_error, rawPrefix <<
"Proc "
191 << gbl_badID <<
", perhaps among other processes, got a bad "
207 Teuchos::Array<size_t> starts (numProcs + 1, 0);
210 size_t numActive = 0;
211 int needSendBuff = 0;
213 for (
size_t i = 0; i < numExports; ++i) {
214 const int exportID = exportProcIDs[i];
229 if (needSendBuff == 0 && starts[exportID] > 1 &&
230 exportID != exportProcIDs[i-1]) {
237 #if defined(HAVE_TPETRA_PRINT_EFFICIENCY_WARNINGS)
239 int global_needSendBuff;
240 reduceAll<int, int> (*comm_, REDUCE_MAX, needSendBuff,
241 outArg (global_needSendBuff));
243 global_needSendBuff != 0,
244 "::createFromSends: Grouping export IDs together by process rank often "
245 "improves performance.");
251 if (starts[myProcID] != 0) {
252 sendMessageToSelf_ =
true;
255 sendMessageToSelf_ =
false;
258 if (! needSendBuff) {
260 numSendsToOtherProcs_ = 0;
263 for (
int i = 0; i < numProcs; ++i) {
265 ++numSendsToOtherProcs_;
271 indicesTo_.resize(0);
274 procIdsToSendTo_.assign(numSendsToOtherProcs_,0);
275 startsTo_.assign(numSendsToOtherProcs_,0);
276 lengthsTo_.assign(numSendsToOtherProcs_,0);
283 size_t procIndex = 0;
284 for (
size_t i = 0; i < numSendsToOtherProcs_; ++i) {
285 while (exportProcIDs[procIndex] < 0) {
288 startsTo_[i] = procIndex;
289 int procID = exportProcIDs[procIndex];
290 procIdsToSendTo_[i] = procID;
291 procIndex += starts[procID];
296 if (numSendsToOtherProcs_ > 0) {
297 sort2(procIdsToSendTo_.begin(), procIdsToSendTo_.end(), startsTo_.begin());
301 for (
size_t i = 0; i < numSendsToOtherProcs_; ++i) {
302 int procID = procIdsToSendTo_[i];
303 lengthsTo_[i] = starts[procID];
304 if ((procID != myProcID) && (lengthsTo_[i] > maxSendLength_)) {
305 maxSendLength_ = lengthsTo_[i];
316 if (starts[0] == 0 ) {
317 numSendsToOtherProcs_ = 0;
320 numSendsToOtherProcs_ = 1;
322 for (Teuchos::Array<size_t>::iterator i=starts.begin()+1,
324 i != starts.end(); ++i)
326 if (*i != 0) ++numSendsToOtherProcs_;
332 for (Teuchos::Array<size_t>::reverse_iterator ip1=starts.rbegin(),
334 i != starts.rend(); ++i)
343 indicesTo_.resize(numActive);
345 for (
size_t i = 0; i < numExports; ++i) {
346 if (exportProcIDs[i] >= 0) {
348 indicesTo_[starts[exportProcIDs[i]]] = i;
350 ++starts[exportProcIDs[i]];
362 for (
int proc = numProcs-1; proc != 0; --proc) {
363 starts[proc] = starts[proc-1];
366 starts[numProcs] = numActive;
373 procIdsToSendTo_.resize(numSendsToOtherProcs_);
374 startsTo_.resize(numSendsToOtherProcs_);
375 lengthsTo_.resize(numSendsToOtherProcs_);
382 for (
int proc = 0; proc < numProcs; ++proc ) {
383 if (starts[proc+1] != starts[proc]) {
384 lengthsTo_[snd] = starts[proc+1] - starts[proc];
385 startsTo_[snd] = starts[proc];
387 if ((proc != myProcID) && (lengthsTo_[snd] > maxSendLength_)) {
388 maxSendLength_ = lengthsTo_[snd];
390 procIdsToSendTo_[snd] = proc;
396 if (sendMessageToSelf_) {
397 --numSendsToOtherProcs_;
405 howInitialized_ = Details::DISTRIBUTOR_INITIALIZED_BY_CREATE_FROM_SENDS;
407 return totalReceiveLength_;
410 void DistributorPlan::createFromRecvs(
const Teuchos::ArrayView<const int>& remoteProcIDs)
412 *
this = *getReversePlan();
413 howInitialized_ = Details::DISTRIBUTOR_INITIALIZED_BY_CREATE_FROM_RECVS;
416 void DistributorPlan::createFromSendsAndRecvs(
const Teuchos::ArrayView<const int>& exportProcIDs,
417 const Teuchos::ArrayView<const int>& remoteProcIDs)
426 howInitialized_ = Tpetra::Details::DISTRIBUTOR_INITIALIZED_BY_CREATE_FROM_SENDS_N_RECVS;
428 int myProcID = comm_->getRank ();
429 int numProcs = comm_->getSize();
431 const size_t numExportIDs = exportProcIDs.size();
432 Teuchos::Array<size_t> starts (numProcs + 1, 0);
434 size_t numActive = 0;
435 int needSendBuff = 0;
437 for(
size_t i = 0; i < numExportIDs; i++ )
439 if( needSendBuff==0 && i && (exportProcIDs[i] < exportProcIDs[i-1]) )
441 if( exportProcIDs[i] >= 0 )
443 ++starts[ exportProcIDs[i] ];
448 sendMessageToSelf_ = ( starts[myProcID] != 0 ) ? 1 : 0;
450 numSendsToOtherProcs_ = 0;
454 if (starts[0] == 0 ) {
455 numSendsToOtherProcs_ = 0;
458 numSendsToOtherProcs_ = 1;
460 for (Teuchos::Array<size_t>::iterator i=starts.begin()+1,
462 i != starts.end(); ++i)
464 if (*i != 0) ++numSendsToOtherProcs_;
470 for (Teuchos::Array<size_t>::reverse_iterator ip1=starts.rbegin(),
472 i != starts.rend(); ++i)
481 indicesTo_.resize(numActive);
483 for (
size_t i = 0; i < numExportIDs; ++i) {
484 if (exportProcIDs[i] >= 0) {
486 indicesTo_[starts[exportProcIDs[i]]] = i;
488 ++starts[exportProcIDs[i]];
491 for (
int proc = numProcs-1; proc != 0; --proc) {
492 starts[proc] = starts[proc-1];
495 starts[numProcs] = numActive;
496 procIdsToSendTo_.resize(numSendsToOtherProcs_);
497 startsTo_.resize(numSendsToOtherProcs_);
498 lengthsTo_.resize(numSendsToOtherProcs_);
501 for (
int proc = 0; proc < numProcs; ++proc ) {
502 if (starts[proc+1] != starts[proc]) {
503 lengthsTo_[snd] = starts[proc+1] - starts[proc];
504 startsTo_[snd] = starts[proc];
506 if ((proc != myProcID) && (lengthsTo_[snd] > maxSendLength_)) {
507 maxSendLength_ = lengthsTo_[snd];
509 procIdsToSendTo_[snd] = proc;
516 numSendsToOtherProcs_ = 0;
519 for (
int i = 0; i < numProcs; ++i) {
521 ++numSendsToOtherProcs_;
527 indicesTo_.resize(0);
530 procIdsToSendTo_.assign(numSendsToOtherProcs_,0);
531 startsTo_.assign(numSendsToOtherProcs_,0);
532 lengthsTo_.assign(numSendsToOtherProcs_,0);
539 size_t procIndex = 0;
540 for (
size_t i = 0; i < numSendsToOtherProcs_; ++i) {
541 while (exportProcIDs[procIndex] < 0) {
544 startsTo_[i] = procIndex;
545 int procID = exportProcIDs[procIndex];
546 procIdsToSendTo_[i] = procID;
547 procIndex += starts[procID];
552 if (numSendsToOtherProcs_ > 0) {
553 sort2(procIdsToSendTo_.begin(), procIdsToSendTo_.end(), startsTo_.begin());
557 for (
size_t i = 0; i < numSendsToOtherProcs_; ++i) {
558 int procID = procIdsToSendTo_[i];
559 lengthsTo_[i] = starts[procID];
560 if ((procID != myProcID) && (lengthsTo_[i] > maxSendLength_)) {
561 maxSendLength_ = lengthsTo_[i];
567 numSendsToOtherProcs_ -= sendMessageToSelf_;
568 std::vector<int> recv_list;
569 recv_list.reserve(numSendsToOtherProcs_);
572 for(
int i=0; i<remoteProcIDs.size(); i++) {
573 if(remoteProcIDs[i]>last_pid) {
574 recv_list.push_back(remoteProcIDs[i]);
575 last_pid = remoteProcIDs[i];
577 else if (remoteProcIDs[i]<last_pid)
578 throw std::runtime_error(
"Tpetra::Distributor:::createFromSendsAndRecvs expected RemotePIDs to be in sorted order");
580 numReceives_ = recv_list.size();
582 procsFrom_.assign(numReceives_,0);
583 lengthsFrom_.assign(numReceives_,0);
584 indicesFrom_.assign(numReceives_,0);
585 startsFrom_.assign(numReceives_,0);
587 for(
size_t i=0,j=0; i<numReceives_; ++i) {
589 procsFrom_[i] = recv_list[i];
591 for( ; j<(size_t)remoteProcIDs.size() &&
592 remoteProcIDs[jlast]==remoteProcIDs[j] ; j++){;}
593 lengthsFrom_[i] = j-jlast;
595 totalReceiveLength_ = remoteProcIDs.size();
596 indicesFrom_.clear ();
597 numReceives_-=sendMessageToSelf_;
600 Teuchos::RCP<DistributorPlan> DistributorPlan::getReversePlan()
const {
601 if (reversePlan_.is_null()) createReversePlan();
605 void DistributorPlan::createReversePlan()
const
607 reversePlan_ = Teuchos::rcp(
new DistributorPlan(comm_));
608 reversePlan_->howInitialized_ = Details::DISTRIBUTOR_INITIALIZED_BY_REVERSE;
609 reversePlan_->sendType_ = sendType_;
614 size_t totalSendLength =
615 std::accumulate(lengthsTo_.begin(), lengthsTo_.end(), 0);
620 size_t maxReceiveLength = 0;
621 const int myProcID = comm_->getRank();
622 for (
size_t i=0; i < numReceives_; ++i) {
623 if (procsFrom_[i] != myProcID) {
625 if (lengthsFrom_[i] > maxReceiveLength) {
626 maxReceiveLength = lengthsFrom_[i];
631 reversePlan_->sendMessageToSelf_ = sendMessageToSelf_;
632 reversePlan_->numSendsToOtherProcs_ = numReceives_;
633 reversePlan_->procIdsToSendTo_ = procsFrom_;
634 reversePlan_->startsTo_ = startsFrom_;
635 reversePlan_->lengthsTo_ = lengthsFrom_;
636 reversePlan_->maxSendLength_ = maxReceiveLength;
637 reversePlan_->indicesTo_ = indicesFrom_;
638 reversePlan_->numReceives_ = numSendsToOtherProcs_;
639 reversePlan_->totalReceiveLength_ = totalSendLength;
640 reversePlan_->lengthsFrom_ = lengthsTo_;
641 reversePlan_->procsFrom_ = procIdsToSendTo_;
642 reversePlan_->startsFrom_ = startsTo_;
643 reversePlan_->indicesFrom_ = indicesTo_;
645 #if defined(HAVE_TPETRACORE_MPI_ADVANCE)
647 reversePlan_->initializeMpiAdvance();
651 void DistributorPlan::computeReceives()
653 using Teuchos::Array;
654 using Teuchos::ArrayRCP;
656 using Teuchos::CommStatus;
657 using Teuchos::CommRequest;
658 using Teuchos::ireceive;
661 using Teuchos::REDUCE_SUM;
662 using Teuchos::receive;
663 using Teuchos::reduce;
664 using Teuchos::scatter;
666 using Teuchos::waitAll;
668 const int myRank = comm_->getRank();
669 const int numProcs = comm_->getSize();
671 const int mpiTag = DEFAULT_MPI_TAG;
679 Array<int> toProcsFromMe (numProcs, 0);
680 #ifdef HAVE_TPETRA_DEBUG
681 bool counting_error =
false;
682 #endif // HAVE_TPETRA_DEBUG
683 for (
size_t i = 0; i < (numSendsToOtherProcs_ + (sendMessageToSelf_ ? 1 : 0)); ++i) {
684 #ifdef HAVE_TPETRA_DEBUG
685 if (toProcsFromMe[procIdsToSendTo_[i]] != 0) {
686 counting_error =
true;
688 #endif // HAVE_TPETRA_DEBUG
689 toProcsFromMe[procIdsToSendTo_[i]] = 1;
691 #ifdef HAVE_TPETRA_DEBUG
694 "Tpetra::Distributor::computeReceives: There was an error on at least "
695 "one process in counting the number of messages send by that process to "
696 "the other processs. Please report this bug to the Tpetra developers.",
698 #endif // HAVE_TPETRA_DEBUG
753 Array<int> numRecvsOnEachProc;
754 if (myRank == root) {
755 numRecvsOnEachProc.resize (numProcs);
757 int numReceivesAsInt = 0;
758 reduce<int, int> (toProcsFromMe.getRawPtr (),
759 numRecvsOnEachProc.getRawPtr (),
760 numProcs, REDUCE_SUM, root, *comm_);
761 scatter<int, int> (numRecvsOnEachProc.getRawPtr (), 1,
762 &numReceivesAsInt, 1, root, *comm_);
763 numReceives_ =
static_cast<size_t> (numReceivesAsInt);
769 lengthsFrom_.assign (numReceives_, 0);
770 procsFrom_.assign (numReceives_, 0);
786 const size_t actualNumReceives = numReceives_ - (sendMessageToSelf_ ? 1 : 0);
792 Array<RCP<CommRequest<int> > > requests (actualNumReceives);
793 Array<ArrayRCP<size_t> > lengthsFromBuffers (actualNumReceives);
794 Array<RCP<CommStatus<int> > > statuses (actualNumReceives);
799 const int anySourceProc = MPI_ANY_SOURCE;
801 const int anySourceProc = -1;
805 for (
size_t i = 0; i < actualNumReceives; ++i) {
810 lengthsFromBuffers[i].resize (1);
811 lengthsFromBuffers[i][0] = as<size_t> (0);
812 requests[i] = ireceive<int, size_t> (lengthsFromBuffers[i], anySourceProc,
824 for (
size_t i = 0; i < numSendsToOtherProcs_ + (sendMessageToSelf_ ? 1 : 0); ++i) {
825 if (procIdsToSendTo_[i] != myRank) {
829 const size_t*
const lengthsTo_i = &lengthsTo_[i];
830 send<int, size_t> (lengthsTo_i, 1, as<int> (procIdsToSendTo_[i]), mpiTag, *comm_);
839 lengthsFrom_[numReceives_-1] = lengthsTo_[i];
840 procsFrom_[numReceives_-1] = myRank;
850 waitAll (*comm_, requests (), statuses ());
851 for (
size_t i = 0; i < actualNumReceives; ++i) {
852 lengthsFrom_[i] = *lengthsFromBuffers[i];
853 procsFrom_[i] = statuses[i]->getSourceRank ();
859 sort2 (procsFrom_.begin(), procsFrom_.end(), lengthsFrom_.begin());
862 totalReceiveLength_ =
863 std::accumulate (lengthsFrom_.begin (), lengthsFrom_.end (), 0);
864 indicesFrom_.clear ();
866 startsFrom_.clear ();
867 startsFrom_.reserve (numReceives_);
868 for (
size_t i = 0, j = 0; i < numReceives_; ++i) {
869 startsFrom_.push_back(j);
870 j += lengthsFrom_[i];
873 if (sendMessageToSelf_) {
878 void DistributorPlan::setParameterList(
const Teuchos::RCP<Teuchos::ParameterList>& plist)
880 using Teuchos::FancyOStream;
881 using Teuchos::getIntegralValue;
882 using Teuchos::ParameterList;
883 using Teuchos::parameterList;
887 if (! plist.is_null()) {
888 RCP<const ParameterList> validParams = getValidParameters ();
889 plist->validateParametersAndSetDefaults (*validParams);
892 getIntegralValue<Details::EDistributorSendType> (*plist,
"Send type");
895 sendType_ = sendType;
897 #if defined(HAVE_TPETRACORE_MPI_ADVANCE)
898 initializeMpiAdvance();
903 this->setMyParamList (plist);
909 Teuchos::Array<std::string> sendTypes;
910 sendTypes.push_back (
"Isend");
911 sendTypes.push_back (
"Send");
912 sendTypes.push_back (
"Alltoall");
913 #if defined(HAVE_TPETRACORE_MPI_ADVANCE)
914 sendTypes.push_back (
"MpiAdvanceAlltoall");
915 sendTypes.push_back (
"MpiAdvanceNbralltoallv");
921 Teuchos::Array<EDistributorSendType> res;
922 res.push_back (DISTRIBUTOR_ISEND);
923 res.push_back (DISTRIBUTOR_SEND);
924 res.push_back (DISTRIBUTOR_ALLTOALL);
925 #if defined(HAVE_TPETRACORE_MPI_ADVANCE)
926 res.push_back (DISTRIBUTOR_MPIADVANCE_ALLTOALL);
927 res.push_back (DISTRIBUTOR_MPIADVANCE_NBRALLTOALLV);
932 Teuchos::RCP<const Teuchos::ParameterList>
933 DistributorPlan::getValidParameters()
const
935 using Teuchos::Array;
936 using Teuchos::ParameterList;
937 using Teuchos::parameterList;
939 using Teuchos::setStringToIntegralParameter;
941 Array<std::string> sendTypes = distributorSendTypes ();
942 const Array<Details::EDistributorSendType> sendTypeEnums = distributorSendTypeEnums ();
944 const std::string validatedSendType = validSendTypeOrThrow(Behavior::defaultSendType());
946 RCP<ParameterList> plist = parameterList (
"Tpetra::Distributor");
948 setStringToIntegralParameter<Details::EDistributorSendType> (
"Send type",
949 validatedSendType,
"When using MPI, the variant of send to use in "
950 "do[Reverse]Posts()", sendTypes(), sendTypeEnums(), plist.getRawPtr());
951 plist->set (
"Timer Label",
"",
"Label for Time Monitor output");
953 return Teuchos::rcp_const_cast<
const ParameterList> (plist);
956 #if defined(HAVE_TPETRACORE_MPI_ADVANCE)
959 struct MpixCommDeallocator {
960 void free(MPIX_Comm **comm)
const {
961 MPIX_Comm_free(comm);
965 void DistributorPlan::initializeMpiAdvance() {
968 TEUCHOS_ASSERT(mpixComm_.is_null());
971 Teuchos::RCP<const Teuchos::MpiComm<int> > mpiComm = Teuchos::rcp_dynamic_cast<
const Teuchos::MpiComm<int> >(comm_);
972 Teuchos::RCP<const Teuchos::OpaqueWrapper<MPI_Comm> > rawComm = mpiComm->getRawMpiComm();
974 if (sendType_ == DISTRIBUTOR_MPIADVANCE_ALLTOALL ||
975 sendType_ == DISTRIBUTOR_MPIADVANCE_NBRALLTOALLV ) {
976 MPIX_Comm **mpixComm =
new(MPIX_Comm*);
977 err = MPIX_Comm_init(mpixComm, (*rawComm)());
978 mpixComm_ = Teuchos::RCP(mpixComm,
979 MpixCommDeallocator(),
984 TEUCHOS_ASSERT(err == 0);
989 DistributorPlan::SubViewLimits DistributorPlan::getImportViewLimits(
size_t numPackets)
const {
990 const size_t actualNumReceives = getNumReceives() + (hasSelfMessage() ? 1 : 0);
992 IndexView importStarts(actualNumReceives);
993 IndexView importLengths(actualNumReceives);
996 for (
size_t i = 0; i < actualNumReceives; ++i) {
997 importStarts[i] = offset;
998 offset += getLengthsFrom()[i] * numPackets;
999 importLengths[i] = getLengthsFrom()[i] * numPackets;
1001 return std::make_pair(importStarts, importLengths);
1004 DistributorPlan::SubViewLimits DistributorPlan::getImportViewLimits(
const Teuchos::ArrayView<const size_t> &numImportPacketsPerLID)
const {
1006 const size_t actualNumReceives = getNumReceives() + (hasSelfMessage() ? 1 : 0);
1008 IndexView importStarts(actualNumReceives);
1009 IndexView importLengths(actualNumReceives);
1012 size_t curLIDoffset = 0;
1013 for (
size_t i = 0; i < actualNumReceives; ++i) {
1014 size_t totalPacketsFrom_i = 0;
1015 for (
size_t j = 0; j < getLengthsFrom()[i]; ++j) {
1016 totalPacketsFrom_i += numImportPacketsPerLID[curLIDoffset + j];
1018 curLIDoffset += getLengthsFrom()[i];
1019 importStarts[i] = offset;
1020 offset += totalPacketsFrom_i;
1021 importLengths[i] = totalPacketsFrom_i;
1023 return std::make_pair(importStarts, importLengths);
1027 DistributorPlan::SubViewLimits DistributorPlan::getExportViewLimits(
size_t numPackets)
const {
1028 if (getIndicesTo().is_null()) {
1030 const size_t actualNumSends = getNumSends() + (hasSelfMessage() ? 1 : 0);
1031 IndexView exportStarts(actualNumSends);
1032 IndexView exportLengths(actualNumSends);
1033 for (
size_t pp = 0; pp < actualNumSends; ++pp) {
1034 exportStarts[pp] = getStartsTo()[pp] * numPackets;
1035 exportLengths[pp] = getLengthsTo()[pp] * numPackets;
1037 return std::make_pair(exportStarts, exportLengths);
1039 const size_t numIndices = getIndicesTo().size();
1040 IndexView exportStarts(numIndices);
1041 IndexView exportLengths(numIndices);
1042 for (
size_t j = 0; j < numIndices; ++j) {
1043 exportStarts[j] = getIndicesTo()[j]*numPackets;
1044 exportLengths[j] = numPackets;
1046 return std::make_pair(exportStarts, exportLengths);
1050 DistributorPlan::SubViewLimits DistributorPlan::getExportViewLimits(
const Teuchos::ArrayView<const size_t> &numExportPacketsPerLID)
const {
1051 if (getIndicesTo().is_null()) {
1052 const size_t actualNumSends = getNumSends() + (hasSelfMessage() ? 1 : 0);
1053 IndexView exportStarts(actualNumSends);
1054 IndexView exportLengths(actualNumSends);
1056 for (
size_t pp = 0; pp < actualNumSends; ++pp) {
1057 size_t numPackets = 0;
1058 for (
size_t j = getStartsTo()[pp];
1059 j < getStartsTo()[pp] + getLengthsTo()[pp]; ++j) {
1060 numPackets += numExportPacketsPerLID[j];
1062 exportStarts[pp] = offset;
1063 offset += numPackets;
1064 exportLengths[pp] = numPackets;
1066 return std::make_pair(exportStarts, exportLengths);
1068 const size_t numIndices = getIndicesTo().size();
1069 IndexView exportStarts(numIndices);
1070 IndexView exportLengths(numIndices);
1072 for (
size_t j = 0; j < numIndices; ++j) {
1073 exportStarts[j] = offset;
1074 offset += numExportPacketsPerLID[j];
1075 exportLengths[j] = numExportPacketsPerLID[j];
1077 return std::make_pair(exportStarts, exportLengths);
EDistributorHowInitialized
Enum indicating how and whether a Distributor was initialized.
static bool debug()
Whether Tpetra is in debug mode.
const std::string & validSendTypeOrThrow(const std::string &s)
Valid enum values of distributor send types.
std::string DistributorSendTypeEnumToString(EDistributorSendType sendType)
Convert an EDistributorSendType enum value to a string.
#define TPETRA_EFFICIENCY_WARNING(throw_exception_test, msg)
Print or throw an efficency warning.
EDistributorSendType DistributorSendTypeStringToEnum(const std::string_view s)
Convert a string to an EDistributorSendType. Throw on error.
Teuchos::Array< EDistributorSendType > distributorSendTypeEnums()
Valid enum values of distributor send types.
void sort2(const IT1 &first1, const IT1 &last1, const IT2 &first2, const bool stableSort=false)
Sort the first array, and apply the resulting permutation to the second array.
std::string DistributorHowInitializedEnumToString(EDistributorHowInitialized how)
Convert an EDistributorHowInitialized enum value to a string.
Stand-alone utility functions and macros.
Teuchos::Array< std::string > distributorSendTypes()
Valid string values for Distributor's "Send type" parameter.
#define SHARED_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg, comm)
Test for exception, with reduction over the given communicator.
EDistributorSendType
The type of MPI send that Distributor should use.
Declaration of Tpetra::Details::Behavior, a class that describes Tpetra's behavior.