47 #ifndef THYRA_MUELU_REFMAXWELL_PRECONDITIONER_FACTORY_DECL_HPP
48 #define THYRA_MUELU_REFMAXWELL_PRECONDITIONER_FACTORY_DECL_HPP
52 #ifdef HAVE_MUELU_STRATIMIKOS
55 #include "Thyra_DefaultPreconditioner.hpp"
56 #include "Thyra_BlockedLinearOpBase.hpp"
57 #include "Thyra_DiagonalLinearOpBase.hpp"
59 #ifdef HAVE_MUELU_TPETRA
60 #include "Thyra_TpetraLinearOp.hpp"
61 #include "Thyra_TpetraThyraWrappers.hpp"
63 #ifdef HAVE_MUELU_EPETRA
64 #include "Thyra_EpetraLinearOp.hpp"
65 #include "Thyra_EpetraThyraWrappers.hpp"
68 #include "Teuchos_Ptr.hpp"
70 #include "Teuchos_Assert.hpp"
78 #include <MueLu_Hierarchy.hpp>
80 #include <MueLu_HierarchyUtils.hpp>
81 #include <MueLu_Utilities.hpp>
82 #include <MueLu_ParameterListInterpreter.hpp>
83 #include <MueLu_MLParameterListInterpreter.hpp>
86 #include <MueLu_RefMaxwell.hpp>
87 #ifdef HAVE_MUELU_TPETRA
88 #include <MueLu_TpetraOperator.hpp>
90 #ifdef HAVE_MUELU_EPETRA
95 #include "Thyra_PreconditionerFactoryBase.hpp"
110 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node = KokkosClassic::DefaultNode::DefaultNodeType>
125 bool isCompatible(
const LinearOpSourceBase<Scalar>& fwdOp)
const;
130 PreconditionerBase<Scalar>* prec,
131 const ESupportSolveUse supportSolveUse
136 ESupportSolveUse* supportSolveUse
172 #ifdef HAVE_MUELU_EPETRA
203 #ifdef HAVE_MUELU_TPETRA
204 if (Xpetra::ThyraUtils<Scalar,LocalOrdinal,GlobalOrdinal,Node>::isTpetra(fwdOp))
return true;
207 #ifdef HAVE_MUELU_EPETRA
208 if (Xpetra::ThyraUtils<Scalar,LocalOrdinal,GlobalOrdinal,Node>::isEpetra(fwdOp))
return true;
211 if (Xpetra::ThyraUtils<Scalar,LocalOrdinal,GlobalOrdinal,Node>::isBlockedOperator(fwdOp))
return true;
223 PreconditionerBase<Scalar>* prec,
224 const ESupportSolveUse
226 using Teuchos::rcp_dynamic_cast;
231 typedef Xpetra::ThyraUtils<Scalar,LocalOrdinal,GlobalOrdinal,Node> XpThyUtils;
238 typedef Thyra::LinearOpBase<Scalar> ThyLinOpBase;
239 typedef Thyra::DiagonalLinearOpBase<Scalar> ThyDiagLinOpBase;
240 #ifdef HAVE_MUELU_TPETRA
242 #if ((defined(EPETRA_HAVE_OMP) && (defined(HAVE_TPETRA_INST_OPENMP) && defined(HAVE_TPETRA_INST_INT_INT))) || \
243 (!defined(EPETRA_HAVE_OMP) && (defined(HAVE_TPETRA_INST_SERIAL) && defined(HAVE_TPETRA_INST_INT_INT))))
244 typedef Thyra::TpetraLinearOp<Scalar,LocalOrdinal,GlobalOrdinal,Node> ThyTpLinOp;
247 #if defined(HAVE_MUELU_EPETRA)
248 typedef Thyra::EpetraLinearOp ThyEpLinOp;
269 bool bIsEpetra = XpThyUtils::isEpetra(fwdOp);
270 bool bIsTpetra = XpThyUtils::isTpetra(fwdOp);
271 bool bIsBlocked = XpThyUtils::isBlockedOperator(fwdOp);
279 Teuchos::rcp_dynamic_cast<
const Thyra::BlockedLinearOpBase<Scalar> >(fwdOp);
282 TEUCHOS_TEST_FOR_EXCEPT(ThyBlockedOp->blockExists(0,0)==
false);
291 RCP<XpCrsMat> xpetraFwdCrsMatNonConst00 = Teuchos::rcp_const_cast<XpCrsMat>(xpetraFwdCrsMat00);
312 RCP<XpCrsMat> xpetraFwdCrsMatNonConst = Teuchos::rcp_const_cast<XpCrsMat>(xpetraFwdCrsMat);
326 thyra_precOp = rcp_dynamic_cast<Thyra::LinearOpBase<Scalar> >(defaultPrec->getNonconstUnspecifiedPrecOp(),
true);
334 const bool startingOver = (thyra_precOp.
is_null() || !paramList.
isParameter(
"reuse: type") || paramList.
get<std::string>(
"reuse: type") ==
"none");
336 if (startingOver ==
true) {
343 #ifdef HAVE_MUELU_TPETRA
345 #if ((defined(EPETRA_HAVE_OMP) && (defined(HAVE_TPETRA_INST_OPENMP) && defined(HAVE_TPETRA_INST_INT_INT))) || \
346 (!defined(EPETRA_HAVE_OMP) && (defined(HAVE_TPETRA_INST_SERIAL) && defined(HAVE_TPETRA_INST_INT_INT))))
347 typedef Tpetra::Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node> tV;
348 typedef Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> tMV;
349 typedef Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> TpCrsMat;
352 paramList.
remove(
"Nullspace");
353 RCP<XpMultVec> nullspace = MueLu::TpetraMultiVector_To_XpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>(tpetra_nullspace);
370 RCP<XpCrsMat> crsM1NonConst = Teuchos::rcp_const_cast<XpCrsMat>(crsM1);
394 RCP<XpCrsMat> crsD0NonConst = Teuchos::rcp_const_cast<XpCrsMat>(crsD0);
409 paramList.
remove(
"M0inv");
410 RCP<XpCrsMat> xM0inv = rcp_dynamic_cast<XpCrsMat>(tM0inv,
true);
414 paramList.
remove(
"M0inv");
417 RCP<const tV> tDiag = Thyra::TpetraOperatorVectorExtraction<Scalar,LocalOrdinal,GlobalOrdinal,Node>::getConstTpetraVector(diag);
422 paramList.
remove(
"M0inv");
426 RCP<XpCrsMat> crsM0invNonConst = Teuchos::rcp_const_cast<XpCrsMat>(crsM0inv);
439 "Thyra::MueLuRefMaxwellPreconditionerFactory: Tpetra does not support GO=int and or EpetraNode.");
443 #ifdef HAVE_MUELU_EPETRA
448 paramList.
remove(
"Nullspace");
451 RCP<XpMultVec> nullspace = rcp_dynamic_cast<XpMultVec>(xpEpNullspaceMult,
true);
460 RCP<XpCrsMat> xCrsM1 = rcp_dynamic_cast<XpCrsMat>(xeM1,
true);
462 RCP<XpMat> xM1 = rcp_dynamic_cast<XpMat>(xwM1);
471 RCP<XpCrsMat> crsM1NonConst = Teuchos::rcp_const_cast<XpCrsMat>(crsM1);
488 RCP<XpCrsMat> xCrsD0 = rcp_dynamic_cast<XpCrsMat>(xeD0,
true);
490 RCP<XpMat> xD0 = rcp_dynamic_cast<XpMat>(xwD0);
499 RCP<XpCrsMat> crsD0NonConst = Teuchos::rcp_const_cast<XpCrsMat>(crsD0);
514 paramList.
remove(
"M0inv");
516 RCP<XpCrsMat> xCrsM0inv = rcp_dynamic_cast<XpCrsMat>(xeM0inv,
true);
518 RCP<XpMat> xM0inv = rcp_dynamic_cast<XpMat>(xwM0inv);
523 paramList.
remove(
"M0inv");
537 paramList.
remove(
"M0inv");
541 RCP<XpCrsMat> crsM0invNonConst = Teuchos::rcp_const_cast<XpCrsMat>(crsM0inv);
555 paramList.
set<
bool>(
"refmaxwell: use as preconditioner",
true);
562 #if defined(HAVE_MUELU_TPETRA)
564 #if ((defined(EPETRA_HAVE_OMP) && (defined(HAVE_TPETRA_INST_OPENMP) && defined(HAVE_TPETRA_INST_INT_INT))) || \
565 (!defined(EPETRA_HAVE_OMP) && (defined(HAVE_TPETRA_INST_SERIAL) && defined(HAVE_TPETRA_INST_INT_INT))))
566 RCP<ThyTpLinOp> tpetr_precOp = rcp_dynamic_cast<ThyTpLinOp>(thyra_precOp);
570 "Thyra::MueLuRefMaxwellPreconditionerFactory: Tpetra does not support GO=int and or EpetraNode.");
574 #if defined(HAVE_MUELU_EPETRA)// && defined(HAVE_MUELU_SERIAL)
576 RCP<ThyEpLinOp> epetr_precOp = rcp_dynamic_cast<ThyEpLinOp>(thyra_precOp);
589 RCP<XpOp> xpOp = Teuchos::rcp_dynamic_cast<XpOp>(preconditioner);
590 thyraPrecOp = Thyra::xpetraLinearOp<Scalar, LocalOrdinal, GlobalOrdinal, Node>(thyraRangeSpace, thyraDomainSpace,xpOp);
594 defaultPrec->initializeUnspecified(thyraPrecOp);
600 ESupportSolveUse* supportSolveUse
610 *fwdOp = Teuchos::null;
613 if (supportSolveUse) {
615 *supportSolveUse = Thyra::SUPPORT_SOLVE_UNSPECIFIED;
618 defaultPrec->uninitialize();
635 return savedParamList;
656 std::string
description()
const {
return "Thyra::MueLuRefMaxwellPreconditionerFactory"; }
666 #endif // HAVE_MUELU_EPETRA
670 #endif // #ifdef HAVE_MUELU_STRATIMIKOS
672 #endif // THYRA_MUELU_REFMAXWELL_PRECONDITIONER_FACTORY_DECL_HPP
void initializePrec(const Teuchos::RCP< const LinearOpSourceBase< Scalar > > &fwdOpSrc, PreconditionerBase< Scalar > *prec, const ESupportSolveUse) const
Teuchos::RCP< Teuchos::ParameterList > unsetParameterList()
T & get(const std::string &name, T def_value)
static RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::magnitudeType, LocalOrdinal, GlobalOrdinal, Node > > ExtractCoordinatesFromParameterList(ParameterList ¶mList)
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
bool nonnull(const std::shared_ptr< T > &p)
MueLuRefMaxwellPreconditionerFactory()
bool is_null(const std::shared_ptr< T > &p)
bool isCompatible(const LinearOpSourceBase< Scalar > &fwdOp) const
Concrete preconditioner factory subclass for Thyra based on MueLu.Add support for MueLu preconditione...
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
const Epetra_CrsGraph & toEpetra(const RCP< const CrsGraph< int, GlobalOrdinal, Node > > &graph)
Teuchos::RCP< Teuchos::ParameterList > paramList_
Teuchos::RCP< PreconditionerBase< Scalar > > createPrec() const
Teuchos::RCP< Teuchos::ParameterList > getNonconstParameterList()
Teuchos::RCP< Teuchos::ParameterList > getNonconstParameterList()
void uninitializePrec(PreconditionerBase< Scalar > *prec, Teuchos::RCP< const LinearOpSourceBase< Scalar > > *fwdOp, ESupportSolveUse *supportSolveUse) const
void uninitializePrec(PreconditionerBase< Scalar > *prec, Teuchos::RCP< const LinearOpSourceBase< Scalar > > *fwdOp, ESupportSolveUse *supportSolveUse) const
Teuchos::RCP< PreconditionerBase< Scalar > > createPrec() const
Preconditioner (wrapped as a Xpetra::Operator) for Maxwell's equations in curl-curl form...
Teuchos::RCP< const Teuchos::ParameterList > getParameterList() const
bool isParameter(const std::string &name) const
bool remove(std::string const &name, bool throwIfNotExists=true)
MueLuRefMaxwellPreconditionerFactory()
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
void initializePrec(const Teuchos::RCP< const LinearOpSourceBase< Scalar > > &fwdOp, PreconditionerBase< Scalar > *prec, const ESupportSolveUse supportSolveUse) const
Teuchos::RCP< Teuchos::ParameterList > paramList_
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > ¶mList)
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > ¶mList)
RCP< const CrsGraph< int, GlobalOrdinal, Node > > toXpetra(const Epetra_CrsGraph &g)
std::string description() const
std::string description() const
bool isType(const std::string &name) const
bool isCompatible(const LinearOpSourceBase< Scalar > &fwdOpSrc) const
Teuchos::RCP< Teuchos::ParameterList > unsetParameterList()
static RCP< Matrix > Build(const RCP< const Map > &rowMap, size_t maxNumEntriesPerRow, Xpetra::ProfileType pftype=Xpetra::DynamicProfile)
Exception throws to report errors in the internal logical of the program.
#define TEUCHOS_ASSERT(assertion_test)
#define TEUCHOS_TEST_FOR_EXCEPT(throw_exception_test)
Teuchos::RCP< const Teuchos::ParameterList > getParameterList() const