44 #ifndef TPETRA_DETAILS_GETDIAGCOPYWITHOUTOFFSETS_DEF_HPP 45 #define TPETRA_DETAILS_GETDIAGCOPYWITHOUTOFFSETS_DEF_HPP 55 #include "Tpetra_Details_gathervPrint.hpp" 56 #include "Tpetra_RowGraph.hpp" 57 #include "Tpetra_CrsGraph.hpp" 58 #include "Tpetra_RowMatrix.hpp" 59 #include "Tpetra_Vector.hpp" 72 template<
class SC,
class LO,
class GO,
class NT>
73 class GetLocalDiagCopyWithoutOffsetsNotFillCompleteFunctor {
75 typedef ::Tpetra::RowMatrix<SC, LO, GO, NT> row_matrix_type;
76 typedef ::Tpetra::Vector<SC, LO, GO, NT> vec_type;
83 typedef typename vec_type::dual_view_type::host_mirror_space::execution_space host_execution_space;
87 graphIsSorted (
const row_matrix_type& A)
90 using Teuchos::rcp_dynamic_cast;
99 RCP<const row_graph_type> G_row = A.getGraph ();
100 if (! G_row.is_null ()) {
101 RCP<const crs_graph_type> G_crs =
102 rcp_dynamic_cast<
const crs_graph_type> (G_row);
103 if (! G_crs.is_null ()) {
104 sorted = G_crs->isSorted ();
113 GetLocalDiagCopyWithoutOffsetsNotFillCompleteFunctor (LO& lclNumErrs,
115 const row_matrix_type& A) :
117 lclRowMap_ (A.getRowMap ()->getLocalMap ()),
118 lclColMap_ (A.getColMap ()->getLocalMap ()),
119 sorted_ (graphIsSorted (A))
121 const LO lclNumRows =
static_cast<LO
> (diag.getLocalLength ());
123 const LO matLclNumRows =
124 static_cast<LO
> (lclRowMap_.getNodeNumElements ());
125 TEUCHOS_TEST_FOR_EXCEPTION
126 (lclNumRows != matLclNumRows, std::invalid_argument,
127 "diag.getLocalLength() = " << lclNumRows <<
" != " 128 "A.getRowMap()->getNodeNumElements() = " << matLclNumRows <<
".");
133 diag.template modify<Kokkos::HostSpace> ();
134 D_lcl_ = diag.template getLocalView<Kokkos::HostSpace> ();
135 D_lcl_1d_ = Kokkos::subview (D_lcl_, Kokkos::ALL (), 0);
137 Kokkos::RangePolicy<host_execution_space, LO> range (0, lclNumRows);
139 Kokkos::parallel_reduce (range, *
this, lclNumErrs);
143 diag.template sync<typename device_type::memory_space> ();
146 void operator () (
const LO& lclRowInd, LO& errCount)
const {
147 using KokkosSparse::findRelOffset;
149 D_lcl_1d_(lclRowInd) = Kokkos::Details::ArithTraits<IST>::zero ();
150 const GO gblInd = lclRowMap_.getGlobalElement (lclRowInd);
151 const LO lclColInd = lclColMap_.getLocalElement (gblInd);
158 const LO* lclColInds;
160 const LO err = A_.getLocalRowViewRaw (lclRowInd, numEnt, lclColInds, curVals);
169 findRelOffset (lclColInds, numEnt, lclColInd, hint, sorted_);
170 if (offset == numEnt) {
174 D_lcl_1d_(lclRowInd) = curVals[offset];
181 const row_matrix_type& A_;
184 typename vec_type::dual_view_type::t_host D_lcl_;
185 decltype (Kokkos::subview (D_lcl_, Kokkos::ALL (), 0)) D_lcl_1d_;
190 template<class SC, class LO, class GO, class NT>
196 using ::Tpetra::Details::gathervPrint;
197 using Teuchos::outArg;
198 using Teuchos::REDUCE_MIN;
199 using Teuchos::reduceAll;
200 typedef GetLocalDiagCopyWithoutOffsetsNotFillCompleteFunctor<SC,
201 LO, GO, NT> functor_type;
211 std::ostringstream errStrm;
212 Teuchos::RCP<const Teuchos::Comm<int> > commPtr = A.getComm ();
213 if (commPtr.is_null ()) {
216 const Teuchos::Comm<int>& comm = *commPtr;
219 functor_type functor (lclNumErrs, diag, A);
221 catch (std::exception& e) {
223 errStrm <<
"Process " << A.getComm ()->getRank () <<
": " 224 << e.what () << std::endl;
226 if (lclNumErrs != 0) {
230 reduceAll<int, int> (comm, REDUCE_MIN, lclSuccess, outArg (gblSuccess));
231 if (gblSuccess == -1) {
232 if (comm.getRank () == 0) {
240 std::cerr <<
"getLocalDiagCopyWithoutOffsetsNotFillComplete threw an " 241 "exception on one or more MPI processes in the matrix's comunicator." 247 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::runtime_error,
"");
249 else if (gblSuccess == 0) {
250 TEUCHOS_TEST_FOR_EXCEPTION
251 (gblSuccess != 1, std::runtime_error,
252 "getLocalDiagCopyWithoutOffsetsNotFillComplete failed on " 253 "one or more MPI processes in the matrix's communicator.");
257 functor_type functor (lclNumErrs, diag, A);
269 #define TPETRA_DETAILS_GETDIAGCOPYWITHOUTOFFSETS_INSTANT( SCALAR, LO, GO, NODE ) \ 271 Details::getLocalDiagCopyWithoutOffsetsNotFillComplete< SCALAR, LO, GO, NODE > \ 272 ( ::Tpetra::Vector< SCALAR, LO, GO, NODE >& diag, \ 273 const ::Tpetra::RowMatrix< SCALAR, LO, GO, NODE >& A, \ 276 #endif // TPETRA_DETAILS_GETDIAGCOPYWITHOUTOFFSETS_DEF_HPP Namespace Tpetra contains the class and methods constituting the Tpetra library.
base_type::impl_scalar_type impl_scalar_type
The type used internally in place of Scalar.
An abstract interface for graphs accessed by rows.
Details::LocalMap< LocalOrdinal, GlobalOrdinal, device_type > local_map_type
Type of the "local" Map.
Implementation details of Tpetra.
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.
Traits class for "invalid" (flag) values of integer types that Tpetra uses as local ordinals or globa...
base_type::map_type map_type
The type of the Map specialization used by this class.
LO getLocalDiagCopyWithoutOffsetsNotFillComplete(::Tpetra::Vector< SC, LO, GO, NT > &diag, const ::Tpetra::RowMatrix< SC, LO, GO, NT > &A, const bool debug=false)
Given a locally indexed, global sparse matrix, extract the matrix's diagonal entries into a Tpetra::V...
A distributed graph accessed by rows (adjacency lists) and stored sparsely.
A read-only, row-oriented interface to a sparse matrix.
A distributed dense vector.
Node::device_type device_type
The Kokkos device type.