40 #ifndef TPETRA_LOCALCRSMATRIXOPERATOR_DECL_HPP
41 #define TPETRA_LOCALCRSMATRIXOPERATOR_DECL_HPP
43 #include "Tpetra_LocalCrsMatrixOperator_fwd.hpp"
44 #include "Tpetra_LocalOperator.hpp"
45 #include "KokkosSparse_CrsMatrix.hpp"
59 template<
class MultiVectorScalar,
class MatrixScalar,
class Device>
63 using mv_scalar_type =
65 using matrix_scalar_type =
66 typename LocalOperator<MatrixScalar, Device>::scalar_type;
71 using local_ordinal_type =
73 using execution_space =
typename Device::execution_space;
75 using local_matrix_device_type =
76 KokkosSparse::CrsMatrix<matrix_scalar_type,
83 using local_cusparse_matrix_type =
84 KokkosSparse::CrsMatrix<matrix_scalar_type,
89 using local_graph_device_type =
typename local_matrix_device_type::StaticCrsGraphType;
92 using ordinal_view_type =
typename local_graph_device_type::entries_type::non_const_type;
95 LocalCrsMatrixOperator (
const std::shared_ptr<local_matrix_device_type>& A,
const ordinal_view_type& A_ordinal_rowptrs);
99 apply (Kokkos::View<
const mv_scalar_type**, array_layout,
100 device_type, Kokkos::MemoryTraits<Kokkos::Unmanaged> > X,
101 Kokkos::View<mv_scalar_type**, array_layout,
102 device_type, Kokkos::MemoryTraits<Kokkos::Unmanaged> > Y,
103 const Teuchos::ETransp mode,
104 const mv_scalar_type alpha,
105 const mv_scalar_type beta)
const override;
109 Kokkos::View<
const mv_scalar_type**, array_layout,
110 device_type, Kokkos::MemoryTraits<Kokkos::Unmanaged> > X,
111 Kokkos::View<mv_scalar_type**, array_layout,
112 device_type, Kokkos::MemoryTraits<Kokkos::Unmanaged> > Y,
113 const Teuchos::ETransp mode,
114 const mv_scalar_type alpha,
115 const mv_scalar_type beta)
const;
117 bool hasTransposeApply ()
const override;
119 const local_matrix_device_type& getLocalMatrixDevice ()
const;
122 std::shared_ptr<local_matrix_device_type> A_;
123 local_cusparse_matrix_type A_cusparse;
124 const bool have_A_cusparse;
129 #endif // TPETRA_LOCALCRSMATRIXOPERATOR_DECL_HPP
Abstract interface for local operators (e.g., matrices and preconditioners).
Abstract interface for local operators (e.g., matrices and preconditioners).
int local_ordinal_type
Default value of Scalar template parameter.
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...