42 #ifndef THYRA_BELOS_TPETRA_PRECONDITIONERFACTORY_DEF_HPP
43 #define THYRA_BELOS_TPETRA_PRECONDITIONERFACTORY_DEF_HPP
47 #include "Thyra_DefaultPreconditioner.hpp"
51 #include "BelosTpetraOperator.hpp"
52 #include "BelosConfigDefs.hpp"
53 #include "BelosLinearProblem.hpp"
54 #include "BelosTpetraAdapter.hpp"
56 #include "Tpetra_MixedScalarMultiplyOp.hpp"
72 #if (!defined(HAVE_TPETRA_INST_DOUBLE) || (defined(HAVE_TPETRA_INST_DOUBLE) && defined(HAVE_TPETRA_INST_FLOAT))) && \
73 (!defined(HAVE_TPETRA_INST_COMPLEX_DOUBLE) || (defined(HAVE_TPETRA_INST_COMPLEX_DOUBLE) && defined(HAVE_TPETRA_INST_COMPLEX_FLOAT)))
74 # define THYRA_BELOS_PREC_ENABLE_HALF_PRECISION
82 template <
typename MatrixType>
90 template <
typename MatrixType>
92 const LinearOpSourceBase<scalar_type> &fwdOpSrc
95 typedef typename MatrixType::local_ordinal_type local_ordinal_type;
96 typedef typename MatrixType::global_ordinal_type global_ordinal_type;
97 typedef typename MatrixType::node_type node_type;
101 const auto tpetraFwdOp = TpetraExtractHelper::getConstTpetraOperator(fwdOp);
109 template <
typename MatrixType>
113 return Teuchos::rcp(
new DefaultPreconditioner<scalar_type>);
117 template <
typename MatrixType>
119 const Teuchos::RCP<
const LinearOpSourceBase<scalar_type> > &fwdOpSrc,
120 PreconditionerBase<scalar_type> *prec,
121 const ESupportSolveUse
134 totalTimer.start(
true);
140 *out <<
"\nEntering Thyra::BelosTpetraPreconditionerFactory::initializePrec(...) ...\n";
148 typedef typename MatrixType::local_ordinal_type local_ordinal_type;
149 typedef typename MatrixType::global_ordinal_type global_ordinal_type;
150 typedef typename MatrixType::node_type node_type;
152 typedef Tpetra::Operator<scalar_type, local_ordinal_type, global_ordinal_type, node_type> TpetraLinOp;
154 const auto tpetraFwdOp = TpetraExtractHelper::getConstTpetraOperator(fwdOp);
158 typedef Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type> TpetraMV;
159 typedef Belos::TpetraOperator<scalar_type, local_ordinal_type, global_ordinal_type, node_type> BelosTpOp;
160 typedef Belos::LinearProblem<scalar_type, TpetraMV, TpetraLinOp> BelosTpLinProb;
162 #ifdef THYRA_BELOS_PREC_ENABLE_HALF_PRECISION
167 typedef Tpetra::Operator<half_scalar_type, local_ordinal_type, global_ordinal_type, node_type> TpetraLinOpHalf;
168 typedef Tpetra::MultiVector<half_scalar_type, local_ordinal_type, global_ordinal_type, node_type> TpetraMVHalf;
169 typedef Belos::TpetraOperator<half_scalar_type, local_ordinal_type, global_ordinal_type, node_type> BelosTpOpHalf;
170 typedef Belos::LinearProblem<half_scalar_type, TpetraMVHalf, TpetraLinOpHalf> BelosTpLinProbHalf;
176 Teuchos::ptr(
dynamic_cast<DefaultPreconditioner<scalar_type> *
>(prec));
181 if (paramList_.is_null ()) {
185 innerParamList = paramList_;
188 bool useHalfPrecision =
false;
190 useHalfPrecision = Teuchos::getParameter<bool>(*innerParamList,
"half precision");
192 const std::string solverType = Teuchos::getParameter<std::string>(*innerParamList,
"BelosPrec Solver Type");
196 std::string solverTypeUpper (solverType);
197 std::transform(solverTypeUpper.begin(), solverTypeUpper.end(),solverTypeUpper.begin(), ::toupper);
201 *out <<
"\nCreating a new BelosTpetra::Preconditioner object...\n";
205 if (useHalfPrecision) {
206 #ifdef THYRA_BELOS_PREC_ENABLE_HALF_PRECISION
208 Teuchos::OSTab(out).o() <<
"> Creating half precision preconditioner\n";
210 const RCP<const MatrixType> tpetraFwdMatrix = Teuchos::rcp_dynamic_cast<
const MatrixType>(tpetraFwdOp);
212 auto tpetraFwdMatrixHalf = tpetraFwdMatrix->template convert<half_scalar_type>();
215 belosLinProbHalf->setOperator(tpetraFwdMatrixHalf);
216 RCP<TpetraLinOpHalf> belosOpRCPHalf =
rcp(
new BelosTpOpHalf(belosLinProbHalf, packageParamList, solverType,
true));
217 RCP<TpetraLinOp> wrappedOp =
rcp(
new Tpetra::MixedScalarMultiplyOp<scalar_type,half_scalar_type,local_ordinal_type,global_ordinal_type,node_type>(belosOpRCPHalf));
226 belosLinProb->setOperator(tpetraFwdOp);
227 RCP<TpetraLinOp> belosOpRCP =
rcp(
new BelosTpOp(belosLinProb, packageParamList, solverType,
true));
230 defaultPrec->initializeUnspecified(thyraPrecOp);
234 *out <<
"\nTotal time in Thyra::BelosTpetraPreconditionerFactory::initializePrec(...) = " << totalTimer.totalElapsedTime() <<
" sec\n";
238 *out <<
"\nLeaving Thyra::BelosTpetraPreconditionerFactory::initializePrec(...) ...\n";
243 template <
typename MatrixType>
245 PreconditionerBase<scalar_type> *prec,
246 Teuchos::RCP<
const LinearOpSourceBase<scalar_type> > *fwdOp,
247 ESupportSolveUse *supportSolveUse
257 Teuchos::ptr(
dynamic_cast<DefaultPreconditioner<scalar_type> *
>(prec));
265 if (supportSolveUse) {
267 *supportSolveUse = Thyra::SUPPORT_SOLVE_UNSPECIFIED;
270 defaultPrec->uninitialize();
277 template <
typename MatrixType>
284 const auto validParamList = this->getValidParameters();
288 paramList_ = paramList;
289 Teuchos::readVerboseObjectSublist(paramList_.getRawPtr(),
this);
293 template <
typename MatrixType>
301 template <
typename MatrixType>
307 return savedParamList;
311 template <
typename MatrixType>
318 template <
typename MatrixType>
327 "BelosPrec Solver Type",
"GMRES",
328 "Name of Belos solver to be used as a preconditioner. (Use valid names for Belos::SolverFactory.)"
331 "half precision",
false,
332 "Whether a half-of-standard-precision Belos-solver-as-preconditioner should be built."
335 "BelosPrec Solver Params",
false,
336 "Belos solver settings that are passed onto the Belos solver itself."
338 Teuchos::setupVerboseObjectSublist(validParamList.
getRawPtr());
341 return validParamList;
347 template <
typename MatrixType>
350 return "Thyra::BelosTpetraPreconditionerFactory";
356 #endif // THYRA_BELOS_TPETRA_PRECONDITIONERFACTORY_DEF_HPP
std::string description() const
void uninitializePrec(PreconditionerBase< scalar_type > *prec, Teuchos::RCP< const LinearOpSourceBase< scalar_type > > *fwdOp, ESupportSolveUse *supportSolveUse) const
bool isCompatible(const LinearOpSourceBase< scalar_type > &fwdOp) const
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)
bool is_null(const std::shared_ptr< T > &p)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
void initializePrec(const Teuchos::RCP< const LinearOpSourceBase< scalar_type > > &fwdOp, PreconditionerBase< scalar_type > *prec, const ESupportSolveUse supportSolveUse) const
RCP< ParameterList > sublist(const RCP< ParameterList > ¶mList, const std::string &name, bool mustAlreadyExist=false, const std::string &docString="")
basic_OSTab< char > OSTab
RCP< LinearOpBase< Scalar > > createLinearOp(const RCP< Tpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraOperator, const RCP< const VectorSpaceBase< Scalar > > rangeSpace=Teuchos::null, const RCP< const VectorSpaceBase< Scalar > > domainSpace=Teuchos::null)
bool isParameter(const std::string &name) const
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > ¶mList)
void validateParametersAndSetDefaults(ParameterList const &validParamList, int const depth=1000)
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
TEUCHOSCORE_LIB_DLL_EXPORT bool includesVerbLevel(const EVerbosityLevel verbLevel, const EVerbosityLevel requestedVerbLevel, const bool isDefaultLevel=false)
void validateParameters(ParameterList const &validParamList, int const depth=1000, EValidateUsed const validateUsed=VALIDATE_USED_ENABLED, EValidateDefaults const validateDefaults=VALIDATE_DEFAULTS_ENABLED) const
Teuchos::RCP< const Teuchos::ParameterList > getParameterList() const
Teuchos::RCP< Teuchos::ParameterList > getNonconstParameterList()
Teuchos::RCP< PreconditionerBase< scalar_type > > createPrec() const
BelosTpetraPreconditionerFactory()
#define TEUCHOS_ASSERT(assertion_test)
#define TEUCHOS_TEST_FOR_EXCEPT(throw_exception_test)
Teuchos::RCP< Teuchos::ParameterList > unsetParameterList()