47 #ifndef PACKAGES_XPETRA_CRSMATRIX_UTILS_HPP_
48 #define PACKAGES_XPETRA_CRSMATRIX_UTILS_HPP_
52 #include "Xpetra_Map.hpp"
54 #ifdef HAVE_XPETRA_EPETRA
55 #include "Epetra_Util.h"
58 #ifdef HAVE_XPETRA_TPETRA
59 #include "Tpetra_Import_Util2.hpp"
71 template <
class Scalar,
74 class Node = Tpetra::KokkosClassic::DefaultNode::DefaultNodeType>
75 class CrsMatrixUtils {
76 #undef XPETRA_CRSMATRIXUTILS_SHORT
82 const Teuchos::ArrayView<LocalOrdinal>& CRS_colind,
83 const Teuchos::ArrayView<Scalar>& CRS_vals,
86 #if defined(HAVE_XPETRA_EPETRA)
87 throw(
Xpetra::Exceptions::RuntimeError(
"Xpetra::CrsMatrixUtils::sortCrsEntries only available for GO=int or GO=long long with EpetraNode (Serial or OpenMP depending on configuration)"));
88 #endif // HAVE_XPETRA_EPETRA
90 #ifdef HAVE_XPETRA_TPETRA
91 Tpetra::Import_Util::sortCrsEntries(CRS_rowptr, CRS_colind, CRS_vals);
92 #endif // HAVE_XPETRA_TPETRA
101 const Teuchos::ArrayView<LocalOrdinal>& CRS_colind,
102 const Teuchos::ArrayView<Scalar>& CRS_vals,
105 #if defined(HAVE_XPETRA_EPETRA)
106 throw(
Xpetra::Exceptions::RuntimeError(
"Xpetra::CrsMatrixUtils::sortAndMergeCrsEntries only available for GO=int or GO=long long with EpetraNode (Serial or OpenMP depending on configuration)"));
107 #endif // HAVE_XPETRA_EPETRA
109 #ifdef HAVE_XPETRA_TPETRA
110 Tpetra::Import_Util::sortAndMergeCrsEntries(CRS_rowptr, CRS_colind, CRS_vals);
111 #endif // HAVE_XPETRA_TPETRA
119 #ifdef HAVE_XPETRA_EPETRA
127 #undef XPETRA_CRSMATRIXUTILS_SHORT
133 const Teuchos::ArrayView<LocalOrdinal>& CRS_colind,
134 const Teuchos::ArrayView<Scalar>& CRS_vals,
137 #if defined(HAVE_XPETRA_EPETRA)
138 int rv = Epetra_Util::SortCrsEntries(Teuchos::as<int>(CRS_rowptr.size() - 1),
139 CRS_rowptr.getRawPtr(),
140 CRS_colind.getRawPtr(),
141 CRS_vals.getRawPtr());
145 #endif // HAVE_XPETRA_EPETRA
147 #ifdef HAVE_XPETRA_TPETRA
148 Tpetra::Import_Util::sortCrsEntries(CRS_rowptr, CRS_colind, CRS_vals);
149 #endif // HAVE_XPETRA_TPETRA
158 const Teuchos::ArrayView<LocalOrdinal>& CRS_colind,
159 const Teuchos::ArrayView<Scalar>& CRS_vals,
162 #if defined(HAVE_XPETRA_EPETRA)
163 int rv = Epetra_Util::SortAndMergeCrsEntries(Teuchos::as<int>(CRS_rowptr.size() - 1),
164 CRS_rowptr.getRawPtr(),
165 CRS_colind.getRawPtr(),
166 CRS_vals.getRawPtr());
170 #endif // HAVE_XPETRA_EPETRA
172 #ifdef HAVE_XPETRA_TPETRA
173 Tpetra::Import_Util::sortAndMergeCrsEntries(CRS_rowptr, CRS_colind, CRS_vals);
174 #endif // HAVE_XPETRA_TPETRA
189 #undef XPETRA_CRSMATRIXUTILS_SHORT
195 const Teuchos::ArrayView<LocalOrdinal>& CRS_colind,
196 const Teuchos::ArrayView<Scalar>& CRS_vals,
199 #if defined(HAVE_XPETRA_EPETRA)
200 int rv = Epetra_Util::SortCrsEntries(Teuchos::as<int>(CRS_rowptr.size() - 1),
201 CRS_rowptr.getRawPtr(),
202 CRS_colind.getRawPtr(),
203 CRS_vals.getRawPtr());
207 #endif // HAVE_XPETRA_EPETRA
209 #ifdef HAVE_XPETRA_TPETRA
210 Tpetra::Import_Util::sortCrsEntries(CRS_rowptr, CRS_colind, CRS_vals);
211 #endif // HAVE_XPETRA_TPETRA
220 const Teuchos::ArrayView<LocalOrdinal>& CRS_colind,
221 const Teuchos::ArrayView<Scalar>& CRS_vals,
224 #if defined(HAVE_XPETRA_EPETRA)
225 int rv = Epetra_Util::SortAndMergeCrsEntries(Teuchos::as<int>(CRS_rowptr.size() - 1),
226 CRS_rowptr.getRawPtr(),
227 CRS_colind.getRawPtr(),
228 CRS_vals.getRawPtr());
232 #endif // HAVE_XPETRA_EPETRA
234 #ifdef HAVE_XPETRA_TPETRA
235 Tpetra::Import_Util::sortAndMergeCrsEntries(CRS_rowptr, CRS_colind, CRS_vals);
236 #endif // HAVE_XPETRA_TPETRA
243 #endif // HAVE_XPETRA_EPETRA for Epetra scpecialization
247 #define XPETRA_CRSMATRIXUTILS_SHORT
249 #endif // PACKAGES_XPETRA_CRSMATRIX_UTILS_HPP_
static void sortAndMergeCrsEntries(const Teuchos::ArrayView< size_t > &CRS_rowptr, const Teuchos::ArrayView< LocalOrdinal > &CRS_colind, const Teuchos::ArrayView< Scalar > &CRS_vals, const UnderlyingLib lib)
Sort and merge the entries of the (raw CSR) matrix by column index within each row.
std::string toString(Xpetra::UnderlyingLib lib)
Convert a Xpetra::UnderlyingLib to a std::string.
Xpetra utility class for CrsMatrix-related routines.
Exception throws to report errors in the internal logical of the program.
static void sortCrsEntries(const Teuchos::ArrayView< size_t > &CRS_rowptr, const Teuchos::ArrayView< LocalOrdinal > &CRS_colind, const Teuchos::ArrayView< Scalar > &CRS_vals, const UnderlyingLib lib)
Sort the entries of the (raw CSR) matrix by column index within each row.
static void sortCrsEntries(const Teuchos::ArrayView< size_t > &CRS_rowptr, const Teuchos::ArrayView< LocalOrdinal > &CRS_colind, const Teuchos::ArrayView< Scalar > &CRS_vals, const UnderlyingLib lib)
Sort the entries of the (raw CSR) matrix by column index within each row.
static void sortAndMergeCrsEntries(const Teuchos::ArrayView< size_t > &CRS_rowptr, const Teuchos::ArrayView< LocalOrdinal > &CRS_colind, const Teuchos::ArrayView< Scalar > &CRS_vals, const UnderlyingLib lib)
Sort and merge the entries of the (raw CSR) matrix by column index within each row.
static void sortCrsEntries(const Teuchos::ArrayView< size_t > &CRS_rowptr, const Teuchos::ArrayView< LocalOrdinal > &CRS_colind, const Teuchos::ArrayView< Scalar > &CRS_vals, const UnderlyingLib lib)
Sort the entries of the (raw CSR) matrix by column index within each row.
Tpetra::KokkosCompat::KokkosSerialWrapperNode EpetraNode
static void sortAndMergeCrsEntries(const Teuchos::ArrayView< size_t > &CRS_rowptr, const Teuchos::ArrayView< LocalOrdinal > &CRS_colind, const Teuchos::ArrayView< Scalar > &CRS_vals, const UnderlyingLib lib)
Sort and merge the entries of the (raw CSR) matrix by column index within each row.