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"
51template<
class MultiVectorScalar,
class MatrixScalar,
class Device>
52LocalCrsMatrixOperator<MultiVectorScalar, MatrixScalar, Device>::
53LocalCrsMatrixOperator (
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.");
62template<
class MultiVectorScalar,
class MatrixScalar,
class Device>
63LocalCrsMatrixOperator<MultiVectorScalar, MatrixScalar, Device>::
64LocalCrsMatrixOperator (
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.");
76template<
class MultiVectorScalar,
class MatrixScalar,
class Device>
78LocalCrsMatrixOperator<MultiVectorScalar, MatrixScalar, Device>::
79hasTransposeApply ()
const
84template<
class MultiVectorScalar,
class MatrixScalar,
class Device>
86LocalCrsMatrixOperator<MultiVectorScalar, MatrixScalar, Device>::
87apply (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.");
112 const auto op = transpose ?
113 (conjugate ? KokkosSparse::ConjugateTranspose :
114 KokkosSparse::Transpose) : KokkosSparse::NoTranspose;
117 if(X.extent(1) == size_t(1) && have_A_cusparse)
119 KokkosSparse::spmv (op, alpha, A_cusparse, Kokkos::subview(X, Kokkos::ALL(), 0),
120 beta, Kokkos::subview(Y, Kokkos::ALL(), 0));
124 KokkosSparse::spmv (op, alpha, *A_, X, beta, Y);
130template<
class MultiVectorScalar,
class MatrixScalar,
class Device>
134 Kokkos::View<
const mv_scalar_type**, array_layout,
135 device_type, Kokkos::MemoryTraits<Kokkos::Unmanaged> > X,
136 Kokkos::View<mv_scalar_type**, array_layout,
137 device_type, Kokkos::MemoryTraits<Kokkos::Unmanaged> > Y,
138 const Teuchos::ETransp mode,
139 const mv_scalar_type alpha,
140 const mv_scalar_type beta)
const
142 const bool conjugate = (mode == Teuchos::CONJ_TRANS);
143 const bool transpose = (mode != Teuchos::NO_TRANS);
145#ifdef HAVE_TPETRA_DEBUG
146 const char tfecfFuncName[] =
"applyLoadBalanced: ";
148 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC
149 (X.extent (1) != Y.extent (1), std::runtime_error,
150 "X.extent(1) = " << X.extent (1) <<
" != Y.extent(1) = "
151 << Y.extent (1) <<
".");
154 TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC
155 (X.data () == Y.data () && X.data () !=
nullptr,
156 std::runtime_error,
"X and Y may not alias one another.");
159 const auto op = transpose ?
160 (conjugate ? KokkosSparse::ConjugateTranspose :
161 KokkosSparse::Transpose) : KokkosSparse::NoTranspose;
166 if((Details::Behavior::useMergePathMultiVector() || X.extent(1) ==
size_t(1)) && have_A_cusparse)
168 KokkosKernels::Experimental::Controls controls;
169 controls.setParameter(
"algorithm",
"merge");
171 for(
size_t vec = 0; vec < X.extent(1); vec++)
173 KokkosSparse::spmv (controls, op,
174 alpha, A_cusparse, Kokkos::subview(X, Kokkos::ALL(), vec),
175 beta, Kokkos::subview(Y, Kokkos::ALL(), vec));
181 KokkosSparse::spmv (op, alpha, *A_, X, beta, Y);
185template<
class MultiVectorScalar,
class MatrixScalar,
class Device>
186const typename LocalCrsMatrixOperator<MultiVectorScalar, MatrixScalar, Device>::local_matrix_device_type&
204#define TPETRA_LOCALCRSMATRIXOPERATOR_INSTANT(SC,NT) \
205 template class LocalCrsMatrixOperator< SC, SC, NT::device_type >;
209#define TPETRA_LOCALCRSMATRIXOPERATOR_MIXED_INSTANT(SC,MATSC,LO,GO,NT) \
210 template class LocalCrsMatrixOperator< SC, MATSC, NT::device_type >;
Declaration of Tpetra::Details::Behavior, a class that describes Tpetra's behavior.
Abstract interface for local operators (e.g., matrices and preconditioners).
Namespace Tpetra contains the class and methods constituting the Tpetra library.