42 #ifndef Thyra_TpetraExplicitAdjointModelEvaluator_hpp
43 #define Thyra_TpetraExplicitAdjointModelEvaluator_hpp
45 #include "Thyra_ModelEvaluatorDelegatorBase.hpp"
47 #include "Tpetra_RowMatrixTransposer.hpp"
61 template <
typename Scalar,
typename LocalOrdinal,
typename GlobalOrdinal=LocalOrdinal,
typename Node=Tpetra::KokkosClassic::DefaultNode::DefaultNodeType>
63 public ModelEvaluatorDelegatorBase<Scalar>{
68 const RCP<
const ModelEvaluator<Scalar> >& model) :
69 ModelEvaluatorDelegatorBase<Scalar>(model) {}
73 const RCP<ModelEvaluator<Scalar> >& model) :
74 ModelEvaluatorDelegatorBase<Scalar>(model) {}
86 ModelEvaluatorBase::InArgsSetup<Scalar> inArgs =
87 this->getUnderlyingModel()->createInArgs();
88 inArgs.setModelEvalDescription(this->description());
94 typedef Tpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node> TO;
95 typedef Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> TCM;
96 typedef Tpetra::RowMatrixTransposer<Scalar,LocalOrdinal,GlobalOrdinal,Node> TRMT;
98 thyra_fwd_op = this->getUnderlyingModel()->create_W_op();
101 RCP<TO> tpetra_fwd_op = thyra_tpetra_fwd_op->getTpetraOperator();
103 Teuchos::rcp_dynamic_cast<TCM>(tpetra_fwd_op,
true);
104 TRMT transposer(tpetra_fwd_mat);
105 RCP<TCM> tpetra_trans_mat = transposer.createTranspose();
106 return tpetraLinearOp(thyra_op->range(), thyra_op->domain(),
116 typedef ModelEvaluatorBase MEB;
117 MEB::OutArgs<Scalar> model_outArgs =
118 this->getUnderlyingModel()->createOutArgs();
119 MEB::OutArgsSetup<Scalar> outArgs;
120 outArgs.setModelEvalDescription(this->description());
121 outArgs.setSupports(MEB::OUT_ARG_f);
122 outArgs.setSupports(MEB::OUT_ARG_W_op);
128 const ModelEvaluatorBase::InArgs<Scalar> &inArgs,
129 const ModelEvaluatorBase::OutArgs<Scalar> &outArgs)
const
131 typedef ModelEvaluatorBase MEB;
133 typedef Tpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node> TO;
134 typedef Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> TCM;
135 typedef Tpetra::RowMatrixTransposer<Scalar,LocalOrdinal,GlobalOrdinal,Node> TRMT;
137 MEB::OutArgs<Scalar> model_outArgs =
138 this->getUnderlyingModel()->createOutArgs();
140 if (model_outArgs.supports(MEB::OUT_ARG_W_op) &&
144 thyra_fwd_op = this->getUnderlyingModel()->create_W_op();
146 this->getUnderlyingModel()->evalModel(inArgs, model_outArgs);
153 RCP<TO> tpetra_fwd_op = thyra_tpetra_fwd_op->getTpetraOperator();
155 Teuchos::rcp_dynamic_cast<TCM>(tpetra_fwd_op,
true);
156 TRMT transposer(tpetra_fwd_mat);
157 RCP<TCM> tpetra_trans_mat = transposer.createTranspose();
160 RCP<LOB> thyra_adj_op = outArgs.get_W_op();
162 Teuchos::rcp_dynamic_cast<TLO>(thyra_adj_op,
true);
163 RCP<TO> tpetra_adj_op = thyra_tpetra_adj_op->getTpetraOperator();
165 Teuchos::rcp_dynamic_cast<TCM>(tpetra_adj_op,
true);
166 *tpetra_adj_mat = *tpetra_trans_mat;
172 template <
typename Scalar>
173 RCP<TpetraExplicitAdjointModelEvaluator<Scalar> >
175 const RCP<
const ModelEvaluator<Scalar> >& model)
180 template <
typename Scalar>
181 RCP<TpetraExplicitAdjointModelEvaluator<Scalar> >
183 const RCP<ModelEvaluator<Scalar> >& model)
void evalModelImpl(const ModelEvaluatorBase::InArgs< Scalar > &inArgs, const ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
RCP< Thyra::LinearOpBase< Scalar > > thyra_fwd_op
ModelEvaluatorBase::InArgs< Scalar > createInArgs() const
RCP< LinearOpBase< Scalar > > create_W_op() const
virtual ~TpetraExplicitAdjointModelEvaluator()=default
Destructor.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
A model evaluator decorator for computing an explicit adjoint.
ModelEvaluatorBase::OutArgs< Scalar > createOutArgsImpl() const
#define TEUCHOS_ASSERT(assertion_test)
Concrete Thyra::LinearOpBase subclass for Tpetra::Operator.
TpetraExplicitAdjointModelEvaluator(const RCP< const ModelEvaluator< Scalar > > &model)
Constructor.
TpetraExplicitAdjointModelEvaluator(const RCP< ModelEvaluator< Scalar > > &model)
Constructor.
RCP< TpetraExplicitAdjointModelEvaluator< Scalar > > tpetraExplicitAdjointModelEvaluator(const RCP< const ModelEvaluator< Scalar > > &model)