9 #ifndef Tempus_CombinedForwardSensitivityModelEvaluator_impl_hpp
10 #define Tempus_CombinedForwardSensitivityModelEvaluator_impl_hpp
12 #include "Thyra_DefaultMultiVectorProductVector.hpp"
16 #include "Thyra_VectorStdOps.hpp"
17 #include "Thyra_MultiVectorStdOps.hpp"
21 template <
typename Scalar>
34 sens_residual_model_(sens_residual_model),
35 sens_solve_model_(sens_solve_model),
36 dxdp_init_(dxdp_init),
37 dx_dotdp_init_(dx_dotdp_init),
38 dx_dotdotdp_init_(dx_dotdotdp_init),
42 xdot_tangent_index_(2),
43 xdotdot_tangent_index_(3),
44 use_dfdp_as_tangent_(false),
45 use_dgdp_as_tangent_(false),
55 if (pList != Teuchos::null) *pl = *pList;
59 p_index_ = pl->
get<
int>(
"Sensitivity Parameter Index");
83 MEB::InArgs<Scalar> me_inArgs =
model_->createInArgs();
87 MEB::InArgsSetup<Scalar> inArgs;
88 inArgs.setModelEvalDescription(this->
description());
89 inArgs.setSupports(MEB::IN_ARG_x);
90 if (me_inArgs.supports(MEB::IN_ARG_x_dot))
91 inArgs.setSupports(MEB::IN_ARG_x_dot);
92 if (me_inArgs.supports(MEB::IN_ARG_t)) inArgs.setSupports(MEB::IN_ARG_t);
93 if (me_inArgs.supports(MEB::IN_ARG_alpha))
94 inArgs.setSupports(MEB::IN_ARG_alpha);
95 if (me_inArgs.supports(MEB::IN_ARG_beta))
96 inArgs.setSupports(MEB::IN_ARG_beta);
97 if (me_inArgs.supports(MEB::IN_ARG_W_x_dot_dot_coeff))
98 inArgs.setSupports(MEB::IN_ARG_W_x_dot_dot_coeff);
101 inArgs.set_Np(me_inArgs.Np());
104 MEB::OutArgs<Scalar> me_outArgs =
model_->createOutArgs();
107 MEB::OutArgsSetup<Scalar> outArgs;
108 outArgs.setModelEvalDescription(this->
description());
109 outArgs.set_Np_Ng(me_outArgs.Np(), me_outArgs.Ng() +
g_offset_);
110 outArgs.setSupports(MEB::OUT_ARG_f);
111 if (sens_mer_outArgs.supports(MEB::OUT_ARG_W_op) ||
112 sens_mes_outArgs.supports(MEB::OUT_ARG_W_op))
113 outArgs.setSupports(MEB::OUT_ARG_W_op);
117 for (
int j = 0; j < me_outArgs.Ng(); ++j) {
118 outArgs.setSupports(MEB::OUT_ARG_DgDx_dot, j +
g_offset_,
119 me_outArgs.supports(MEB::OUT_ARG_DgDx_dot, j));
120 outArgs.setSupports(MEB::OUT_ARG_DgDx, j +
g_offset_,
121 me_outArgs.supports(MEB::OUT_ARG_DgDx, j));
122 for (
int l = 0; l < me_outArgs.Np(); ++l) {
123 outArgs.setSupports(MEB::OUT_ARG_DgDp, j +
g_offset_, l,
124 me_outArgs.supports(MEB::OUT_ARG_DgDp, j, l));
134 if (me_inArgs.supports(MEB::IN_ARG_x_dot)) {
139 .supports(MEB::DERIV_MV_JACOBIAN_FORM));
142 template <
typename Scalar>
146 return model_->get_p_space(p);
149 template <
typename Scalar>
153 return model_->get_p_names(p);
156 template <
typename Scalar>
160 return x_dxdp_space_;
163 template <
typename Scalar>
167 return f_dfdp_space_;
170 template <
typename Scalar>
178 return model_->get_g_space(g_index_);
180 return model_->get_g_space(j - g_offset_);
183 template <
typename Scalar>
190 for (
int i = 0; i < names.
size(); ++i)
191 names[i] = names[i] +
"_reduced sensitivity";
195 return model_->get_g_names(g_index_);
197 return model_->get_g_names(j - g_offset_);
200 template <
typename Scalar>
205 sens_solve_model_->create_W_op();
206 return Thyra::nonconstMultiVectorLinearOp(op, num_param_ + 1);
209 template <
typename Scalar>
214 return model_->create_DgDx_dot_op(j - g_offset_);
217 template <
typename Scalar>
221 return model_->create_DgDx_op(j - g_offset_);
224 template <
typename Scalar>
229 return model_->create_DgDp_op(j - g_offset_, l);
232 template <
typename Scalar>
237 sens_solve_model_->get_W_factory();
238 if (factory == Teuchos::null)
239 return Teuchos::null;
241 return Thyra::multiVectorLinearOpWithSolveFactory(factory, f_dfdp_space_,
245 template <
typename Scalar>
249 return prototypeInArgs_;
252 template <
typename Scalar>
260 using Teuchos::rcp_dynamic_cast;
262 MEB::InArgs<Scalar> me_nominal = model_->getNominalValues();
263 MEB::InArgs<Scalar> nominal = this->createInArgs();
269 RCP<const Thyra::VectorBase<Scalar> > me_x = me_nominal.get_x();
270 if (me_x != Teuchos::null) {
271 RCP<Thyra::VectorBase<Scalar> > x = Thyra::createMember(*x_dxdp_space_);
272 RCP<DMVPV> x_dxdp = rcp_dynamic_cast<DMVPV>(x,
true);
273 Thyra::assign(x_dxdp->getNonconstMultiVector()->col(0).ptr(), *me_x);
274 if (dxdp_init_ == Teuchos::null)
275 Thyra::assign(x_dxdp->getNonconstMultiVector()->subView(rng).ptr(), zero);
277 Thyra::assign(x_dxdp->getNonconstMultiVector()->subView(rng).ptr(),
283 RCP<const Thyra::VectorBase<Scalar> > me_xdot;
284 if (me_nominal.supports(MEB::IN_ARG_x_dot)) me_xdot = me_nominal.get_x_dot();
285 if (me_xdot != Teuchos::null) {
286 RCP<Thyra::VectorBase<Scalar> > xdot = Thyra::createMember(*x_dxdp_space_);
287 RCP<DMVPV> xdot_dxdp = rcp_dynamic_cast<DMVPV>(xdot,
true);
288 Thyra::assign(xdot_dxdp->getNonconstMultiVector()->col(0).ptr(), *me_xdot);
289 if (dx_dotdp_init_ == Teuchos::null)
290 Thyra::assign(xdot_dxdp->getNonconstMultiVector()->subView(rng).ptr(),
293 Thyra::assign(xdot_dxdp->getNonconstMultiVector()->subView(rng).ptr(),
295 nominal.set_x_dot(xdot);
300 RCP<const Thyra::VectorBase<Scalar> > me_xdotdot;
301 if (me_nominal.supports(MEB::IN_ARG_x_dot_dot))
302 me_xdotdot = me_nominal.get_x_dot_dot();
303 if (me_xdotdot != Teuchos::null) {
304 RCP<Thyra::VectorBase<Scalar> > xdotdot =
305 Thyra::createMember(*x_dxdp_space_);
306 RCP<DMVPV> xdotdot_dxdp = rcp_dynamic_cast<DMVPV>(xdotdot,
true);
307 Thyra::assign(xdotdot_dxdp->getNonconstMultiVector()->col(0).ptr(),
309 if (dx_dotdotdp_init_ == Teuchos::null)
310 Thyra::assign(xdotdot_dxdp->getNonconstMultiVector()->subView(rng).ptr(),
313 Thyra::assign(xdotdot_dxdp->getNonconstMultiVector()->subView(rng).ptr(),
315 nominal.set_x_dot_dot(xdotdot);
318 const int np = model_->Np();
319 for (
int i = 0; i < np; ++i) nominal.set_p(i, me_nominal.get_p(i));
324 template <
typename Scalar>
328 return prototypeOutArgs_;
331 template <
typename Scalar>
340 using Teuchos::rcp_dynamic_cast;
342 const bool use_tangent = use_dfdp_as_tangent_ || use_dgdp_as_tangent_;
345 RCP<const Thyra::VectorBase<Scalar> > x, xdot, xdotdot;
346 RCP<const Thyra::MultiVectorBase<Scalar> > dxdp, dxdotdp, dxdotdotdp;
347 RCP<const Thyra::VectorBase<Scalar> > dxdp_vec, dxdotdp_vec, dxdotdotdp_vec;
348 MEB::InArgs<Scalar> me_inArgs = model_->getNominalValues();
349 RCP<const DMVPV> x_dxdp = rcp_dynamic_cast<
const DMVPV>(inArgs.
get_x(),
true);
350 x = x_dxdp->getMultiVector()->col(0);
351 dxdp = x_dxdp->getMultiVector()->subView(
Range1D(1, num_param_));
354 dxdp_vec = Thyra::multiVectorProductVector(dxdp_space_, dxdp);
355 if (use_dfdp_as_tangent_) me_inArgs.set_p(x_tangent_index_, dxdp_vec);
357 if (me_inArgs.supports(MEB::IN_ARG_x_dot)) {
358 if (inArgs.
get_x_dot() != Teuchos::null) {
359 RCP<const DMVPV> xdot_dxdotdp =
360 rcp_dynamic_cast<
const DMVPV>(inArgs.
get_x_dot(),
true);
361 xdot = xdot_dxdotdp->getMultiVector()->col(0);
362 dxdotdp = xdot_dxdotdp->getMultiVector()->subView(
Range1D(1, num_param_));
363 me_inArgs.set_x_dot(xdot);
365 dxdotdp_vec = Thyra::multiVectorProductVector(dxdp_space_, dxdotdp);
366 if (use_dfdp_as_tangent_)
367 me_inArgs.set_p(xdot_tangent_index_, dxdotdp_vec);
371 me_inArgs.set_x_dot(Teuchos::null);
373 if (me_inArgs.supports(MEB::IN_ARG_x_dot_dot)) {
375 RCP<const DMVPV> xdotdot_dxdotdotdp =
377 xdotdot = xdotdot_dxdotdotdp->getMultiVector()->col(0);
379 xdotdot_dxdotdotdp->getMultiVector()->subView(
Range1D(1, num_param_));
380 me_inArgs.set_x_dot_dot(xdotdot);
383 Thyra::multiVectorProductVector(dxdp_space_, dxdotdotdp);
384 if (use_dfdp_as_tangent_)
385 me_inArgs.set_p(xdotdot_tangent_index_, dxdotdotdp_vec);
389 me_inArgs.set_x_dot_dot(Teuchos::null);
391 if (me_inArgs.supports(MEB::IN_ARG_t)) me_inArgs.set_t(inArgs.
get_t());
392 if (me_inArgs.supports(MEB::IN_ARG_alpha))
394 if (me_inArgs.supports(MEB::IN_ARG_beta))
395 me_inArgs.set_beta(inArgs.
get_beta());
396 if (me_inArgs.supports(MEB::IN_ARG_W_x_dot_dot_coeff))
402 const int np = me_inArgs.Np();
403 for (
int i = 0; i < np; ++i) {
404 if (inArgs.
get_p(i) != Teuchos::null)
406 (use_tangent && i != x_tangent_index_ && i != xdot_tangent_index_ &&
407 i != xdotdot_tangent_index_))
408 me_inArgs.set_p(i, inArgs.
get_p(i));
412 RCP<Thyra::VectorBase<Scalar> > f;
413 RCP<Thyra::MultiVectorBase<Scalar> > dfdp;
414 MEB::OutArgs<Scalar> me_outArgs = model_->createOutArgs();
415 if (outArgs.
get_f() != Teuchos::null) {
416 RCP<DMVPV> f_dfdp = rcp_dynamic_cast<DMVPV>(outArgs.
get_f(),
true);
417 f = f_dfdp->getNonconstMultiVector()->col(0);
418 dfdp = f_dfdp->getNonconstMultiVector()->subView(
Range1D(1, num_param_));
420 me_outArgs.set_DfDp(p_index_, dfdp);
422 if (outArgs.
supports(MEB::OUT_ARG_W_op) &&
423 outArgs.
get_W_op() != Teuchos::null &&
424 model_.ptr() == sens_solve_model_.ptr()) {
425 RCP<Thyra::LinearOpBase<Scalar> > op = outArgs.
get_W_op();
426 RCP<Thyra::MultiVectorLinearOp<Scalar> > mv_op =
428 me_outArgs.set_W_op(mv_op->getNonconstLinearOp());
430 for (
int j = g_offset_; j < outArgs.
Ng(); ++j) {
431 if (!me_outArgs.supports(MEB::OUT_ARG_DgDx_dot, j - g_offset_).none())
432 me_outArgs.set_DgDx_dot(j - g_offset_, outArgs.
get_DgDx_dot(j));
433 if (!me_outArgs.supports(MEB::OUT_ARG_DgDx, j - g_offset_).none())
434 me_outArgs.set_DgDx(j - g_offset_, outArgs.
get_DgDx(j));
435 for (
int l = 0; l < outArgs.
Np(); ++l)
436 if (!me_outArgs.supports(MEB::OUT_ARG_DgDp, j - g_offset_, l).none())
437 me_outArgs.set_DgDp(j - g_offset_, l, outArgs.
get_DgDp(j, l));
439 if (g_index_ >= 0 && outArgs.
get_g(1) != Teuchos::null)
440 me_outArgs.set_g(g_index_, outArgs.
get_g(1));
443 model_->evalModel(me_inArgs, me_outArgs);
446 if (outArgs.
supports(MEB::OUT_ARG_W_op) &&
447 outArgs.
get_W_op() != Teuchos::null &&
448 model_.ptr() != sens_solve_model_.ptr()) {
449 MEB::OutArgs<Scalar> sens_me_outArgs = sens_solve_model_->createOutArgs();
450 RCP<Thyra::LinearOpBase<Scalar> > op = outArgs.
get_W_op();
451 RCP<Thyra::MultiVectorLinearOp<Scalar> > mv_op =
453 sens_me_outArgs.set_W_op(mv_op->getNonconstLinearOp());
454 sens_solve_model_->evalModel(me_inArgs, sens_me_outArgs);
460 if (!use_dfdp_as_tangent_) {
461 if (dxdp != Teuchos::null && dfdp != Teuchos::null) {
462 if (my_dfdx_ == Teuchos::null)
463 my_dfdx_ = sens_residual_model_->create_W_op();
464 MEB::OutArgs<Scalar> meo = sens_residual_model_->createOutArgs();
465 meo.set_W_op(my_dfdx_);
466 if (me_inArgs.supports(MEB::IN_ARG_alpha)) me_inArgs.set_alpha(0.0);
467 if (me_inArgs.supports(MEB::IN_ARG_beta)) me_inArgs.set_beta(1.0);
468 if (me_inArgs.supports(MEB::IN_ARG_W_x_dot_dot_coeff))
469 me_inArgs.set_W_x_dot_dot_coeff(0.0);
470 sens_residual_model_->evalModel(me_inArgs, meo);
474 if (dxdotdp != Teuchos::null && dfdp != Teuchos::null) {
475 if (my_dfdxdot_ == Teuchos::null)
476 my_dfdxdot_ = sens_residual_model_->create_W_op();
477 MEB::OutArgs<Scalar> meo = sens_residual_model_->createOutArgs();
478 meo.set_W_op(my_dfdxdot_);
479 if (me_inArgs.supports(MEB::IN_ARG_alpha)) me_inArgs.set_alpha(1.0);
480 if (me_inArgs.supports(MEB::IN_ARG_beta)) me_inArgs.set_beta(0.0);
481 if (me_inArgs.supports(MEB::IN_ARG_W_x_dot_dot_coeff))
482 me_inArgs.set_W_x_dot_dot_coeff(0.0);
483 sens_residual_model_->evalModel(me_inArgs, meo);
484 my_dfdxdot_->apply(
Thyra::NOTRANS, *dxdotdp, dfdp.ptr(), Scalar(1.0),
487 if (dxdotdotdp != Teuchos::null && dfdp != Teuchos::null) {
488 if (my_dfdxdotdot_ == Teuchos::null)
489 my_dfdxdotdot_ = sens_residual_model_->create_W_op();
490 MEB::OutArgs<Scalar> meo = sens_residual_model_->createOutArgs();
491 meo.set_W_op(my_dfdxdotdot_);
492 if (me_inArgs.supports(MEB::IN_ARG_alpha)) me_inArgs.set_alpha(0.0);
493 if (me_inArgs.supports(MEB::IN_ARG_beta)) me_inArgs.set_beta(0.0);
494 if (me_inArgs.supports(MEB::IN_ARG_W_x_dot_dot_coeff))
495 me_inArgs.set_W_x_dot_dot_coeff(1.0);
496 sens_residual_model_->evalModel(me_inArgs, meo);
498 Scalar(1.0), Scalar(1.0));
503 if (g_index_ >= 0 && outArgs.
get_g(0) != Teuchos::null) {
504 MEB::OutArgs<Scalar> meo = model_->createOutArgs();
505 RCP<Thyra::MultiVectorBase<Scalar> > dgdp =
506 rcp_dynamic_cast<DMVPV>(outArgs.
get_g(0),
true)
507 ->getNonconstMultiVector();
508 RCP<Thyra::MultiVectorBase<Scalar> > dgdp_trans;
509 MEB::DerivativeSupport dgdp_support =
510 meo.supports(MEB::OUT_ARG_DgDp, g_index_, p_index_);
511 if (dgdp_support.supports(MEB::DERIV_MV_JACOBIAN_FORM))
512 meo.set_DgDp(g_index_, p_index_,
513 MEB::Derivative<Scalar>(dgdp, MEB::DERIV_MV_JACOBIAN_FORM));
514 else if (dgdp_support.supports(MEB::DERIV_MV_GRADIENT_FORM)) {
515 dgdp_trans = createMembers(model_->get_p_space(p_index_), num_response_);
518 MEB::Derivative<Scalar>(dgdp_trans, MEB::DERIV_MV_GRADIENT_FORM));
522 true, std::logic_error,
523 "Operator form of dg/dp not supported for reduced sensitivity");
525 if (!use_dgdp_as_tangent_ || dxdp == Teuchos::null) {
526 MEB::DerivativeSupport dgdx_support =
527 meo.supports(MEB::OUT_ARG_DgDx, g_index_);
528 if (dgdx_support.supports(MEB::DERIV_LINEAR_OP)) {
529 if (my_dgdx_ == Teuchos::null)
530 my_dgdx_ = model_->create_DgDx_op(g_index_);
531 meo.set_DgDx(g_index_, my_dgdx_);
533 else if (dgdx_support.supports(MEB::DERIV_MV_GRADIENT_FORM)) {
534 if (my_dgdx_mv_ == Teuchos::null)
536 createMembers(model_->get_g_space(g_index_), num_param_);
537 meo.set_DgDx(g_index_, MEB::Derivative<Scalar>(
538 my_dgdx_mv_, MEB::DERIV_MV_GRADIENT_FORM));
542 true, std::logic_error,
543 "Jacobian form of dg/dx not supported for reduced sensitivity");
548 if (!use_dgdp_as_tangent_ && use_dfdp_as_tangent_ &&
549 dxdp_vec != Teuchos::null)
550 me_inArgs.set_p(x_tangent_index_, Teuchos::null);
551 if (!use_dgdp_as_tangent_ && use_dfdp_as_tangent_ &&
552 dxdotdp_vec != Teuchos::null)
553 me_inArgs.set_p(xdot_tangent_index_, Teuchos::null);
554 if (!use_dgdp_as_tangent_ && use_dfdp_as_tangent_ &&
555 dxdotdotdp_vec != Teuchos::null)
556 me_inArgs.set_p(xdotdot_tangent_index_, Teuchos::null);
559 if (use_dgdp_as_tangent_ && !use_dfdp_as_tangent_ &&
560 dxdp_vec != Teuchos::null)
561 me_inArgs.set_p(x_tangent_index_, dxdp_vec);
562 if (use_dgdp_as_tangent_ && !use_dfdp_as_tangent_ &&
563 dxdotdp_vec != Teuchos::null)
564 me_inArgs.set_p(xdot_tangent_index_, dxdotdp_vec);
565 if (use_dgdp_as_tangent_ && !use_dfdp_as_tangent_ &&
566 dxdotdotdp_vec != Teuchos::null)
567 me_inArgs.set_p(xdotdot_tangent_index_, dxdotdotdp_vec);
569 model_->evalModel(me_inArgs, meo);
572 if (dgdp_trans != Teuchos::null) {
575 for (
int i = 0; i < num_param_; ++i)
576 for (
int j = 0; j < num_response_; ++j)
577 dgdp_view(j, i) = dgdp_trans_view(i, j);
582 if (!use_dgdp_as_tangent_ && dxdp != Teuchos::null) {
583 MEB::DerivativeSupport dgdx_support =
584 meo.supports(MEB::OUT_ARG_DgDx, g_index_);
585 if (dgdx_support.supports(MEB::DERIV_LINEAR_OP))
588 else if (dgdx_support.supports(MEB::DERIV_MV_GRADIENT_FORM))
589 my_dgdx_mv_->apply(
Thyra::TRANS, *dxdp, dgdp.ptr(), Scalar(1.0),
593 true, std::logic_error,
594 "Jacobian form of dg/dx not supported for reduced sensitivity");
599 template <
class Scalar>
604 pl->
set<
bool>(
"Use DfDp as Tangent",
false);
605 pl->
set<
bool>(
"Use DgDp as Tangent",
false);
606 pl->
set<
int>(
"Sensitivity Parameter Index", 0);
607 pl->
set<
int>(
"Response Function Index", -1);
608 pl->
set<
int>(
"Sensitivity X Tangent Index", 1);
609 pl->
set<
int>(
"Sensitivity X-Dot Tangent Index", 2);
610 pl->
set<
int>(
"Sensitivity X-Dot-Dot Tangent Index", 3);
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_p_space(int p) const
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_g_space(int j) const
Teuchos::RCP< Thyra::LinearOpBase< Scalar > > create_DgDx_dot_op(int j) const
RCP< const VectorBase< Scalar > > get_x_dot() const
T & get(const std::string &name, T def_value)
Thyra::ModelEvaluatorBase::InArgs< Scalar > createInArgs() const
Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > sens_residual_model_
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
Teuchos::RCP< const DMVPVS > dxdp_space_
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Teuchos::RCP< const Teuchos::Array< std::string > > get_p_names(int p) const
Teuchos::RCP< Thyra::LinearOpBase< Scalar > > create_DgDp_op(int j, int l) const
Teuchos::RCP< const DMVPVS > f_dfdp_space_
Evaluation< VectorBase< Scalar > > get_g(int j) const
Thyra::ModelEvaluatorBase::InArgs< Scalar > prototypeInArgs_
static Teuchos::RCP< const Teuchos::ParameterList > getValidParameters()
Evaluation< VectorBase< Scalar > > get_f() const
Teuchos::RCP< const DMVPVS > dgdp_space_
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_x_space() const
Teuchos::RCP< Thyra::LinearOpBase< Scalar > > create_DgDx_op(int j) const
bool supports(EOutArgsMembers arg) const
Derivative< Scalar > get_DgDp(int j, int l) const
Thyra::ModelEvaluatorBase::OutArgs< Scalar > createOutArgsImpl() const
Thyra::ModelEvaluatorBase::InArgs< Scalar > getNominalValues() const
bool use_dgdp_as_tangent_
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
CombinedForwardSensitivityModelEvaluator(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &sens_residual_model, const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &sens_solve_model, const Teuchos::RCP< const Teuchos::ParameterList > &pList=Teuchos::null, const Teuchos::RCP< MultiVector > &dxdp_init=Teuchos::null, const Teuchos::RCP< MultiVector > &dx_dotdp_init=Teuchos::null, const Teuchos::RCP< MultiVector > &dx_dotdot_dp_init=Teuchos::null)
Constructor.
Teuchos::ArrayView< const std::string > get_g_names(int j) const
void evalModelImpl(const Thyra::ModelEvaluatorBase::InArgs< Scalar > &inArgs, const Thyra::ModelEvaluatorBase::OutArgs< Scalar > &outArgs) const
Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > sens_solve_model_
virtual std::string description() const
void validateParametersAndSetDefaults(ParameterList const &validParamList, int const depth=1000)
Teuchos::RCP< const DMVPVS > dfdp_space_
Scalar get_W_x_dot_dot_coeff() const
Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > model_
int xdotdot_tangent_index_
Teuchos::RCP< const Thyra::LinearOpWithSolveFactoryBase< Scalar > > get_W_factory() const
Implicit concrete LinearOpBase subclass that takes a flattended out multi-vector and performs a multi...
bool use_dfdp_as_tangent_
Teuchos::RCP< const DMVPVS > x_dxdp_space_
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_f_space() const
Teuchos::RCP< Thyra::LinearOpBase< Scalar > > create_W_op() const
Thyra::ModelEvaluatorBase::OutArgs< Scalar > prototypeOutArgs_
Derivative< Scalar > get_DgDx_dot(int j) const
RCP< const VectorBase< Scalar > > get_x_dot_dot() const
#define TEUCHOS_ASSERT(assertion_test)
RCP< LinearOpBase< Scalar > > get_W_op() const
RCP< const VectorBase< Scalar > > get_x() const
RCP< const VectorBase< Scalar > > get_p(int l) const
Derivative< Scalar > get_DgDx(int j) const