43#ifndef PANZER_EXPLICIT_MODEL_EVALUATOR_IMPL_HPP
44#define PANZER_EXPLICIT_MODEL_EVALUATOR_IMPL_HPP
46#include "Thyra_DefaultDiagonalLinearOp.hpp"
47#include "Thyra_LinearOpWithSolveFactoryBase.hpp"
49#include "PanzerDiscFE_config.hpp"
50#ifdef PANZER_HAVE_EPETRA_STACK
51#include "Thyra_EpetraModelEvaluator.hpp"
52#include "Thyra_get_Epetra_Operator.hpp"
58template<
typename Scalar>
61 bool constantMassMatrix,
63 bool applyMassInverse)
64 :
Thyra::ModelEvaluatorDelegatorBase<Scalar>(model)
65 , constantMassMatrix_(constantMassMatrix)
66 , massLumping_(useLumpedMass)
69 using Teuchos::rcp_dynamic_cast;
74 panzerModel_ = rcp_dynamic_cast<panzer::ModelEvaluator<Scalar> >(model);
76#ifdef PANZER_HAVE_EPETRA_STACK
78 RCP<Thyra::EpetraModelEvaluator> epME = rcp_dynamic_cast<Thyra::EpetraModelEvaluator>(model);
79 if(epME!=Teuchos::null)
80 panzerEpetraModel_ = rcp_dynamic_cast<const panzer::ModelEvaluator_Epetra>(epME->getEpetraModel());
87template<
typename Scalar>
91 typedef Thyra::ModelEvaluatorBase MEB;
93 MEB::InArgs<Scalar> nomVals = createInArgs();
94 nomVals.setArgs(this->getUnderlyingModel()->getNominalValues(),
true);
99template<
typename Scalar>
103 return prototypeInArgs_;
106template<
typename Scalar>
110 return prototypeOutArgs_;
113template<
typename Scalar>
117 return Teuchos::rcp_dynamic_cast<panzer::ModelEvaluator<Scalar> >(this->getNonconstUnderlyingModel());
120template<
typename Scalar>
122evalModelImpl(
const Thyra::ModelEvaluatorBase::InArgs<Scalar> &inArgs,
123 const Thyra::ModelEvaluatorBase::OutArgs<Scalar> &outArgs)
const
125 typedef Thyra::ModelEvaluatorBase MEB;
127 RCP<const Thyra::ModelEvaluator<Scalar> > under_me = this->getUnderlyingModel();
129 MEB::InArgs<Scalar> under_inArgs = under_me->createInArgs();
130 under_inArgs.setArgs(inArgs);
133 under_inArgs.set_x_dot(inArgs.get_x_dot());
134 under_inArgs.set_alpha(0.0);
136 Teuchos::RCP<Thyra::VectorBase<Scalar> > f = outArgs.get_f();
138 MEB::OutArgs<Scalar> under_outArgs = under_me->createOutArgs();
139 under_outArgs.setArgs(outArgs);
140 if(f!=Teuchos::null) {
142 if(scrap_f_==Teuchos::null)
143 scrap_f_ = Thyra::createMember(*under_me->get_f_space());
145 Thyra::assign(scrap_f_.ptr(),0.0);
146 under_outArgs.set_f(scrap_f_);
149 under_me->evalModel(under_inArgs,under_outArgs);
152 if(invMassMatrix_==Teuchos::null || constantMassMatrix_==
false)
153 buildInverseMassMatrix(inArgs);
156 Thyra::Vt_S(scrap_f_.ptr(),-1.0);
159 if(f!=Teuchos::null && this->applyMassInverse_) {
160 Thyra::apply(*invMassMatrix_,Thyra::NOTRANS,*scrap_f_,f.ptr());
162 else if(f!=Teuchos::null){
163 Thyra::V_V(f.ptr(),*scrap_f_);
167template<
typename Scalar>
171 typedef Thyra::ModelEvaluatorBase MEB;
173 using Thyra::createMember;
175 RCP<const Thyra::ModelEvaluator<Scalar> > me = this->getUnderlyingModel();
178 mass_ = me->create_W_op();
181 if(zero_==Teuchos::null) {
182 zero_ = Thyra::createMember(*me->get_x_space());
183 Thyra::assign(zero_.ptr(),0.0);
188 MEB::InArgs<Scalar> inArgs_new = me->createInArgs();
189 inArgs_new.setArgs(inArgs);
190 inArgs_new.set_x_dot(inArgs.get_x_dot());
191 inArgs_new.set_alpha(1.0);
192 inArgs_new.set_beta(0.0);
197 if(panzerModel_!=Teuchos::null)
198 panzerModel_->setOneTimeDirichletBeta(1.0);
199#ifdef PANZER_HAVE_EPETRA_STACK
200 else if(panzerEpetraModel_!=Teuchos::null)
201 panzerEpetraModel_->setOneTimeDirichletBeta(1.0);
207 setOneTimeDirichletBeta(1.0,*this->getUnderlyingModel());
211 MEB::OutArgs<Scalar> outArgs = me->createOutArgs();
212 outArgs.set_W_op(mass_);
215 me->evalModel(inArgs_new,outArgs);
221 invMassMatrix_ = Thyra::inverse<Scalar>(*me->get_W_factory(),mass_);
225 Teuchos::RCP<Thyra::VectorBase<Scalar> > ones = Thyra::createMember(*mass_->domain());
226 Thyra::assign(ones.ptr(),1.0);
228 RCP<Thyra::VectorBase<Scalar> > invLumpMass = Thyra::createMember(*mass_->range());
229 Thyra::apply(*mass_,Thyra::NOTRANS,*ones,invLumpMass.ptr());
230 Thyra::reciprocal(*invLumpMass,invLumpMass.ptr());
232 invMassMatrix_ = Thyra::diagonal(invLumpMass);
236template<
typename Scalar>
240 typedef Thyra::ModelEvaluatorBase MEB;
242 MEB::InArgsSetup<Scalar> inArgs(this->getUnderlyingModel()->createInArgs());
243 inArgs.setModelEvalDescription(this->description());
244 inArgs.setSupports(MEB::IN_ARG_alpha,
true);
245 inArgs.setSupports(MEB::IN_ARG_beta,
true);
246 inArgs.setSupports(MEB::IN_ARG_x_dot,
true);
247 prototypeInArgs_ = inArgs;
249 MEB::OutArgsSetup<Scalar> outArgs(this->getUnderlyingModel()->createOutArgs());
250 outArgs.setModelEvalDescription(this->description());
251 outArgs.setSupports(MEB::OUT_ARG_W,
false);
252 outArgs.setSupports(MEB::OUT_ARG_W_op,
false);
253 prototypeOutArgs_ = outArgs;
256template<
typename Scalar>
261 using Teuchos::ptrFromRef;
262 using Teuchos::ptr_dynamic_cast;
265 Ptr<const panzer::ModelEvaluator<Scalar> > panzerModel = ptr_dynamic_cast<const panzer::ModelEvaluator<Scalar> >(ptrFromRef(me));
266 if(panzerModel!=Teuchos::null) {
267 panzerModel->setOneTimeDirichletBeta(beta);
271#ifdef PANZER_HAVE_EPETRA_STACK
272 Ptr<const Thyra::EpetraModelEvaluator> epModel = ptr_dynamic_cast<const Thyra::EpetraModelEvaluator>(ptrFromRef(me));
273 if(epModel!=Teuchos::null) {
274 Ptr<const panzer::ModelEvaluator_Epetra> panzerEpetraModel
275 = ptr_dynamic_cast<const panzer::ModelEvaluator_Epetra>(epModel->getEpetraModel().ptr());
277 if(panzerEpetraModel!=Teuchos::null) {
278 panzerEpetraModel->setOneTimeDirichletBeta(beta);
288 Ptr<const Thyra::ModelEvaluatorDelegatorBase<Scalar> > delegator
289 = ptr_dynamic_cast<const Thyra::ModelEvaluatorDelegatorBase<Scalar> >(ptrFromRef(me));
290 if(delegator!=Teuchos::null) {
291 setOneTimeDirichletBeta(beta,*delegator->getUnderlyingModel());
295 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::logic_error,
296 "panzer::ExplicitModelEvaluator::setOneTimeDirichletBeta can't find a panzer::ME or a panzer::EpetraME. "
297 "The deepest model is also not a delegator. Thus the recursion failed and an exception was generated.");
Teuchos::RCP< const panzer::ModelEvaluator< Scalar > > panzerModel_
Access to the panzer model evaluator pointer (thyra version)
void evalModelImpl(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
Evaluation of model, applies the mass matrix to the RHS.
void setOneTimeDirichletBeta(double beta, const Thyra::ModelEvaluator< Scalar > &me) const
Thyra::ModelEvaluatorBase::InArgs< Scalar > getNominalValues() const
Build the nominal values, modifies the underlying models in args slightly.
void buildInverseMassMatrix(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs) const
Thyra::ModelEvaluatorBase::OutArgs< Scalar > createOutArgs() const
Build the out args, modifies the underlying models in args slightly.
Teuchos::RCP< panzer::ModelEvaluator< Scalar > > getPanzerUnderlyingModel()
Get the underlying panzer::ModelEvaluator.
void buildArgsPrototypes()
Build prototype in/out args.
Thyra::ModelEvaluatorBase::InArgs< Scalar > createInArgs() const
Build the in args, modifies the underlying models in args slightly.
bool applyMassInverse_
Apply mass matrix inverse within the evaluator.