46 #ifndef BELOS_MUELU_ADAPTER_HPP 47 #define BELOS_MUELU_ADAPTER_HPP 49 #ifdef HAVE_MUELU_EPETRA 50 #include <BelosOperator.hpp> 53 #ifdef HAVE_MUELU_AMGX 54 #include "MueLu_AMGXOperator.hpp" 60 #include "MueLu_Hierarchy.hpp" 64 using Teuchos::rcpFromRef;
79 template <
class Scalar,
80 class LocalOrdinal = int,
81 class GlobalOrdinal = LocalOrdinal,
82 class Node = KokkosClassic::DefaultNode::DefaultNodeType>
84 public OperatorT<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
85 #ifdef HAVE_MUELU_TPETRA 86 ,
public OperatorT<Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
96 #ifdef HAVE_MUELU_AMGX 111 void Apply(
const Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>& x, Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>& y, ETrans trans = NOTRANS )
const {
114 "Belos::MueLuOp::Apply, transpose mode != NOTRANS not supported by MueLu preconditionners.");
119 #ifdef HAVE_MUELU_AMGX 120 if (!AMGX_.is_null()) {
121 Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> tX = Xpetra::toTpetra(x);
122 Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> tY = Xpetra::toTpetra(y);
124 AMGX_->apply(tX, tY);
128 if (!Hierarchy_.is_null())
129 Hierarchy_->Iterate(x, y, 1,
true);
133 #ifdef HAVE_MUELU_TPETRA 140 void Apply(
const Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>& x, Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>& y, ETrans trans = NOTRANS )
const {
143 "Belos::MueLuOp::Apply, transpose mode != NOTRANS not supported by MueLu preconditionners.");
148 #ifdef HAVE_MUELU_AMGX 149 if (!AMGX_.is_null())
153 if (!Hierarchy_.is_null()) {
154 Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> & temp_x =
const_cast<Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> &
>(x);
156 const Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> tX(rcpFromRef(temp_x));
157 Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> tY(rcpFromRef(y));
158 Hierarchy_->Iterate(tX, tY, 1,
true);
164 RCP<MueLu::Hierarchy<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
Hierarchy_;
165 #ifdef HAVE_MUELU_AMGX 166 RCP<MueLu::AMGXOperator<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
AMGX_;
172 public OperatorT<Xpetra::MultiVector<double, int, int, Xpetra::EpetraNode> >
173 #ifdef HAVE_MUELU_TPETRA 175 ,
public OperatorT<Tpetra::MultiVector<double, int, int, Xpetra::EpetraNode> >
177 #ifdef HAVE_MUELU_EPETRA 179 ,
public Belos::Operator<double>
190 #ifdef HAVE_MUELU_AMGX 195 void Apply(
const Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>& x, Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>& y, ETrans trans = NOTRANS )
const {
198 "Belos::MueLuOp::Apply, transpose mode != NOTRANS not supported by MueLu preconditionners.");
203 #ifdef HAVE_MUELU_AMGX 204 if (!AMGX_.is_null()) {
205 Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> tX = Xpetra::toTpetra(x);
206 Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> tY = Xpetra::toTpetra(y);
208 AMGX_->apply(tX, tY);
211 if (!Hierarchy_.is_null())
212 Hierarchy_->Iterate(x, y, 1,
true);
215 #ifdef HAVE_MUELU_TPETRA 216 void Apply (
const Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>& x, Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>& y, ETrans trans=NOTRANS )
const {
218 "Belos::MueLuOp::Apply, transpose mode != NOTRANS not supported by MueLu preconditionners.");
223 #ifdef HAVE_MUELU_AMGX 224 if (!AMGX_.is_null())
228 if (!Hierarchy_.is_null()) {
229 Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> & temp_x =
const_cast<Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> &
>(x);
231 const Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> tX(rcpFromRef(temp_x));
232 Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> tY(rcpFromRef(y));
236 Hierarchy_->Iterate(tX, tY, 1,
true);
241 #ifdef HAVE_MUELU_EPETRA 248 void Apply(
const Epetra_MultiVector& x, Epetra_MultiVector& y, ETrans trans = NOTRANS)
const {
250 "Belos::MueLuOp::Apply, transpose mode != NOTRANS not supported by MueLu preconditionners.");
252 Epetra_MultiVector& temp_x =
const_cast<Epetra_MultiVector&
>(x);
254 const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> tX(rcpFromRef(temp_x));
255 Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> tY(rcpFromRef(y));
260 Hierarchy_->Iterate(tX, tY, 1,
true);
268 void Apply(
const Belos::MultiVec<double>& x, Belos::MultiVec<double>& y, ETrans trans = NOTRANS )
const {
269 const Epetra_MultiVector* vec_x =
dynamic_cast<const Epetra_MultiVector*
>(&x);
270 Epetra_MultiVector* vec_y =
dynamic_cast<Epetra_MultiVector*
>(&y);
272 TEUCHOS_TEST_FOR_EXCEPTION(vec_x==NULL || vec_y==NULL,
MueLuOpFailure,
273 "Belos::MueLuOp::Apply, x and/or y cannot be dynamic cast to an Epetra_MultiVector.");
275 Apply(*vec_x, *vec_y, trans);
280 RCP<MueLu::Hierarchy<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
Hierarchy_;
281 #ifdef HAVE_MUELU_AMGX 282 RCP<MueLu::AMGXOperator<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
AMGX_;
288 #endif // BELOS_MUELU_ADAPTER_HPP MueLuOp(const RCP< MueLu::Hierarchy< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &H)
void Apply(const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x, Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &y, ETrans trans=NOTRANS) const
void Apply(const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x, Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &y, ETrans trans=NOTRANS) const
This routine takes the Tpetra::MultiVector x and applies the operator to it resulting in the Tpetra::...
MueLuOpFailure is thrown when a return value from an MueLu call on an Xpetra::Operator or MueLu::Hier...
RCP< MueLu::AMGXOperator< Scalar, LocalOrdinal, GlobalOrdinal, Node > > AMGX_
void Apply(const Epetra_MultiVector &x, Epetra_MultiVector &y, ETrans trans=NOTRANS) const
This routine takes the Tpetra::MultiVector x and applies the operator to it resulting in the Tpetra::...
void Apply(const Belos::MultiVec< double > &x, Belos::MultiVec< double > &y, ETrans trans=NOTRANS) const
This routine takes the Belos::MultiVec x and applies the operator to it resulting in the Belos::Multi...
void Apply(const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x, Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &y, ETrans trans=NOTRANS) const
This routine takes the Xpetra::MultiVector x and applies the operator to it resulting in the Xpetra::...
RCP< MueLu::AMGXOperator< Scalar, LocalOrdinal, GlobalOrdinal, Node > > AMGX_
MueLuOpFailure(const std::string &what_arg)
RCP< MueLu::Hierarchy< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Hierarchy_
RCP< MueLu::Hierarchy< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Hierarchy_
MueLuOp(const RCP< MueLu::Hierarchy< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &H)
Default constructor.
MueLuOp(const RCP< MueLu::AMGXOperator< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A)
MueLuOp(const RCP< MueLu::AMGXOperator< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A)
Kokkos::Compat::KokkosSerialWrapperNode EpetraNode
void Apply(const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x, Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &y, ETrans trans=NOTRANS) const
Provides methods to build a multigrid hierarchy and apply multigrid cycles.
Adapter for AmgX library from Nvidia.