42 #ifndef TPETRA_DETAILS_COOMATRIX_HPP
43 #define TPETRA_DETAILS_COOMATRIX_HPP
50 #include "TpetraCore_config.h"
51 #include "Tpetra_Details_PackTriples.hpp"
52 #include "Tpetra_DistObject.hpp"
55 #include "Teuchos_TypeNameTraits.hpp"
57 #include <initializer_list>
61 #include <type_traits>
74 template<
class IndexType>
84 template<
class IndexType>
90 return (x.row < y.row) || (x.row == y.row && x.col < y.col);
96 template<
class SC,
class GO>
106 typedef std::map<CooGraphEntry<GO>, SC,
111 entries_type entries_;
135 auto result = this->entries_.insert ({ij, val});
136 if (! result.second) {
137 result.first->second += val;
151 const GO gblColInds[],
153 const std::size_t numEnt)
155 for (std::size_t k = 0; k < numEnt; ++k) {
160 const SC val = vals[k];
161 auto result = this->entries_.insert ({ij, val});
162 if (! result.second) {
163 result.first->second += val;
172 return this->entries_.size ();
181 for (
auto iter = this->entries_.begin ();
182 iter != this->entries_.end (); ++iter) {
183 f (iter->first.row, iter->first.col, iter->second);
200 entries_type& tgtEntries = this->entries_;
201 const entries_type& srcEntries = src.entries_;
207 auto srcBeg = srcEntries.lower_bound ({srcGblRow, std::numeric_limits<GO>::min ()});
208 auto srcEnd = srcEntries.upper_bound ({srcGblRow+1, std::numeric_limits<GO>::min ()});
209 auto tgtBeg = tgtEntries.lower_bound ({tgtGblRow, std::numeric_limits<GO>::min ()});
210 auto tgtEnd = tgtEntries.upper_bound ({tgtGblRow+1, std::numeric_limits<GO>::min ()});
214 if (srcBeg != srcEnd) {
215 for (
auto tgtCur = tgtBeg; tgtCur != tgtEnd; ++tgtCur) {
216 auto srcCur = srcBeg;
217 while (srcCur != srcEnd && srcCur->first.col < tgtCur->first.col) {
229 if (srcCur != srcEnd) {
230 if (srcCur->first.col == tgtCur->first.col) {
231 tgtCur->second += srcCur->second;
238 (void) tgtEntries.insert (tgtCur, *srcCur);
257 const ::Teuchos::Comm<int>& comm,
258 std::ostream* errStrm = NULL)
const
260 using ::Tpetra::Details::countPackTriples;
261 using ::Tpetra::Details::countPackTriplesCount;
263 const char prefix[] =
"Tpetra::Details::Impl::CooMatrixImpl::countPackRow: ";
264 #ifdef HAVE_TPETRACORE_MPI
265 int errCode = MPI_SUCCESS;
268 #endif // HAVE_TPETRACORE_MPI
271 const GO minGO = std::numeric_limits<GO>::min ();
272 auto beg = this->entries_.lower_bound ({gblRow, minGO});
273 auto end = this->entries_.upper_bound ({gblRow+1, minGO});
275 for (
auto iter = beg; iter != end; ++iter) {
277 if (numEnt == std::numeric_limits<int>::max ()) {
278 if (errStrm != NULL) {
279 *errStrm << prefix <<
"In (global) row " << gblRow <<
", the number "
280 "of entries thus far, numEnt = " << numEnt <<
", has reached the "
281 "maximum int value. We need to store this count as int for MPI's "
288 int numPacketsForCount = 0;
293 if (errStrm != NULL) {
294 *errStrm << prefix <<
"countPackTriplesCount "
295 "returned errCode = " << errCode <<
" != 0." << endl;
299 if (numPacketsForCount < 0) {
300 if (errStrm != NULL) {
301 *errStrm << prefix <<
"countPackTriplesCount returned "
302 "numPacketsForCount = " << numPacketsForCount <<
" < 0." << endl;
308 int numPacketsForTriples = 0;
310 errCode = countPackTriples<SC, GO> (numEnt, comm, numPacketsForTriples);
311 TEUCHOS_TEST_FOR_EXCEPTION
312 (errCode != 0, std::runtime_error, prefix <<
"countPackTriples "
313 "returned errCode = " << errCode <<
" != 0.");
314 TEUCHOS_TEST_FOR_EXCEPTION
315 (numPacketsForTriples < 0, std::logic_error, prefix <<
"countPackTriples "
316 "returned numPacketsForTriples = " << numPacketsForTriples <<
" < 0.");
319 numPackets = numPacketsForCount + numPacketsForTriples;
339 const int outBufSize,
341 const ::Teuchos::Comm<int>& comm,
342 std::vector<GO>& gblRowInds,
343 std::vector<GO>& gblColInds,
344 std::vector<SC>& vals,
345 const GO gblRow)
const
347 using ::Tpetra::Details::packTriples;
348 using ::Tpetra::Details::packTriplesCount;
349 const char prefix[] =
"Tpetra::Details::Impl::CooMatrixImpl::packRow: ";
351 const GO minGO = std::numeric_limits<GO>::min ();
352 auto beg = this->entries_.lower_bound ({gblRow, minGO});
353 auto end = this->entries_.upper_bound ({gblRow+1, minGO});
357 gblRowInds.resize (0);
358 gblColInds.resize (0);
362 for (
auto iter = beg; iter != end; ++iter) {
363 gblRowInds.push_back (iter->first.row);
364 gblColInds.push_back (iter->first.col);
365 vals.push_back (iter->second);
367 TEUCHOS_TEST_FOR_EXCEPTION
368 (numEnt == std::numeric_limits<int>::max (), std::runtime_error, prefix
369 <<
"In (global) row " << gblRow <<
", the number of entries thus far, "
370 "numEnt = " << numEnt <<
", has reached the maximum int value. "
371 "We need to store this count as int for MPI's sake.");
377 TEUCHOS_TEST_FOR_EXCEPTION
378 (errCode != 0, std::runtime_error, prefix
379 <<
"In (global) row " << gblRow <<
", packTriplesCount returned "
380 "errCode = " << errCode <<
" != 0.");
392 TEUCHOS_TEST_FOR_EXCEPTION
393 (errCode != 0, std::runtime_error, prefix <<
"In (global) row "
394 << gblRow <<
", packTriples returned errCode = " << errCode
411 GO lclMinRowInd = std::numeric_limits<GO>::max ();
412 for (
typename entries_type::const_iterator iter = this->entries_.begin ();
413 iter != this->entries_.end (); ++iter) {
414 const GO rowInd = iter->first.row;
415 if (rowInd < lclMinRowInd) {
416 lclMinRowInd = rowInd;
418 if (rowInds.empty () || rowInds.back () != rowInd) {
419 rowInds.push_back (rowInd);
425 template<
class OffsetType>
427 buildCrs (std::vector<OffsetType>& rowOffsets,
431 static_assert (std::is_integral<OffsetType>::value,
432 "OffsetType must be a built-in integer type.");
440 rowOffsets.push_back (0);
443 typename entries_type::const_iterator iter = this->entries_.begin ();
444 GO prevGblRowInd = iter->first.row;
447 for ( ; iter != this->entries_.end (); ++iter, ++k) {
448 const GO gblRowInd = iter->first.row;
449 if (k == 0 || gblRowInd != prevGblRowInd) {
453 rowOffsets.push_back (k);
454 prevGblRowInd = gblRowInd;
456 gblColInds[k] = iter->first.col;
458 static_assert (std::is_same<
typename std::decay<decltype (iter->second)>::type, SC>::value,
459 "The type of iter->second != SC.");
460 vals[k] = iter->second;
462 rowOffsets.push_back (static_cast<OffsetType> (numEnt));
478 template<
class OffsetType,
class LO>
483 std::function<LO (
const GO)> gblToLcl)
const
485 static_assert (std::is_integral<OffsetType>::value,
486 "OffsetType must be a built-in integer type.");
487 static_assert (std::is_integral<LO>::value,
488 "LO must be a built-in integer type.");
496 rowOffsets.push_back (0);
499 typename entries_type::const_iterator iter = this->entries_.begin ();
500 GO prevGblRowInd = iter->first.row;
503 for ( ; iter != this->entries_.end (); ++iter, ++k) {
504 const GO gblRowInd = iter->first.row;
505 if (k == 0 || gblRowInd != prevGblRowInd) {
509 rowOffsets.push_back (k);
510 prevGblRowInd = gblRowInd;
512 lclColInds[k] = gblToLcl (iter->first.col);
513 vals[k] = iter->second;
515 rowOffsets.push_back (static_cast<OffsetType> (numEnt));
580 typedef LO local_ordinal_type;
581 typedef GO global_ordinal_type;
582 typedef NT node_type;
583 typedef typename NT::device_type device_type;
585 typedef ::Tpetra::Map<local_ordinal_type, global_ordinal_type, node_type>
map_type;
589 typedef ::Tpetra::DistObject<
packet_type, local_ordinal_type,
590 global_ordinal_type, node_type>
base_type;
603 localError_ (new bool (false)),
604 errs_ (new std::shared_ptr<std::ostringstream> ())
616 localError_ (new bool (false)),
617 errs_ (new std::shared_ptr<std::ostringstream> ())
634 this->impl_.sumIntoGlobalValue (gblRowInd, gblColInd, val);
647 const GO gblColInds[],
649 const std::size_t numEnt)
651 this->impl_.sumIntoGlobalValues (gblRowInds, gblColInds, vals, numEnt);
657 std::initializer_list<GO> gblColInds,
658 std::initializer_list<SC> vals,
659 const std::size_t numEnt)
661 this->impl_.sumIntoGlobalValues (gblRowInds.begin (), gblColInds.begin (),
662 vals.begin (), numEnt);
669 return this->impl_.getLclNumEntries ();
672 template<
class OffsetType>
674 buildCrs (::Kokkos::View<OffsetType*, device_type>& rowOffsets,
675 ::Kokkos::View<GO*, device_type>& gblColInds,
676 ::Kokkos::View<typename ::Kokkos::ArithTraits<SC>::val_type*, device_type>& vals)
678 static_assert (std::is_integral<OffsetType>::value,
679 "OffsetType must be a built-in integer type.");
680 using ::Kokkos::create_mirror_view;
682 using ::Kokkos::View;
683 typedef typename ::Kokkos::ArithTraits<SC>::val_type ISC;
687 gblColInds = View<GO*, device_type> (
"gblColInds", numEnt);
688 vals = View<ISC*, device_type> (
"vals", numEnt);
689 auto gblColInds_h = create_mirror_view (gblColInds);
690 auto vals_h = create_mirror_view (vals);
692 std::vector<std::size_t> rowOffsetsSV;
693 this->impl_.buildCrs (rowOffsetsSV,
694 gblColInds_h.data (),
697 View<OffsetType*, device_type> (
"rowOffsets", rowOffsetsSV.size ());
698 typename View<OffsetType*, device_type>::HostMirror
699 rowOffsets_h (rowOffsetsSV.data (), rowOffsetsSV.size ());
719 if (comm.is_null ()) {
720 this->
map_ = ::Teuchos::null;
723 this->
map_ = this->buildMap (comm);
738 TEUCHOS_TEST_FOR_EXCEPTION
739 (this->
getMap ().is_null (), std::runtime_error,
"Tpetra::Details::"
740 "CooMatrix::fillComplete: This object does not yet have a Map. "
741 "You must call the version of fillComplete "
742 "that takes an input communicator.");
782 return ((*errs_).get () == NULL) ? std::string (
"") : (*errs_)->str ();
798 std::shared_ptr<bool> localError_;
807 std::shared_ptr<std::shared_ptr<std::ostringstream> > errs_;
811 markLocalErrorAndGetStream ()
813 * (this->localError_) =
true;
814 if ((*errs_).get () == NULL) {
815 *errs_ = std::shared_ptr<std::ostringstream> (
new std::ostringstream ());
824 using Teuchos::TypeNameTraits;
826 std::ostringstream os;
827 os <<
"\"Tpetra::Details::CooMatrix\": {"
828 <<
"SC: " << TypeNameTraits<SC>::name ()
829 <<
", LO: " << TypeNameTraits<LO>::name ()
830 <<
", GO: " << TypeNameTraits<GO>::name ()
831 <<
", NT: " << TypeNameTraits<NT>::name ();
832 if (this->getObjectLabel () !=
"") {
833 os <<
", Label: \"" << this->getObjectLabel () <<
"\"";
835 os <<
", Has Map: " << (this->
map_.is_null () ?
"false" :
"true")
844 const Teuchos::EVerbosityLevel verbLevel =
845 Teuchos::Describable::verbLevel_default)
const
847 using ::Tpetra::Details::gathervPrint;
848 using ::Teuchos::EVerbosityLevel;
849 using ::Teuchos::OSTab;
850 using ::Teuchos::TypeNameTraits;
851 using ::Teuchos::VERB_DEFAULT;
852 using ::Teuchos::VERB_LOW;
853 using ::Teuchos::VERB_MEDIUM;
856 const auto vl = (verbLevel == VERB_DEFAULT) ? VERB_LOW : verbLevel;
858 auto comm = this->
getMap ().is_null () ?
859 Teuchos::null : this->
getMap ()->getComm ();
860 const int myRank = comm.is_null () ? 0 : comm->getRank ();
863 if (vl != Teuchos::VERB_NONE) {
867 out <<
"\"Tpetra::Details::CooMatrix\":" << endl;
871 out <<
"Template parameters:" << endl;
874 out <<
"SC: " << TypeNameTraits<SC>::name () << endl
875 <<
"LO: " << TypeNameTraits<LO>::name () << endl
876 <<
"GO: " << TypeNameTraits<GO>::name () << endl
877 <<
"NT: " << TypeNameTraits<NT>::name () << endl;
879 if (this->getObjectLabel () !=
"") {
880 out <<
"Label: \"" << this->getObjectLabel () <<
"\"" << endl;
882 out <<
"Has Map: " << (this->
map_.is_null () ?
"false" :
"true") << endl;
886 if (! this->
map_.is_null ()) {
888 out <<
"Map:" << endl;
891 this->
map_->describe (out, vl);
896 std::ostringstream os;
897 os <<
"Process " << myRank <<
":" << endl;
900 const std::size_t numEnt = this->impl_.getLclNumEntries ();
901 os <<
" Local number of entries: " << numEnt << endl;
902 if (vl > VERB_MEDIUM) {
903 os <<
" Local entries:" << endl;
905 this->impl_.forAllEntries ([&os] (
const GO row,
const GO col,
const SC& val) {
927 Teuchos::RCP<const map_type>
928 buildMap (const ::Teuchos::RCP<const ::Teuchos::Comm<int> >& comm)
930 using ::Teuchos::outArg;
931 using ::Teuchos::rcp;
932 using ::Teuchos::REDUCE_MIN;
933 using ::Teuchos::reduceAll;
934 typedef ::Tpetra::global_size_t GST;
938 if (comm.is_null ()) {
939 return ::Teuchos::null;
950 std::vector<GO> rowInds;
951 const GO lclMinGblRowInd = this->impl_.getMyGlobalRowIndices (rowInds);
954 GO gblMinGblRowInd = 0;
955 reduceAll<int, GO> (*comm, REDUCE_MIN, lclMinGblRowInd,
956 outArg (gblMinGblRowInd));
957 const GO indexBase = gblMinGblRowInd;
958 const GST INV = Tpetra::Details::OrdinalTraits<GST>::invalid ();
959 return rcp (
new map_type (INV, rowInds.data (), rowInds.size (),
968 return static_cast<size_t> (0);
979 const char prefix[] =
"Tpetra::Details::CooMatrix::checkSizes: ";
981 const this_COO_type* src =
dynamic_cast<const this_COO_type*
> (&source);
984 std::ostream& err = markLocalErrorAndGetStream ();
985 err << prefix <<
"The source object of the Import or Export "
986 "must be a CooMatrix with the same template parameters as the "
987 "target object." << endl;
989 else if (this->
map_.is_null ()) {
990 std::ostream& err = markLocalErrorAndGetStream ();
991 err << prefix <<
"The target object of the Import or Export "
992 "must be a CooMatrix with a nonnull Map." << endl;
994 return ! (* (this->localError_));
999 typename ::Tpetra::DistObject<char, LO, GO, NT>::buffer_device_type;
1003 (const ::Tpetra::SrcDistObject& sourceObject,
1004 const size_t numSameIDs,
1005 const Kokkos::DualView<
const LO*,
1007 const Kokkos::DualView<
const LO*,
1013 const char prefix[] =
"Tpetra::Details::CooMatrix::copyAndPermute: ";
1018 if (* (this->localError_)) {
1019 std::ostream& err = this->markLocalErrorAndGetStream ();
1020 err << prefix <<
"The target object of the Import or Export is "
1021 "already in an error state." << endl;
1025 const this_COO_type* src =
dynamic_cast<const this_COO_type*
> (&sourceObject);
1026 if (src ==
nullptr) {
1027 std::ostream& err = this->markLocalErrorAndGetStream ();
1028 err << prefix <<
"Input argument 'sourceObject' is not a CooMatrix."
1033 const size_t numPermuteIDs =
1034 static_cast<size_t> (permuteToLIDs.extent (0));
1035 if (numPermuteIDs != static_cast<size_t> (permuteFromLIDs.extent (0))) {
1036 std::ostream& err = this->markLocalErrorAndGetStream ();
1037 err << prefix <<
"permuteToLIDs.extent(0) = "
1038 << numPermuteIDs <<
" != permuteFromLIDs.extent(0) = "
1039 << permuteFromLIDs.extent (0) <<
"." << endl;
1042 if (
sizeof (
int) <=
sizeof (
size_t) &&
1043 numPermuteIDs > static_cast<size_t> (std::numeric_limits<int>::max ())) {
1044 std::ostream& err = this->markLocalErrorAndGetStream ();
1045 err << prefix <<
"numPermuteIDs = " << numPermuteIDs
1046 <<
", a size_t, overflows int." << endl;
1053 if (
sizeof (
int) <=
sizeof (
size_t) &&
1054 numSameIDs > static_cast<size_t> (std::numeric_limits<int>::max ())) {
1055 std::ostream& err = this->markLocalErrorAndGetStream ();
1056 err << prefix <<
"numSameIDs = " << numSameIDs
1057 <<
", a size_t, overflows int." << endl;
1064 const LO numSame =
static_cast<int> (numSameIDs);
1069 LO numInvalidSameRows = 0;
1070 for (LO lclRow = 0; lclRow < numSame; ++lclRow) {
1074 const GO gblRow = this->
map_->getGlobalElement (lclRow);
1075 if (gblRow == ::Tpetra::Details::OrdinalTraits<GO>::invalid ()) {
1076 ++numInvalidSameRows;
1080 this->impl_.mergeIntoRow (gblRow, src->impl_, gblRow);
1089 const LO numPermute =
static_cast<int> (numPermuteIDs);
1094 LO numInvalidRowsFrom = 0;
1099 LO numInvalidRowsTo = 0;
1101 TEUCHOS_ASSERT( ! permuteFromLIDs.need_sync_host () );
1102 TEUCHOS_ASSERT( ! permuteToLIDs.need_sync_host () );
1103 auto permuteFromLIDs_h = permuteFromLIDs.view_host ();
1104 auto permuteToLIDs_h = permuteToLIDs.view_host ();
1106 for (LO k = 0; k < numPermute; ++k) {
1107 const LO lclRowFrom = permuteFromLIDs_h[k];
1108 const LO lclRowTo = permuteToLIDs_h[k];
1109 const GO gblRowFrom = src->map_->getGlobalElement (lclRowFrom);
1110 const GO gblRowTo = this->
map_->getGlobalElement (lclRowTo);
1112 bool bothConversionsValid =
true;
1113 if (gblRowFrom == ::Tpetra::Details::OrdinalTraits<GO>::invalid ()) {
1114 ++numInvalidRowsFrom;
1115 bothConversionsValid =
false;
1117 if (gblRowTo == ::Tpetra::Details::OrdinalTraits<GO>::invalid ()) {
1119 bothConversionsValid =
false;
1121 if (bothConversionsValid) {
1122 this->impl_.mergeIntoRow (gblRowTo, src->impl_, gblRowFrom);
1127 if (numInvalidSameRows != 0 || numInvalidRowsFrom != 0 ||
1128 numInvalidRowsTo != 0) {
1131 std::vector<std::pair<LO, GO> > invalidSameRows;
1132 invalidSameRows.reserve (numInvalidSameRows);
1133 std::vector<std::pair<LO, GO> > invalidRowsFrom;
1134 invalidRowsFrom.reserve (numInvalidRowsFrom);
1135 std::vector<std::pair<LO, GO> > invalidRowsTo;
1136 invalidRowsTo.reserve (numInvalidRowsTo);
1138 for (LO lclRow = 0; lclRow < numSame; ++lclRow) {
1142 const GO gblRow = this->
map_->getGlobalElement (lclRow);
1143 if (gblRow == ::Tpetra::Details::OrdinalTraits<GO>::invalid ()) {
1144 invalidSameRows.push_back ({lclRow, gblRow});
1148 for (LO k = 0; k < numPermute; ++k) {
1149 const LO lclRowFrom = permuteFromLIDs_h[k];
1150 const LO lclRowTo = permuteToLIDs_h[k];
1151 const GO gblRowFrom = src->map_->getGlobalElement (lclRowFrom);
1152 const GO gblRowTo = this->
map_->getGlobalElement (lclRowTo);
1154 if (gblRowFrom == ::Tpetra::Details::OrdinalTraits<GO>::invalid ()) {
1155 invalidRowsFrom.push_back ({lclRowFrom, gblRowFrom});
1157 if (gblRowTo == ::Tpetra::Details::OrdinalTraits<GO>::invalid ()) {
1158 invalidRowsTo.push_back ({lclRowTo, gblRowTo});
1162 std::ostringstream os;
1163 if (numInvalidSameRows != 0) {
1164 os <<
"Invalid permute \"same\" (local, global) index pairs: [";
1165 for (std::size_t k = 0; k < invalidSameRows.size (); ++k) {
1166 const auto& p = invalidSameRows[k];
1167 os <<
"(" << p.first <<
"," << p.second <<
")";
1168 if (k + 1 < invalidSameRows.size ()) {
1172 os <<
"]" << std::endl;
1174 if (numInvalidRowsFrom != 0) {
1175 os <<
"Invalid permute \"from\" (local, global) index pairs: [";
1176 for (std::size_t k = 0; k < invalidRowsFrom.size (); ++k) {
1177 const auto& p = invalidRowsFrom[k];
1178 os <<
"(" << p.first <<
"," << p.second <<
")";
1179 if (k + 1 < invalidRowsFrom.size ()) {
1183 os <<
"]" << std::endl;
1185 if (numInvalidRowsTo != 0) {
1186 os <<
"Invalid permute \"to\" (local, global) index pairs: [";
1187 for (std::size_t k = 0; k < invalidRowsTo.size (); ++k) {
1188 const auto& p = invalidRowsTo[k];
1189 os <<
"(" << p.first <<
"," << p.second <<
")";
1190 if (k + 1 < invalidRowsTo.size ()) {
1194 os <<
"]" << std::endl;
1197 std::ostream& err = this->markLocalErrorAndGetStream ();
1198 err << prefix << os.str ();
1205 (const ::Tpetra::SrcDistObject& sourceObject,
1206 const Kokkos::DualView<
const local_ordinal_type*,
1210 Kokkos::DualView<
size_t*,
1212 size_t& constantNumPackets)
1214 using Teuchos::Comm;
1218 const char prefix[] =
"Tpetra::Details::CooMatrix::packAndPrepare: ";
1219 const char suffix[] =
" This should never happen. "
1220 "Please report this bug to the Tpetra developers.";
1224 constantNumPackets = 0;
1226 const this_COO_type* src =
dynamic_cast<const this_COO_type*
> (&sourceObject);
1227 if (src ==
nullptr) {
1228 std::ostream& err = this->markLocalErrorAndGetStream ();
1229 err << prefix <<
"Input argument 'sourceObject' is not a CooMatrix."
1232 else if (* (src->localError_)) {
1233 std::ostream& err = this->markLocalErrorAndGetStream ();
1234 err << prefix <<
"The source (input) object of the Import or Export "
1235 "is already in an error state on this process."
1238 else if (* (this->localError_)) {
1239 std::ostream& err = this->markLocalErrorAndGetStream ();
1240 err << prefix <<
"The target (output, \"this\") object of the Import "
1241 "or Export is already in an error state on this process." << endl;
1246 if (* (this->localError_)) {
1251 auto numPacketsPerLID_tmp = numPacketsPerLID;
1252 numPacketsPerLID_tmp.sync_host ();
1253 numPacketsPerLID_tmp.modify_host ();
1260 const size_t numExports = exportLIDs.extent (0);
1261 if (numExports == 0) {
1265 RCP<const Comm<int> > comm = src->getMap ().is_null () ?
1266 Teuchos::null : src->getMap ()->getComm ();
1267 if (comm.is_null () || comm->getSize () == 1) {
1268 if (numExports != static_cast<size_t> (0)) {
1269 std::ostream& err = this->markLocalErrorAndGetStream ();
1270 err << prefix <<
"The input communicator is either null or "
1271 "has only one process, but numExports = " << numExports <<
" != 0."
1282 numPacketsPerLID.sync_host ();
1283 numPacketsPerLID.modify_host ();
1285 TEUCHOS_ASSERT( ! exportLIDs.need_sync_host () );
1286 auto exportLIDs_h = exportLIDs.view_host ();
1288 int totalNumPackets = 0;
1289 size_t numInvalidRowInds = 0;
1290 std::ostringstream errStrm;
1291 for (
size_t k = 0; k < numExports; ++k) {
1292 const LO lclRow = exportLIDs_h[k];
1295 const GO gblRow = src->map_->getGlobalElement (lclRow);
1296 if (gblRow == ::Tpetra::Details::OrdinalTraits<GO>::invalid ()) {
1298 ++numInvalidRowInds;
1299 numPacketsPerLID.h_view[k] = 0;
1307 src->impl_.countPackRow (numPackets, gblRow, *comm, &errStrm);
1309 std::ostream& err = this->markLocalErrorAndGetStream ();
1310 err << prefix << errStrm.str () << endl;
1311 numPacketsPerLID.h_view[k] = 0;
1317 const long long newTotalNumPackets =
1318 static_cast<long long> (totalNumPackets) +
1319 static_cast<long long> (numPackets);
1320 if (newTotalNumPackets >
1321 static_cast<long long> (std::numeric_limits<int>::max ())) {
1322 std::ostream& err = this->markLocalErrorAndGetStream ();
1323 err << prefix <<
"The new total number of packets "
1324 << newTotalNumPackets <<
" does not fit in int." << endl;
1328 for (
size_t k2 = k; k2 < numExports; ++k2) {
1329 numPacketsPerLID.h_view[k2] = 0;
1333 numPacketsPerLID.h_view[k] =
static_cast<size_t> (numPackets);
1334 totalNumPackets =
static_cast<int> (newTotalNumPackets);
1339 if (numInvalidRowInds != 0) {
1340 std::vector<std::pair<LO, GO> > invalidRowInds;
1341 for (
size_t k = 0; k < numExports; ++k) {
1342 const LO lclRow = exportLIDs_h[k];
1346 const GO gblRow = src->map_->getGlobalElement (lclRow);
1347 if (gblRow == ::Tpetra::Details::OrdinalTraits<GO>::invalid ()) {
1348 invalidRowInds.push_back ({lclRow, gblRow});
1351 std::ostringstream os;
1352 os << prefix <<
"We found " << numInvalidRowInds <<
" invalid row ind"
1353 << (numInvalidRowInds !=
static_cast<size_t> (1) ?
"ices" :
"ex")
1354 <<
" out of " << numExports <<
" in exportLIDs. Here is the list "
1355 <<
" of invalid row indices: [";
1356 for (
size_t k = 0; k < invalidRowInds.size (); ++k) {
1357 os <<
"(LID: " << invalidRowInds[k].first <<
", GID: "
1358 << invalidRowInds[k].second <<
")";
1359 if (k + 1 < invalidRowInds.size ()) {
1365 std::ostream& err = this->markLocalErrorAndGetStream ();
1366 err << prefix << os.str () << std::endl;
1371 const bool reallocated =
1373 "CooMatrix exports");
1375 exports.sync_host ();
1378 exports.modify_host ();
1382 std::vector<GO> gblRowInds;
1383 std::vector<GO> gblColInds;
1384 std::vector<SC> vals;
1386 int outBufCurPos = 0;
1388 for (
size_t k = 0; k < numExports; ++k) {
1389 const LO lclRow = exportLIDs.h_view[k];
1392 const GO gblRow = src->map_->getGlobalElement (lclRow);
1394 src->impl_.packRow (outBuf, totalNumPackets, outBufCurPos, *comm,
1395 gblRowInds, gblColInds, vals, gblRow);
1401 (
const Kokkos::DualView<
const local_ordinal_type*,
1405 Kokkos::DualView<
size_t*,
1410 using Teuchos::Comm;
1413 const char prefix[] =
"Tpetra::Details::CooMatrix::unpackAndCombine: ";
1414 const char suffix[] =
" This should never happen. "
1415 "Please report this bug to the Tpetra developers.";
1417 TEUCHOS_ASSERT( ! importLIDs.need_sync_host () );
1418 auto importLIDs_h = importLIDs.view_host ();
1420 const std::size_t numImports = importLIDs.extent (0);
1421 if (numImports == 0) {
1424 else if (imports.extent (0) == 0) {
1425 std::ostream& err = this->markLocalErrorAndGetStream ();
1426 err << prefix <<
"importLIDs.extent(0) = " << numImports <<
" != 0, "
1427 <<
"but imports.extent(0) = 0. This doesn't make sense, because "
1428 <<
"for every incoming LID, CooMatrix packs at least the count of "
1429 <<
"triples associated with that LID, even if the count is zero. "
1430 <<
"importLIDs = [";
1431 for (std::size_t k = 0; k < numImports; ++k) {
1432 err << importLIDs_h[k];
1433 if (k + 1 < numImports) {
1437 err <<
"]. " << suffix << endl;
1441 RCP<const Comm<int> > comm = this->
getMap ().is_null () ?
1442 Teuchos::null : this->
getMap ()->getComm ();
1443 if (comm.is_null () || comm->getSize () == 1) {
1444 if (numImports != static_cast<size_t> (0)) {
1445 std::ostream& err = this->markLocalErrorAndGetStream ();
1446 err << prefix <<
"The input communicator is either null or "
1447 "has only one process, but numImports = " << numImports <<
" != 0."
1460 if (static_cast<size_t> (imports.extent (0)) >
1461 static_cast<size_t> (std::numeric_limits<int>::max ())) {
1462 std::ostream& err = this->markLocalErrorAndGetStream ();
1463 err << prefix <<
"imports.extent(0) = "
1464 << imports.extent (0) <<
" does not fit in int." << endl;
1467 const int inBufSize =
static_cast<int> (imports.extent (0));
1469 if (imports.need_sync_host ()) {
1470 imports.sync_host ();
1472 if (numPacketsPerLID.need_sync_host ()) {
1473 numPacketsPerLID.sync_host ();
1475 auto imports_h = imports.view_host ();
1476 auto numPacketsPerLID_h = numPacketsPerLID.view_host ();
1480 std::vector<GO> gblRowInds;
1481 std::vector<GO> gblColInds;
1482 std::vector<SC> vals;
1485 int inBufCurPos = 0;
1486 size_t numInvalidRowInds = 0;
1488 std::ostringstream errStrm;
1489 for (
size_t k = 0; k < numImports; ++k) {
1490 const LO lclRow = importLIDs_h(k);
1491 const GO gblRow = this->
map_->getGlobalElement (lclRow);
1492 if (gblRow == ::Tpetra::Details::OrdinalTraits<GO>::invalid ()) {
1493 ++numInvalidRowInds;
1499 const int origInBufCurPos = inBufCurPos;
1503 numEnt, *comm, &errStrm);
1504 if (errCode != 0 || numEnt < 0 || inBufCurPos < origInBufCurPos) {
1505 std::ostream& err = this->markLocalErrorAndGetStream ();
1507 err << prefix <<
"In unpack loop, k=" << k <<
": ";
1509 err <<
" unpackTriplesCount returned errCode = " << errCode
1510 <<
" != 0." << endl;
1513 err <<
" unpackTriplesCount returned errCode = 0, but numEnt = "
1514 << numEnt <<
" < 0." << endl;
1516 if (inBufCurPos < origInBufCurPos) {
1517 err <<
" After unpackTriplesCount, inBufCurPos = " << inBufCurPos
1518 <<
" < origInBufCurPos = " << origInBufCurPos <<
"." << endl;
1520 err <<
" unpackTriplesCount report: " << errStrm.str () << endl;
1521 err << suffix << endl;
1528 #ifdef HAVE_TPETRA_DEBUG
1535 #endif // HAVE_TPETRA_DEBUG
1540 gblRowInds.resize (numEnt);
1541 gblColInds.resize (numEnt);
1542 vals.resize (numEnt);
1545 gblRowInds.data (), gblColInds.data (),
1546 vals.data (), numEnt, *comm, &errStrm);
1548 std::ostream& err = this->markLocalErrorAndGetStream ();
1549 err << prefix <<
"unpackTriples returned errCode = "
1550 << errCode <<
" != 0. It reports: " << errStrm.str ()
1557 #ifdef HAVE_TPETRA_DEBUG
1564 #endif // HAVE_TPETRA_DEBUG
1567 vals.data (), numEnt);
1572 if (numInvalidRowInds != 0) {
1576 std::ostream& err = this->markLocalErrorAndGetStream ();
1578 std::vector<std::pair<LO, GO> > invalidRowInds;
1579 for (
size_t k = 0; k < numImports; ++k) {
1580 const LO lclRow = importLIDs_h(k);
1581 const GO gblRow = this->
map_->getGlobalElement (lclRow);
1582 if (gblRow == ::Tpetra::Details::OrdinalTraits<GO>::invalid ()) {
1583 invalidRowInds.push_back ({lclRow, gblRow});
1587 err << prefix <<
"We found " << numInvalidRowInds <<
" invalid row ind"
1588 << (numInvalidRowInds !=
static_cast<size_t> (1) ?
"ices" :
"ex")
1589 <<
" out of " << numImports <<
" in importLIDs. Here is the list "
1590 <<
" of invalid row indices: [";
1591 for (
size_t k = 0; k < invalidRowInds.size (); ++k) {
1592 err <<
"(LID: " << invalidRowInds[k].first <<
", GID: "
1593 << invalidRowInds[k].second <<
")";
1594 if (k + 1 < invalidRowInds.size ()) {
1607 #endif // TPETRA_DETAILS_COOMATRIX_HPP
char packet_type
Type for packing and unpacking data.
int packTriplesCount(const int, char[], const int, int &, const ::Teuchos::Comm< int > &, std::ostream *errStrm)
Pack the count (number) of matrix triples.
void sumIntoGlobalValues(std::initializer_list< GO > gblRowInds, std::initializer_list< GO > gblColInds, std::initializer_list< SC > vals, const std::size_t numEnt)
Initializer-list overload of the above method (which see).
Node node_type
The Node type. If you don't know what this is, don't use it.
CooMatrixImpl()=default
Default constructor.
std::string errorMessages() const
The current stream of error messages.
Declaration of a function that prints strings from each process.
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print a descriptiion of this object to the given output stream; overrides Teuchos::Describable method...
void packRow(packet_type outBuf[], const int outBufSize, int &outBufCurPos, const ::Teuchos::Comm< int > &comm, std::vector< GO > &gblRowInds, std::vector< GO > &gblColInds, std::vector< SC > &vals, const GO gblRow) const
Pack the given row of the matrix.
Function comparing two CooGraphEntry structs, lexicographically, first by row index, then by column index.
void sumIntoGlobalValue(const GO gblRowInd, const GO gblColInd, const SC &val)
Insert one entry locally into the sparse matrix, if it does not exist there yet. If it does exist...
Sparse matrix used only for file input / output.
void buildLocallyIndexedCrs(std::vector< OffsetType > &rowOffsets, LO lclColInds[], SC vals[], std::function< LO(const GO)> gblToLcl) const
Build a locally indexed version of CRS storage.
void fillComplete()
Special version of fillComplete that assumes that the matrix already has a Map, and reuses its commun...
GlobalOrdinal global_ordinal_type
The type of global indices.
virtual ~CooMatrix()
Destructor (virtual for memory safety of derived classes).
void gathervPrint(std::ostream &out, const std::string &s, const Teuchos::Comm< int > &comm)
On Process 0 in the given communicator, print strings from each process in that communicator, in rank order.
void fillComplete(const ::Teuchos::RCP< const ::Teuchos::Comm< int > > &comm)
Tell the matrix that you are done inserting entries locally, and that the matrix should build its Map...
void deep_copy(MultiVector< DS, DL, DG, DN > &dst, const MultiVector< SS, SL, SG, SN > &src)
Copy the contents of the MultiVector src into dst.
std::size_t getLclNumEntries() const
Number of entries in the sparse matrix on the calling process.
void sumIntoGlobalValue(const GO gblRowInd, const GO gblColInd, const SC &val)
Insert one entry locally into the sparse matrix, if it does not exist there yet. If it does exist...
Type of each (row index, column index) pair in the Tpetra::Details::CooMatrix (see below)...
bool localError() const
Whether this object had an error on the calling process.
int unpackTriplesCount(const char[], const int, int &, int &, const ::Teuchos::Comm< int > &, std::ostream *errStrm)
Unpack just the count of triples from the given input buffer.
SC scalar_type
Type of each entry (value) in the sparse matrix.
CooMatrix()
Default constructor.
virtual size_t constantNumberOfPackets() const
By returning 0, tell DistObject that this class may not necessarily have a constant number of "packet...
void sumIntoGlobalValues(const GO gblRowInds[], const GO gblColInds[], const SC vals[], const std::size_t numEnt)
Insert multiple entries locally into the sparse matrix.
CooMatrix(const ::Teuchos::RCP< const map_type > &map)
Constructor that takes a Map.
virtual bool checkSizes(const ::Tpetra::SrcDistObject &source)
Compare the source and target (this) objects for compatibility.
void sumIntoGlobalValues(const GO gblRowInds[], const GO gblColInds[], const SC vals[], const std::size_t numEnt)
Insert multiple entries locally into the sparse matrix.
CombineMode
Rule for combining data in an Import or Export.
Teuchos::RCP< const map_type > map_
The Map over which this object is distributed.
bool reallocDualViewIfNeeded(Kokkos::DualView< ValueType *, DeviceType > &dv, const size_t newSize, const char newLabel[], const size_t tooBigFactor=2, const bool needFenceBeforeRealloc=true)
Reallocate the DualView in/out argument, if needed.
int countPackTriplesCount(const ::Teuchos::Comm< int > &, int &size, std::ostream *errStrm)
Compute the buffer size required by packTriples for packing the number of matrix entries ("triples")...
::Tpetra::Map< local_ordinal_type, global_ordinal_type, node_type > map_type
Type of the Map specialization to give to the constructor.
Declaration and definition of Tpetra::Details::reallocDualViewIfNeeded, an implementation detail of T...
void forAllEntries(std::function< void(const GO, const GO, const SC &)> f) const
Execute the given function for all entries of the sparse matrix, sequentially (no thread parallelism)...
LocalOrdinal local_ordinal_type
The type of local indices.
int countPackRow(int &numPackets, const GO gblRow, const ::Teuchos::Comm< int > &comm, std::ostream *errStrm=NULL) const
Count the number of packets (bytes, in this case) needed to pack the given row of the matrix...
std::size_t getLclNumEntries() const
Number of entries in the sparse matrix on the calling process.
GO getMyGlobalRowIndices(std::vector< GO > &rowInds) const
Get the global row indices on this process, sorted and made unique, and return the minimum global row...
typename::Tpetra::DistObject< char, LO, GO, NT >::buffer_device_type buffer_device_type
Kokkos::Device specialization for DistObject communication buffers.
void mergeIntoRow(const GO tgtGblRow, const CooMatrixImpl< SC, GO > &src, const GO srcGblRow)
Into global row tgtGblRow of *this, merge global row srcGblRow of src.
char packet_type
This class transfers data as bytes (MPI_BYTE).
virtual std::string description() const
One-line descriptiion of this object; overrides Teuchos::Describable method.
virtual Teuchos::RCP< const map_type > getMap() const
The Map describing the parallel distribution of this object.
int unpackTriples(const char[], const int, int &, OrdinalType[], OrdinalType[], ScalarType[], const int, const ::Teuchos::Comm< int > &, std::ostream *errStrm=NULL)
Unpack matrix entries ("triples" (i, j, A(i,j))) from the given input buffer.
Base class for distributed Tpetra objects that support data redistribution.
Implementation detail of Tpetra::Details::CooMatrix (which see below).
int packTriples(const OrdinalType[], const OrdinalType[], const ScalarType[], const int, char[], const int, int &, const ::Teuchos::Comm< int > &, std::ostream *errStrm=NULL)
Pack matrix entries ("triples" (i, j, A(i,j))) into the given output buffer.