40 #ifndef TPETRA_LOCALCRSMATRIXOPERATOR_DEF_HPP
41 #define TPETRA_LOCALCRSMATRIXOPERATOR_DEF_HPP
43 #include "Tpetra_LocalOperator.hpp"
45 #include "KokkosSparse.hpp"
46 #include "Teuchos_TestForException.hpp"
47 #include "Teuchos_OrdinalTraits.hpp"
51 template<
class MultiVectorScalar,
class MatrixScalar,
class Device>
52 LocalCrsMatrixOperator<MultiVectorScalar, MatrixScalar, Device>::
53 LocalCrsMatrixOperator (
const std::shared_ptr<local_matrix_device_type>& A)
54 : A_ (A), have_A_cusparse(false)
56 const char tfecfFuncName[] =
"LocalCrsMatrixOperator: ";
57 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC
58 (A_.get () ==
nullptr, std::invalid_argument,
59 "Input matrix A is null.");
62 template<
class MultiVectorScalar,
class MatrixScalar,
class Device>
63 LocalCrsMatrixOperator<MultiVectorScalar, MatrixScalar, Device>::
64 LocalCrsMatrixOperator (
const std::shared_ptr<local_matrix_device_type>& A,
const ordinal_view_type& A_ordinal_rowptrs) :
66 A_cusparse(
"LocalCrsMatrixOperator_cuSPARSE", A->numRows(), A->numCols(), A->nnz(),
67 A->values, A_ordinal_rowptrs, A->graph.entries),
70 const char tfecfFuncName[] =
"LocalCrsMatrixOperator: ";
71 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC
72 (A_.get () ==
nullptr, std::invalid_argument,
73 "Input matrix A is null.");
76 template<
class MultiVectorScalar,
class MatrixScalar,
class Device>
78 LocalCrsMatrixOperator<MultiVectorScalar, MatrixScalar, Device>::
79 hasTransposeApply ()
const
84 template<
class MultiVectorScalar,
class MatrixScalar,
class Device>
86 LocalCrsMatrixOperator<MultiVectorScalar, MatrixScalar, Device>::
87 apply (Kokkos::View<
const mv_scalar_type**, array_layout,
88 device_type, Kokkos::MemoryTraits<Kokkos::Unmanaged> > X,
89 Kokkos::View<mv_scalar_type**, array_layout,
90 device_type, Kokkos::MemoryTraits<Kokkos::Unmanaged> > Y,
91 const Teuchos::ETransp mode,
92 const mv_scalar_type alpha,
93 const mv_scalar_type beta)
const
95 const bool conjugate = (mode == Teuchos::CONJ_TRANS);
96 const bool transpose = (mode != Teuchos::NO_TRANS);
98 #ifdef HAVE_TPETRA_DEBUG
99 const char tfecfFuncName[] =
"apply: ";
101 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC
102 (X.extent (1) != Y.extent (1), std::runtime_error,
103 "X.extent(1) = " << X.extent (1) <<
" != Y.extent(1) = "
104 << Y.extent (1) <<
".");
107 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC
108 (X.data () == Y.data () && X.data () !=
nullptr,
109 std::runtime_error,
"X and Y may not alias one another.");
110 #endif // HAVE_TPETRA_DEBUG
112 const auto op = transpose ?
113 (conjugate ? KokkosSparse::ConjugateTranspose :
114 KokkosSparse::Transpose) : KokkosSparse::NoTranspose;
117 KokkosSparse::spmv (op, alpha, A_cusparse, X, beta, Y);
121 KokkosSparse::spmv (op, alpha, *A_, X, beta, Y);
127 template<
class MultiVectorScalar,
class MatrixScalar,
class Device>
131 Kokkos::View<
const mv_scalar_type**, array_layout,
132 device_type, Kokkos::MemoryTraits<Kokkos::Unmanaged> > X,
133 Kokkos::View<mv_scalar_type**, array_layout,
134 device_type, Kokkos::MemoryTraits<Kokkos::Unmanaged> > Y,
135 const Teuchos::ETransp mode,
136 const mv_scalar_type alpha,
137 const mv_scalar_type beta)
const
139 apply(X, Y, mode, alpha, beta);
142 template<
class MultiVectorScalar,
class MatrixScalar,
class Device>
143 const typename LocalCrsMatrixOperator<MultiVectorScalar, MatrixScalar, Device>::local_matrix_device_type&
161 #define TPETRA_LOCALCRSMATRIXOPERATOR_INSTANT(SC,NT) \
162 template class LocalCrsMatrixOperator< SC, SC, NT::device_type >;
166 #define TPETRA_LOCALCRSMATRIXOPERATOR_MIXED_INSTANT(SC,MATSC,LO,GO,NT) \
167 template class LocalCrsMatrixOperator< SC, MATSC, NT::device_type >;
169 #endif // TPETRA_LOCALCRSMATRIXOPERATOR_DEF_HPP
Abstract interface for local operators (e.g., matrices and preconditioners).
void applyImbalancedRows(Kokkos::View< const mv_scalar_type **, array_layout, device_type, Kokkos::MemoryTraits< Kokkos::Unmanaged > > X, Kokkos::View< mv_scalar_type **, array_layout, device_type, Kokkos::MemoryTraits< Kokkos::Unmanaged > > Y, const Teuchos::ETransp mode, const mv_scalar_type alpha, const mv_scalar_type beta) const
Same behavior as apply() above, except give KokkosKernels a hint to use an SPMV algorithm that can ef...
Declaration of Tpetra::Details::Behavior, a class that describes Tpetra's behavior.