47 #ifndef MUELU_TPETRAOPERATOR_DEF_HPP
48 #define MUELU_TPETRAOPERATOR_DEF_HPP
52 #include <Xpetra_BlockedMap.hpp>
54 #include <Xpetra_CrsMatrixWrap.hpp>
57 #include <Xpetra_TpetraMultiVector.hpp>
60 #include "MueLu_Hierarchy.hpp"
61 #include "MueLu_Utilities.hpp"
65 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
73 if (!Hierarchy_.is_null())
74 domainMap = Hierarchy_->GetLevel(0)->template Get<RCP<Matrix> >(
"A")->getDomainMap();
76 domainMap = Operator_->getDomainMap();
79 if (bDomainMap.
is_null() ==
false) {
85 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
92 if (!Hierarchy_.is_null())
93 rangeMap = Hierarchy_->GetLevel(0)->template Get<RCP<Matrix> >(
"A")->getRangeMap();
95 rangeMap = Operator_->getRangeMap();
98 if (bRangeMap.
is_null() ==
false) {
104 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
114 TMV& temp_x =
const_cast<TMV&
>(X);
115 const XTMV tX(rcpFromRef(temp_x));
116 XTMV tY(rcpFromRef(Y));
118 if (!Hierarchy_.is_null())
119 Hierarchy_->Iterate(tX, tY, 1,
true);
121 Operator_->apply(tX, tY);
123 }
catch (std::exception& e) {
124 std::cerr <<
"MueLu::TpetraOperator::apply : detected an exception" << std::endl
125 << e.what() << std::endl;
130 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
135 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
141 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
149 #endif // ifdef MUELU_TPETRAOPERATOR_DEF_HPP
bool hasTransposeApply() const
Indicates whether this operator supports applying the adjoint operator.
const RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > toTpetraNonZero(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Teuchos::RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const
Returns the Tpetra::Map object associated with the range of this operator.
MueLu::DefaultScalar Scalar
RCP< Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > > GetOperator() const
Direct access to the underlying MueLu::Operator.
Teuchos::RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const
Returns the Tpetra::Map object associated with the domain of this operator.
RCP< MueLu::Hierarchy< Scalar, LocalOrdinal, GlobalOrdinal, Node > > GetHierarchy() const
Direct access to the underlying MueLu::Hierarchy.
void apply(const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=Teuchos::ScalarTraits< Scalar >::one(), Scalar beta=Teuchos::ScalarTraits< Scalar >::one()) const
Returns in Y the result of a Tpetra::Operator applied to a Tpetra::MultiVector X. ...