42#ifndef TPETRA_CRSMATRIXMULTIPLYOP_HPP
43#define TPETRA_CRSMATRIXMULTIPLYOP_HPP
51#include "Tpetra_CrsMatrix.hpp"
55#include "Tpetra_LocalCrsMatrixOperator.hpp"
95 template <
class Scalar,
101 public Operator<Scalar, LocalOrdinal, GlobalOrdinal, Node>
111 using local_matrix_device_type =
112 typename crs_matrix_type::local_matrix_device_type;
125 A->getLocalMatrixDevice ()))
140 Teuchos::ETransp
mode = Teuchos::NO_TRANS,
141 Scalar alpha = Teuchos::ScalarTraits<Scalar>::one (),
142 Scalar beta = Teuchos::ScalarTraits<Scalar>::zero ())
const override
145 (!
matrix_->isFillComplete (), std::runtime_error,
146 Teuchos::typeName (*
this) <<
"::apply(): underlying matrix is not fill-complete.");
148 (
X.getNumVectors () !=
Y.getNumVectors (), std::runtime_error,
149 Teuchos::typeName (*
this) <<
"::apply(X,Y): X and Y must have the same number of vectors.");
151 (Teuchos::ScalarTraits<Scalar>::isComplex &&
mode == Teuchos::TRANS, std::logic_error,
152 Teuchos::typeName (*
this) <<
"::apply() does not currently support transposed multiplications for complex scalar types.");
153 if (
mode == Teuchos::NO_TRANS) {
172 return matrix_->getDomainMap ();
177 return matrix_->getRangeMap ();
186 const Teuchos::RCP<const crs_matrix_type>
matrix_;
204 mutable Teuchos::RCP<MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
importMV_;
218 mutable Teuchos::RCP<MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
exportMV_;
225 Teuchos::ETransp
mode,
234 using STS = Teuchos::ScalarTraits<Scalar>;
257 X = Teuchos::rcp (
new MV (
X_in, Teuchos::Copy));
261 X = Teuchos::rcpFromRef (
X_in);
288 auto X_lcl =
X->getLocalViewDevice(Access::ReadOnly);
298 Y_in.putScalar (STS::zero ());
308 if (!
Y_in.isConstantStride () ||
X.getRawPtr () == &
Y_in) {
313 else Y_lcl =
Y.getLocalViewDevice(Access::ReadWrite);
321 else Y_lcl =
Y_in.getLocalViewDevice(Access::ReadWrite);
340 using Teuchos::NO_TRANS;
343 using Teuchos::rcp_const_cast;
344 using Teuchos::rcpFromRef;
347 typedef Teuchos::ScalarTraits<Scalar> STS;
350 if (
alpha == STS::zero ()) {
351 if (
beta == STS::zero ()) {
352 Y_in.putScalar (STS::zero ());
353 }
else if (
beta != STS::one ()) {
379 (!
Y_in.isDistributed () &&
matrix_->getComm ()->getSize () != 1);
397 if (!
X_in.isConstantStride ()) {
415 ProfilingRegion
regionImport (
"Tpetra::CrsMatrixMultiplyOp::apply: Import");
442 auto Y_lcl =
Y_rowMap->getLocalViewDevice(Access::OverwriteAll);
445 alpha, STS::zero ());
447 ProfilingRegion
regionExport (
"Tpetra::CrsMatrixMultiplyOp::apply: Export");
454 Y_in.putScalar (STS::zero());
484 if (
beta != STS::zero ()) {
489 else Y_lcl =
Y_rowMap->getLocalViewDevice(Access::ReadWrite);
497 else Y_lcl =
Y_in.getLocalViewDevice(Access::ReadWrite);
508 ProfilingRegion
regionReduce (
"Tpetra::CrsMatrixMultiplyOp::apply: Reduce Y");
535 const bool force =
false)
const
598 getRowMapMultiVector (
const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Y_rangeMap,
599 const bool force =
false)
const
604 typedef Export<LocalOrdinal,GlobalOrdinal,Node> export_type;
606 const size_t numVecs = Y_rangeMap.getNumVectors ();
607 RCP<const export_type> exporter =
matrix_->getGraph ()->getExporter ();
608 RCP<const map_type> rowMap =
matrix_->getRowMap ();
620 if (! exporter.is_null () || force) {
622 Y_rowMap = rcp (
new MV (rowMap, numVecs));
640 template <
class OpScalar,
646 CrsMatrixMultiplyOp<OpScalar, MatScalar, LocalOrdinal, GlobalOrdinal, Node> >
652 return Teuchos::rcp (
new op_type (
A));
Forward declaration of Tpetra::CrsMatrixMultiplyOp.
Declaration of Tpetra::Details::Behavior, a class that describes Tpetra's behavior.
Declaration of Tpetra::Details::Profiling, a scope guard for Kokkos Profiling.
Stand-alone utility functions and macros.
A class for wrapping a CrsMatrix multiply in a Operator.
Teuchos::RCP< MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > importMV_
Column Map MultiVector used in apply().
bool hasTransposeApply() const override
Whether this Operator's apply() method can apply the transpose or conjugate transpose.
void applyNonTranspose(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X_in, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y_in, Scalar alpha, Scalar beta) const
Apply the matrix (not its transpose) to X_in, producing Y_in.
void apply(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=Teuchos::ScalarTraits< Scalar >::one(), Scalar beta=Teuchos::ScalarTraits< Scalar >::zero()) const override
Compute Y = beta*Y + alpha*Op(A)*X, where Op(A) is either A, , or .
Teuchos::RCP< const map_type > getRangeMap() const override
The range Map of this Operator.
void applyTranspose(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X_in, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y_in, Teuchos::ETransp mode, Scalar alpha, Scalar beta) const
Apply the transpose or conjugate transpose of the matrix to X_in, producing Y_in.
CrsMatrixMultiplyOp(const Teuchos::RCP< const crs_matrix_type > &A)
Constructor.
~CrsMatrixMultiplyOp() override=default
Destructor (virtual for memory safety of derived classes).
Teuchos::RCP< CrsMatrixMultiplyOp< OpScalar, MatScalar, LocalOrdinal, GlobalOrdinal, Node > > createCrsMatrixMultiplyOp(const Teuchos::RCP< const CrsMatrix< MatScalar, LocalOrdinal, GlobalOrdinal, Node > > &A)
Non-member function to create a CrsMatrixMultiplyOp.
const Teuchos::RCP< const crs_matrix_type > matrix_
The underlying CrsMatrix object.
Teuchos::RCP< MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > exportMV_
Row Map MultiVector used in apply().
LocalCrsMatrixOperator< Scalar, MatScalar, typename crs_matrix_type::device_type > localMultiply_
Implementation of local sparse matrix-vector multiply.
Teuchos::RCP< const map_type > getDomainMap() const override
The domain Map of this Operator.
Struct that holds views of the contents of a CrsMatrix.
typename Node::device_type device_type
The Kokkos device type.
Abstract interface for local operators (e.g., matrices and preconditioners).
Abstract interface for operators (e.g., matrices and preconditioners).
Namespace Tpetra contains the class and methods constituting the Tpetra library.
void deep_copy(MultiVector< DS, DL, DG, DN > &dst, const MultiVector< SS, SL, SG, SN > &src)
Copy the contents of the MultiVector src into dst.
@ ADD_ASSIGN
Accumulate new values into existing values (may not be supported in all classes)
@ INSERT
Insert new values that don't currently exist.