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>
32 sens_residual_model_(sens_residual_model),
33 sens_solve_model_(sens_solve_model),
34 dxdp_init_(dxdp_init),
35 dx_dotdp_init_(dx_dotdp_init),
36 dx_dotdotdp_init_(dx_dotdotdp_init),
40 xdot_tangent_index_(2),
41 xdotdot_tangent_index_(3),
42 use_dfdp_as_tangent_(false),
43 use_dgdp_as_tangent_(false),
53 if (pList != Teuchos::null)
58 p_index_ = pl->
get<
int>(
"Sensitivity Parameter Index");
72 Thyra::multiVectorProductVectorSpace(
model_->get_x_space(), num_param_+1);
76 Thyra::multiVectorProductVectorSpace(
model_->get_f_space(), num_param_+1);
79 Thyra::multiVectorProductVectorSpace(
model_->get_g_space(g_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))
93 inArgs.setSupports(MEB::IN_ARG_t);
94 if (me_inArgs.supports(MEB::IN_ARG_alpha))
95 inArgs.setSupports(MEB::IN_ARG_alpha);
96 if (me_inArgs.supports(MEB::IN_ARG_beta))
97 inArgs.setSupports(MEB::IN_ARG_beta);
98 if (me_inArgs.supports(MEB::IN_ARG_W_x_dot_dot_coeff))
99 inArgs.setSupports(MEB::IN_ARG_W_x_dot_dot_coeff);
102 inArgs.set_Np(me_inArgs.Np());
105 MEB::OutArgs<Scalar> me_outArgs =
model_->createOutArgs();
108 MEB::OutArgsSetup<Scalar> outArgs;
109 outArgs.setModelEvalDescription(this->
description());
110 outArgs.set_Np_Ng(me_outArgs.Np(),me_outArgs.Ng()+
g_offset_);
111 outArgs.setSupports(MEB::OUT_ARG_f);
112 if (sens_mer_outArgs.supports(MEB::OUT_ARG_W_op) ||
113 sens_mes_outArgs.supports(MEB::OUT_ARG_W_op))
114 outArgs.setSupports(MEB::OUT_ARG_W_op);
118 for (
int j=0; j<me_outArgs.Ng(); ++j) {
119 outArgs.setSupports(MEB::OUT_ARG_DgDx_dot, j+
g_offset_,
120 me_outArgs.supports(MEB::OUT_ARG_DgDx_dot, j));
121 outArgs.setSupports(MEB::OUT_ARG_DgDx, j+
g_offset_,
122 me_outArgs.supports(MEB::OUT_ARG_DgDx, j));
123 for (
int l=0; l<me_outArgs.Np(); ++l) {
124 outArgs.setSupports(MEB::OUT_ARG_DgDp, j+
g_offset_, l,
125 me_outArgs.supports(MEB::OUT_ARG_DgDp, j, l));
135 if (me_inArgs.supports(MEB::IN_ARG_x_dot)) {
139 TEUCHOS_ASSERT(me_outArgs.supports(MEB::OUT_ARG_DfDp, p_index_).supports(MEB::DERIV_MV_JACOBIAN_FORM));
142 template <
typename Scalar>
147 return model_->get_p_space(p);
150 template <
typename Scalar>
155 return model_->get_p_names(p);
158 template <
typename Scalar>
163 return x_dxdp_space_;
166 template <
typename Scalar>
171 return f_dfdp_space_;
174 template <
typename Scalar>
183 return model_->get_g_space(g_index_);
185 return model_->get_g_space(j-g_offset_);
188 template <
typename Scalar>
196 for (
int i=0; i<names.
size(); ++i)
197 names[i] = names[i] +
"_reduced sensitivity";
201 return model_->get_g_names(g_index_);
203 return model_->get_g_names(j-g_offset_);
206 template <
typename Scalar>
212 sens_solve_model_->create_W_op();
213 return Thyra::nonconstMultiVectorLinearOp(op, num_param_+1);
216 template <
typename Scalar>
221 return model_->create_DgDx_dot_op(j-g_offset_);
224 template <
typename Scalar>
229 return model_->create_DgDx_op(j-g_offset_);
232 template <
typename Scalar>
237 return model_->create_DgDp_op(j-g_offset_,l);
240 template <
typename Scalar>
246 sens_solve_model_->get_W_factory();
247 if (factory == Teuchos::null)
248 return Teuchos::null;
250 return Thyra::multiVectorLinearOpWithSolveFactory(factory,
255 template <
typename Scalar>
260 return prototypeInArgs_;
263 template <
typename Scalar>
271 using Teuchos::rcp_dynamic_cast;
274 MEB::InArgs<Scalar> me_nominal = model_->getNominalValues();
275 MEB::InArgs<Scalar> nominal = this->createInArgs();
281 RCP< const Thyra::VectorBase<Scalar> > me_x = me_nominal.get_x();
282 if (me_x != Teuchos::null) {
283 RCP< Thyra::VectorBase<Scalar> > x = Thyra::createMember(*x_dxdp_space_);
284 RCP<DMVPV> x_dxdp = rcp_dynamic_cast<DMVPV>(x,
true);
285 Thyra::assign(x_dxdp->getNonconstMultiVector()->col(0).ptr(), *me_x);
286 if (dxdp_init_ == Teuchos::null)
287 Thyra::assign(x_dxdp->getNonconstMultiVector()->subView(rng).ptr(),
290 Thyra::assign(x_dxdp->getNonconstMultiVector()->subView(rng).ptr(),
296 RCP< const Thyra::VectorBase<Scalar> > me_xdot;
297 if (me_nominal.supports(MEB::IN_ARG_x_dot))
298 me_xdot = me_nominal.get_x_dot();
299 if (me_xdot != Teuchos::null) {
300 RCP< Thyra::VectorBase<Scalar> > xdot = Thyra::createMember(*x_dxdp_space_);
301 RCP<DMVPV> xdot_dxdp = rcp_dynamic_cast<DMVPV>(xdot,
true);
302 Thyra::assign(xdot_dxdp->getNonconstMultiVector()->col(0).ptr(), *me_xdot);
303 if (dx_dotdp_init_ == Teuchos::null)
304 Thyra::assign(xdot_dxdp->getNonconstMultiVector()->subView(rng).ptr(),
307 Thyra::assign(xdot_dxdp->getNonconstMultiVector()->subView(rng).ptr(),
309 nominal.set_x_dot(xdot);
313 RCP< const Thyra::VectorBase<Scalar> > me_xdotdot;
314 if (me_nominal.supports(MEB::IN_ARG_x_dot_dot))
315 me_xdotdot = me_nominal.get_x_dot_dot();
316 if (me_xdotdot != Teuchos::null) {
317 RCP< Thyra::VectorBase<Scalar> > xdotdot =
318 Thyra::createMember(*x_dxdp_space_);
319 RCP<DMVPV> xdotdot_dxdp = rcp_dynamic_cast<DMVPV>(xdotdot,
true);
320 Thyra::assign(xdotdot_dxdp->getNonconstMultiVector()->col(0).ptr(), *me_xdotdot);
321 if (dx_dotdotdp_init_ == Teuchos::null)
322 Thyra::assign(xdotdot_dxdp->getNonconstMultiVector()->subView(rng).ptr(),
325 Thyra::assign(xdotdot_dxdp->getNonconstMultiVector()->subView(rng).ptr(),
327 nominal.set_x_dot_dot(xdotdot);
330 const int np = model_->Np();
331 for (
int i=0; i<np; ++i)
332 nominal.set_p(i, me_nominal.get_p(i));
337 template <
typename Scalar>
342 return prototypeOutArgs_;
345 template <
typename Scalar>
354 using Teuchos::rcp_dynamic_cast;
357 const bool use_tangent = use_dfdp_as_tangent_ || use_dgdp_as_tangent_;
360 RCP< const Thyra::VectorBase<Scalar> > x, xdot, xdotdot;
361 RCP< const Thyra::MultiVectorBase<Scalar> > dxdp, dxdotdp, dxdotdotdp;
362 RCP< const Thyra::VectorBase<Scalar> > dxdp_vec, dxdotdp_vec, dxdotdotdp_vec;
363 MEB::InArgs<Scalar> me_inArgs = model_->getNominalValues();
364 RCP<const DMVPV> x_dxdp = rcp_dynamic_cast<
const DMVPV>(inArgs.
get_x(),
true);
365 x = x_dxdp->getMultiVector()->col(0);
366 dxdp = x_dxdp->getMultiVector()->subView(
Range1D(1,num_param_));
369 dxdp_vec = Thyra::multiVectorProductVector(dxdp_space_, dxdp);
370 if (use_dfdp_as_tangent_)
371 me_inArgs.set_p(x_tangent_index_, dxdp_vec);
373 if (me_inArgs.supports(MEB::IN_ARG_x_dot)) {
374 if (inArgs.
get_x_dot() != Teuchos::null) {
375 RCP<const DMVPV> xdot_dxdotdp =
376 rcp_dynamic_cast<
const DMVPV>(inArgs.
get_x_dot(),
true);
377 xdot = xdot_dxdotdp->getMultiVector()->col(0);
378 dxdotdp = xdot_dxdotdp->getMultiVector()->subView(
Range1D(1,num_param_));
379 me_inArgs.set_x_dot(xdot);
381 dxdotdp_vec = Thyra::multiVectorProductVector(dxdp_space_, dxdotdp);
382 if (use_dfdp_as_tangent_)
383 me_inArgs.set_p(xdot_tangent_index_, dxdotdp_vec);
387 me_inArgs.set_x_dot(Teuchos::null);
389 if (me_inArgs.supports(MEB::IN_ARG_x_dot_dot)) {
391 RCP<const DMVPV> xdotdot_dxdotdotdp =
393 xdotdot = xdotdot_dxdotdotdp->getMultiVector()->col(0);
394 dxdotdotdp = xdotdot_dxdotdotdp->getMultiVector()->subView(
Range1D(1,num_param_));
395 me_inArgs.set_x_dot_dot(xdotdot);
398 Thyra::multiVectorProductVector(dxdp_space_, dxdotdotdp);
399 if (use_dfdp_as_tangent_)
400 me_inArgs.set_p(xdotdot_tangent_index_, dxdotdotdp_vec);
404 me_inArgs.set_x_dot_dot(Teuchos::null);
406 if (me_inArgs.supports(MEB::IN_ARG_t))
407 me_inArgs.set_t(inArgs.
get_t());
408 if (me_inArgs.supports(MEB::IN_ARG_alpha))
410 if (me_inArgs.supports(MEB::IN_ARG_beta))
411 me_inArgs.set_beta(inArgs.
get_beta());
412 if (me_inArgs.supports(MEB::IN_ARG_W_x_dot_dot_coeff))
418 const int np = me_inArgs.Np();
419 for (
int i=0; i<np; ++i) {
420 if (inArgs.
get_p(i) != Teuchos::null)
422 (use_tangent && i != x_tangent_index_ &&
423 i != xdot_tangent_index_ &&
424 i != xdotdot_tangent_index_ ))
425 me_inArgs.set_p(i, inArgs.
get_p(i));
429 RCP< Thyra::VectorBase<Scalar> > f;
430 RCP< Thyra::MultiVectorBase<Scalar> > dfdp;
431 MEB::OutArgs<Scalar> me_outArgs = model_->createOutArgs();
432 if (outArgs.
get_f() != Teuchos::null) {
433 RCP<DMVPV> f_dfdp = rcp_dynamic_cast<DMVPV>(outArgs.
get_f(),
true);
434 f = f_dfdp->getNonconstMultiVector()->col(0);
435 dfdp = f_dfdp->getNonconstMultiVector()->subView(
Range1D(1,num_param_));
437 me_outArgs.set_DfDp(p_index_, dfdp);
439 if (outArgs.
supports(MEB::OUT_ARG_W_op) &&
440 outArgs.
get_W_op() != Teuchos::null &&
441 model_.ptr() == sens_solve_model_.ptr()) {
442 RCP< Thyra::LinearOpBase<Scalar> > op = outArgs.
get_W_op();
443 RCP< Thyra::MultiVectorLinearOp<Scalar> > mv_op =
445 me_outArgs.set_W_op(mv_op->getNonconstLinearOp());
447 for (
int j=g_offset_; j<outArgs.
Ng(); ++j) {
448 if (!me_outArgs.supports(MEB::OUT_ARG_DgDx_dot,j-g_offset_).none())
449 me_outArgs.set_DgDx_dot(j-g_offset_, outArgs.
get_DgDx_dot(j));
450 if (!me_outArgs.supports(MEB::OUT_ARG_DgDx,j-g_offset_).none())
451 me_outArgs.set_DgDx(j-g_offset_, outArgs.
get_DgDx(j));
452 for (
int l=0; l<outArgs.
Np(); ++l)
453 if (!me_outArgs.supports(MEB::OUT_ARG_DgDp,j-g_offset_,l).none())
454 me_outArgs.set_DgDp(j-g_offset_, l, outArgs.
get_DgDp(j,l));
456 if (g_index_ >= 0 && outArgs.
get_g(1) != Teuchos::null)
457 me_outArgs.set_g(g_index_, outArgs.
get_g(1));
460 model_->evalModel(me_inArgs, me_outArgs);
463 if (outArgs.
supports(MEB::OUT_ARG_W_op) &&
464 outArgs.
get_W_op() != Teuchos::null &&
465 model_.ptr() != sens_solve_model_.ptr()) {
466 MEB::OutArgs<Scalar> sens_me_outArgs = sens_solve_model_->createOutArgs();
467 RCP< Thyra::LinearOpBase<Scalar> > op = outArgs.
get_W_op();
468 RCP< Thyra::MultiVectorLinearOp<Scalar> > mv_op =
470 sens_me_outArgs.set_W_op(mv_op->getNonconstLinearOp());
471 sens_solve_model_->evalModel(me_inArgs, sens_me_outArgs);
477 if (!use_dfdp_as_tangent_) {
478 if (dxdp != Teuchos::null && dfdp != Teuchos::null) {
479 if (my_dfdx_ == Teuchos::null)
480 my_dfdx_ = sens_residual_model_->create_W_op();
481 MEB::OutArgs<Scalar> meo = sens_residual_model_->createOutArgs();
482 meo.set_W_op(my_dfdx_);
483 if (me_inArgs.supports(MEB::IN_ARG_alpha))
484 me_inArgs.set_alpha(0.0);
485 if (me_inArgs.supports(MEB::IN_ARG_beta))
486 me_inArgs.set_beta(1.0);
487 if (me_inArgs.supports(MEB::IN_ARG_W_x_dot_dot_coeff))
488 me_inArgs.set_W_x_dot_dot_coeff(0.0);
489 sens_residual_model_->evalModel(me_inArgs, meo);
490 my_dfdx_->apply(
Thyra::NOTRANS, *dxdp, dfdp.ptr(), Scalar(1.0), Scalar(1.0));
492 if (dxdotdp != Teuchos::null && dfdp != Teuchos::null) {
493 if (my_dfdxdot_ == Teuchos::null)
494 my_dfdxdot_ = sens_residual_model_->create_W_op();
495 MEB::OutArgs<Scalar> meo = sens_residual_model_->createOutArgs();
496 meo.set_W_op(my_dfdxdot_);
497 if (me_inArgs.supports(MEB::IN_ARG_alpha))
498 me_inArgs.set_alpha(1.0);
499 if (me_inArgs.supports(MEB::IN_ARG_beta))
500 me_inArgs.set_beta(0.0);
501 if (me_inArgs.supports(MEB::IN_ARG_W_x_dot_dot_coeff))
502 me_inArgs.set_W_x_dot_dot_coeff(0.0);
503 sens_residual_model_->evalModel(me_inArgs, meo);
504 my_dfdxdot_->apply(
Thyra::NOTRANS, *dxdotdp, dfdp.ptr(), Scalar(1.0), Scalar(1.0));
506 if (dxdotdotdp != Teuchos::null && dfdp != Teuchos::null) {
507 if (my_dfdxdotdot_ == Teuchos::null)
508 my_dfdxdotdot_ = sens_residual_model_->create_W_op();
509 MEB::OutArgs<Scalar> meo = sens_residual_model_->createOutArgs();
510 meo.set_W_op(my_dfdxdotdot_);
511 if (me_inArgs.supports(MEB::IN_ARG_alpha))
512 me_inArgs.set_alpha(0.0);
513 if (me_inArgs.supports(MEB::IN_ARG_beta))
514 me_inArgs.set_beta(0.0);
515 if (me_inArgs.supports(MEB::IN_ARG_W_x_dot_dot_coeff))
516 me_inArgs.set_W_x_dot_dot_coeff(1.0);
517 sens_residual_model_->evalModel(me_inArgs, meo);
518 my_dfdxdotdot_->apply(
Thyra::NOTRANS, *dxdotdotdp, dfdp.ptr(), Scalar(1.0), Scalar(1.0));
523 if (g_index_ >= 0 && outArgs.
get_g(0) != Teuchos::null) {
524 MEB::OutArgs<Scalar> meo = model_->createOutArgs();
525 RCP<Thyra::MultiVectorBase<Scalar> > dgdp =
526 rcp_dynamic_cast<DMVPV>(outArgs.
get_g(0),
true)->getNonconstMultiVector();
527 RCP<Thyra::MultiVectorBase<Scalar> > dgdp_trans;
528 MEB::DerivativeSupport dgdp_support =
529 meo.supports(MEB::OUT_ARG_DgDp, g_index_, p_index_);
530 if (dgdp_support.supports(MEB::DERIV_MV_JACOBIAN_FORM))
531 meo.set_DgDp(g_index_, p_index_,
532 MEB::Derivative<Scalar>(dgdp, MEB::DERIV_MV_JACOBIAN_FORM));
533 else if (dgdp_support.supports(MEB::DERIV_MV_GRADIENT_FORM)) {
534 dgdp_trans = createMembers(model_->get_p_space(p_index_),
536 meo.set_DgDp(g_index_, p_index_,
537 MEB::Derivative<Scalar>(dgdp_trans, MEB::DERIV_MV_GRADIENT_FORM));
541 true, std::logic_error,
542 "Operator form of dg/dp not supported for reduced sensitivity");
544 if (!use_dgdp_as_tangent_ || dxdp == Teuchos::null) {
545 MEB::DerivativeSupport dgdx_support =
546 meo.supports(MEB::OUT_ARG_DgDx, g_index_);
547 if (dgdx_support.supports(MEB::DERIV_LINEAR_OP)) {
548 if (my_dgdx_ == Teuchos::null)
549 my_dgdx_ = model_->create_DgDx_op(g_index_);
550 meo.set_DgDx(g_index_, my_dgdx_);
552 else if (dgdx_support.supports(MEB::DERIV_MV_GRADIENT_FORM)) {
553 if (my_dgdx_mv_ == Teuchos::null)
554 my_dgdx_mv_ = createMembers(model_->get_g_space(g_index_), num_param_);
555 meo.set_DgDx(g_index_,
556 MEB::Derivative<Scalar>(my_dgdx_mv_,
557 MEB::DERIV_MV_GRADIENT_FORM));
561 true, std::logic_error,
562 "Jacobian form of dg/dx not supported for reduced sensitivity");
567 if (!use_dgdp_as_tangent_ && use_dfdp_as_tangent_ &&
568 dxdp_vec != Teuchos::null)
569 me_inArgs.set_p(x_tangent_index_, Teuchos::null);
570 if (!use_dgdp_as_tangent_ && use_dfdp_as_tangent_ &&
571 dxdotdp_vec != Teuchos::null)
572 me_inArgs.set_p(xdot_tangent_index_, Teuchos::null);
573 if (!use_dgdp_as_tangent_ && use_dfdp_as_tangent_ &&
574 dxdotdotdp_vec != Teuchos::null)
575 me_inArgs.set_p(xdotdot_tangent_index_, Teuchos::null);
578 if (use_dgdp_as_tangent_ && !use_dfdp_as_tangent_ &&
579 dxdp_vec != Teuchos::null)
580 me_inArgs.set_p(x_tangent_index_, dxdp_vec);
581 if (use_dgdp_as_tangent_ && !use_dfdp_as_tangent_ &&
582 dxdotdp_vec != Teuchos::null)
583 me_inArgs.set_p(xdot_tangent_index_, dxdotdp_vec);
584 if (use_dgdp_as_tangent_ && !use_dfdp_as_tangent_ &&
585 dxdotdotdp_vec != Teuchos::null)
586 me_inArgs.set_p(xdotdot_tangent_index_, dxdotdotdp_vec);
588 model_->evalModel(me_inArgs, meo);
591 if (dgdp_trans != Teuchos::null) {
594 for (
int i=0; i<num_param_; ++i)
595 for (
int j=0; j<num_response_; ++j)
596 dgdp_view(j,i) = dgdp_trans_view(i,j);
601 if (!use_dgdp_as_tangent_ && dxdp != Teuchos::null) {
602 MEB::DerivativeSupport dgdx_support =
603 meo.supports(MEB::OUT_ARG_DgDx, g_index_);
604 if (dgdx_support.supports(MEB::DERIV_LINEAR_OP))
605 my_dgdx_->apply(
Thyra::NOTRANS, *dxdp, dgdp.ptr(), Scalar(1.0), Scalar(1.0));
606 else if (dgdx_support.supports(MEB::DERIV_MV_GRADIENT_FORM))
607 my_dgdx_mv_->apply(
Thyra::TRANS, *dxdp, dgdp.ptr(), Scalar(1.0), Scalar(1.0));
610 true, std::logic_error,
611 "Jacobian form of dg/dx not supported for reduced sensitivity");
616 template<
class Scalar>
622 pl->
set<
bool>(
"Use DfDp as Tangent",
false);
623 pl->
set<
bool>(
"Use DgDp as Tangent",
false);
624 pl->
set<
int>(
"Sensitivity Parameter Index", 0);
625 pl->
set<
int>(
"Response Function Index", -1);
626 pl->
set<
int>(
"Sensitivity X Tangent Index", 1);
627 pl->
set<
int>(
"Sensitivity X-Dot Tangent Index", 2);
628 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