9 #ifndef Tempus_IntegratorAdjointSensitivity_impl_hpp
10 #define Tempus_IntegratorAdjointSensitivity_impl_hpp
13 #include "Thyra_DefaultMultiVectorProductVectorSpace.hpp"
14 #include "Thyra_DefaultMultiVectorProductVector.hpp"
15 #include "Thyra_DefaultProductVectorSpace.hpp"
16 #include "Thyra_DefaultProductVector.hpp"
17 #include "Thyra_VectorStdOps.hpp"
18 #include "Thyra_MultiVectorStdOps.hpp"
20 #include "Thyra_LinearOpWithSolveFactoryHelpers.hpp"
24 template <
class Scalar>
34 const bool g_depends_on_p,
35 const bool f_depends_on_p,
36 const bool ic_depends_on_p,
37 const bool mass_matrix_is_identity)
39 , state_integrator_(state_integrator)
40 , adjoint_model_(adjoint_model)
41 , adjoint_aux_model_(adjoint_aux_model)
42 , adjoint_integrator_(adjoint_integrator)
43 , solutionHistory_(solutionHistory)
46 , g_depends_on_p_(g_depends_on_p)
47 , f_depends_on_p_(f_depends_on_p)
48 , ic_depends_on_p_(ic_depends_on_p)
49 , mass_matrix_is_identity_(mass_matrix_is_identity)
55 "Error - IntegratorAdjointSensitivity(): Cannot use FSAL with\n"
56 " IntegratorAdjointSensitivity, because the state and adjoint\n"
57 " integrators require ModelEvaluator evaluation in the\n"
58 " constructor to make the initial conditions consistent.\n"
59 " For the adjoint integrator, this requires special construction\n"
60 " which has not been implemented yet.\n");
63 template<
class Scalar>
67 state_integrator_ = createIntegratorBasic<Scalar>();
68 adjoint_integrator_ = createIntegratorBasic<Scalar>();
72 template<
class Scalar>
78 state_integrator_->getTimeStepControl()->getFinalTime();
79 return advanceTime(tfinal);
82 template<
class Scalar>
87 TEMPUS_FUNC_TIME_MONITOR_DIFF(
"Tempus::IntegratorAdjointSensitivity::advanceTime()", TEMPUS_AS_AT);
90 using Teuchos::rcp_dynamic_cast;
99 using Thyra::createMember;
100 using Thyra::createMembers;
107 RCP<const SolutionHistory<Scalar> > state_solution_history =
108 state_integrator_->getSolutionHistory();
109 RCP<const SolutionState<Scalar> > initial_state =
110 (*state_solution_history)[0];
113 bool state_status =
true;
115 TEMPUS_FUNC_TIME_MONITOR_DIFF(
"Tempus::IntegratorAdjointSensitivity::advanceTime::state", TEMPUS_AS_AT_FWD);
117 state_status = state_integrator_->advanceTime(timeFinal);
124 adjoint_aux_model_->setFinalTime(state_integrator_->getTime());
127 adjoint_aux_model_->setForwardSolutionHistory(state_solution_history);
130 RCP<const VectorSpaceBase<Scalar> > g_space = model_->get_g_space(g_index_);
131 RCP<const VectorSpaceBase<Scalar> > x_space = model_->get_x_space();
132 const int num_g = g_space->dim();
133 RCP<MultiVectorBase<Scalar> > dgdx = createMembers(x_space, num_g);
134 MEB::InArgs<Scalar> inargs = model_->getNominalValues();
135 RCP<const SolutionState<Scalar> > state =
136 state_solution_history->getCurrentState();
137 inargs.set_t(state->getTime());
138 inargs.set_x(state->getX());
139 inargs.set_x_dot(state->getXDot());
140 MEB::OutArgs<Scalar> outargs = model_->createOutArgs();
141 MEB::OutArgs<Scalar> adj_outargs = adjoint_model_->createOutArgs();
142 outargs.set_DgDx(g_index_,
143 MEB::Derivative<Scalar>(dgdx, MEB::DERIV_MV_GRADIENT_FORM));
144 model_->evalModel(inargs, outargs);
145 outargs.set_DgDx(g_index_, MEB::Derivative<Scalar>());
151 RCP<DPV> adjoint_init =
152 rcp_dynamic_cast<DPV>(Thyra::createMember(adjoint_aux_model_->get_x_space()));
153 RCP<MultiVectorBase<Scalar> > adjoint_init_mv =
154 rcp_dynamic_cast<DMVPV>(adjoint_init->getNonconstVectorBlock(0))->getNonconstMultiVector();
155 assign(adjoint_init->getNonconstVectorBlock(1).ptr(),
157 if (mass_matrix_is_identity_)
158 assign(adjoint_init_mv.ptr(), *dgdx);
160 inargs.set_alpha(1.0);
161 inargs.set_beta(0.0);
162 RCP<LinearOpWithSolveBase<Scalar> > W;
163 if (adj_outargs.supports(MEB::OUT_ARG_W)) {
165 W = adjoint_model_->create_W();
166 adj_outargs.set_W(W);
167 adjoint_model_->evalModel(inargs, adj_outargs);
168 adj_outargs.set_W(Teuchos::null);
172 RCP<const LinearOpWithSolveFactoryBase<Scalar> > lowsfb =
173 adjoint_model_->get_W_factory();
175 lowsfb == Teuchos::null, std::logic_error,
176 "Adjoint ME must support W out-arg or provide a W_factory for non-identity mass matrix");
179 RCP<LinearOpBase<Scalar> > W_op = adjoint_model_->create_W_op();
180 adj_outargs.set_W_op(W_op);
181 RCP<PreconditionerFactoryBase<Scalar> > prec_factory =
182 lowsfb->getPreconditionerFactory();
183 RCP<PreconditionerBase<Scalar> > W_prec;
184 if (prec_factory != Teuchos::null)
185 W_prec = prec_factory->createPrec();
186 else if (adj_outargs.supports(MEB::OUT_ARG_W_prec)) {
187 W_prec = adjoint_model_->create_W_prec();
188 adj_outargs.set_W_prec(W_prec);
190 adjoint_model_->evalModel(inargs, adj_outargs);
191 adj_outargs.set_W_op(Teuchos::null);
192 if (adj_outargs.supports(MEB::OUT_ARG_W_prec))
193 adj_outargs.set_W_prec(Teuchos::null);
196 W = lowsfb->createOp();
197 if (W_prec != Teuchos::null) {
198 if (prec_factory != Teuchos::null)
199 prec_factory->initializePrec(
200 Thyra::defaultLinearOpSource<Scalar>(W_op), W_prec.get());
201 Thyra::initializePreconditionedOp<Scalar>(
202 *lowsfb, W_op, W_prec, W.ptr());
205 Thyra::initializeOp<Scalar>(*lowsfb, W_op, W.ptr());
208 W == Teuchos::null, std::logic_error,
209 "A null W has been encountered in Tempus::IntegratorAdjointSensitivity::advanceTime!\n");
217 bool sens_status =
true;
219 TEMPUS_FUNC_TIME_MONITOR_DIFF(
"Tempus::IntegratorAdjointSensitivity::advanceTime::adjoint", TEMPUS_AS_AT_ADJ);
221 const Scalar tinit = adjoint_integrator_->getTimeStepControl()->getInitTime();
222 adjoint_integrator_->initializeSolutionHistory(tinit, adjoint_init);
223 sens_status = adjoint_integrator_->advanceTime(timeFinal);
225 RCP<const SolutionHistory<Scalar> > adjoint_solution_history =
226 adjoint_integrator_->getSolutionHistory();
229 RCP<const VectorSpaceBase<Scalar> > p_space = model_->get_p_space(p_index_);
230 dgdp_ = createMembers(p_space, num_g);
231 if (g_depends_on_p_) {
232 MEB::DerivativeSupport dgdp_support =
233 outargs.supports(MEB::OUT_ARG_DgDp, g_index_, p_index_);
234 if (dgdp_support.supports(MEB::DERIV_MV_GRADIENT_FORM)) {
235 outargs.set_DgDp(g_index_, p_index_,
236 MEB::Derivative<Scalar>(dgdp_,
237 MEB::DERIV_MV_GRADIENT_FORM));
238 model_->evalModel(inargs, outargs);
240 else if (dgdp_support.supports(MEB::DERIV_MV_JACOBIAN_FORM)) {
241 const int num_p = p_space->dim();
242 RCP<MultiVectorBase<Scalar> > dgdp_trans =
243 createMembers(g_space, num_p);
244 outargs.set_DgDp(g_index_, p_index_,
245 MEB::Derivative<Scalar>(dgdp_trans,
246 MEB::DERIV_MV_JACOBIAN_FORM));
247 model_->evalModel(inargs, outargs);
250 for (
int i=0; i<num_p; ++i)
251 for (
int j=0; j<num_g; ++j)
252 dgdp_view(i,j) = dgdp_trans_view(j,i);
256 "Invalid dg/dp support");
257 outargs.set_DgDp(g_index_, p_index_, MEB::Derivative<Scalar>());
260 assign(dgdp_.ptr(), Scalar(0.0));
264 if (ic_depends_on_p_ && dxdp_init_ != Teuchos::null) {
265 RCP<const SolutionState<Scalar> > adjoint_state =
266 adjoint_solution_history->getCurrentState();
267 RCP<const VectorBase<Scalar> > adjoint_x =
268 rcp_dynamic_cast<
const DPV>(adjoint_state->getX())->getVectorBlock(0);
269 RCP<const MultiVectorBase<Scalar> > adjoint_mv =
270 rcp_dynamic_cast<
const DMVPV>(adjoint_x)->getMultiVector();
271 if (mass_matrix_is_identity_)
275 inargs.set_t(initial_state->getTime());
276 inargs.set_x(initial_state->getX());
277 inargs.set_x_dot(initial_state->getXDot());
278 inargs.set_alpha(1.0);
279 inargs.set_beta(0.0);
280 RCP<LinearOpBase<Scalar> > W_op = adjoint_model_->create_W_op();
281 adj_outargs.set_W_op(W_op);
282 adjoint_model_->evalModel(inargs, adj_outargs);
283 adj_outargs.set_W_op(Teuchos::null);
284 RCP<MultiVectorBase<Scalar> > tmp = createMembers(x_space, num_g);
295 if (f_depends_on_p_) {
296 RCP<const SolutionState<Scalar> > adjoint_state =
297 adjoint_solution_history->getCurrentState();
298 RCP<const VectorBase<Scalar> > z =
299 rcp_dynamic_cast<
const DPV>(adjoint_state->getX())->getVectorBlock(1);
300 RCP<const MultiVectorBase<Scalar> > z_mv =
301 rcp_dynamic_cast<
const DMVPV>(z)->getMultiVector();
302 Thyra::V_VmV(dgdp_.ptr(), *dgdp_, *z_mv);
305 buildSolutionHistory(state_solution_history, adjoint_solution_history);
307 return state_status && sens_status;
310 template<
class Scalar>
315 return solutionHistory_->getCurrentTime();
318 template<
class Scalar>
323 return solutionHistory_->getCurrentIndex();
326 template<
class Scalar>
331 Status state_status = state_integrator_->getStatus();
332 Status sens_status = adjoint_integrator_->getStatus();
340 template <
class Scalar>
342 state_integrator_->setStatus(st);
343 adjoint_integrator_->setStatus(st);
346 template<
class Scalar>
351 return state_integrator_->getStepper();
354 template<
class Scalar>
359 return solutionHistory_;
362 template<
class Scalar>
367 return state_integrator_->getSolutionHistory();
370 template<
class Scalar>
375 return adjoint_integrator_->getSolutionHistory();
378 template<
class Scalar>
383 return solutionHistory_;
386 template<
class Scalar>
391 return state_integrator_->getTimeStepControl();
394 template<
class Scalar>
399 return state_integrator_->getNonConstTimeStepControl();
402 template<
class Scalar>
407 return state_integrator_->getNonConstTimeStepControl();
410 template<
class Scalar>
415 return adjoint_integrator_->getNonConstTimeStepControl();
418 template<
class Scalar>
428 state_integrator_->initializeSolutionHistory(t0, x0, xdot0, xdotdot0);
432 template<
class Scalar>
437 return state_integrator_->getObserver();
440 template<
class Scalar>
444 state_integrator_->setObserver(obs);
450 template<
class Scalar>
454 state_integrator_->initialize();
455 adjoint_integrator_->initialize();
458 template<
class Scalar>
463 return state_integrator_->getX();
466 template<
class Scalar>
471 return state_integrator_->getXDot();
474 template<
class Scalar>
479 return state_integrator_->getXDotDot();
482 template<
class Scalar>
488 using Teuchos::rcp_dynamic_cast;
492 rcp_dynamic_cast<
const DPV>(adjoint_integrator_->getX());
493 RCP<const DMVPV> mvpv =
494 rcp_dynamic_cast<
const DMVPV>(pv->getVectorBlock(0));
498 template<
class Scalar>
504 using Teuchos::rcp_dynamic_cast;
508 rcp_dynamic_cast<
const DPV>(adjoint_integrator_->getXDot());
509 RCP<const DMVPV> mvpv =
510 rcp_dynamic_cast<
const DMVPV>(pv->getVectorBlock(0));
514 template<
class Scalar>
520 using Teuchos::rcp_dynamic_cast;
524 rcp_dynamic_cast<
const DPV>(adjoint_integrator_->getXDotDot());
525 RCP<const DMVPV> mvpv =
526 rcp_dynamic_cast<
const DMVPV>(pv->getVectorBlock(0));
530 template<
class Scalar>
538 template<
class Scalar>
543 std::string name =
"Tempus::IntegratorAdjointSensitivity";
547 template<
class Scalar>
554 auto l_out = Teuchos::fancyOStream( out.
getOStream() );
556 l_out->setOutputToRootOnly(0);
558 *l_out << description() <<
"::describe" << std::endl;
559 state_integrator_->describe(*l_out, verbLevel);
560 adjoint_integrator_->describe(*l_out, verbLevel);
563 template<
class Scalar>
571 template <
class Scalar>
583 RCP<ParameterList> spl = Teuchos::parameterList();
584 if (inputPL != Teuchos::null)
586 *spl = inputPL->
sublist(
"Sensitivities");
588 if (spl->isParameter(
"Response Depends on Parameters"))
589 spl->
remove(
"Response Depends on Parameters");
590 if (spl->isParameter(
"Residual Depends on Parameters"))
591 spl->remove(
"Residual Depends on Parameters");
592 if (spl->isParameter(
"IC Depends on Parameters"))
593 spl->remove(
"IC Depends on Parameters");
595 const Scalar tinit = state_integrator_->getTimeStepControl()->getInitTime();
596 const Scalar tfinal = state_integrator_->getTimeStepControl()->getFinalTime();
598 model, adjoint_model, tinit, tfinal, spl));
601 template<
class Scalar>
610 using Teuchos::rcp_dynamic_cast;
615 using Thyra::createMembers;
616 using Thyra::multiVectorProductVector;
621 RCP<const VectorSpaceBase<Scalar> > x_space = model_->get_x_space();
622 RCP<const VectorSpaceBase<Scalar> > adjoint_space =
623 rcp_dynamic_cast<
const DPVS>(adjoint_aux_model_->get_x_space())->getBlock(0);
626 spaces[1] = adjoint_space;
627 RCP<const DPVS > prod_space = Thyra::productVectorSpace(spaces());
629 int num_states = state_solution_history->getNumStates();
630 const Scalar t_init = state_integrator_->getTimeStepControl()->getInitTime();
631 const Scalar t_final = state_integrator_->getTime();
632 for (
int i=0; i<num_states; ++i) {
633 RCP<const SolutionState<Scalar> > forward_state =
634 (*state_solution_history)[i];
635 RCP<const SolutionState<Scalar> > adjoint_state =
636 adjoint_solution_history->findState(t_final+t_init-forward_state->getTime());
639 RCP<DPV> x = Thyra::defaultProductVector(prod_space);
640 RCP<const VectorBase<Scalar> > adjoint_x =
641 rcp_dynamic_cast<
const DPV>(adjoint_state->getX())->getVectorBlock(0);
642 assign(x->getNonconstVectorBlock(0).ptr(), *(forward_state->getX()));
643 assign(x->getNonconstVectorBlock(1).ptr(), *(adjoint_x));
644 RCP<VectorBase<Scalar> > x_b = x;
647 RCP<DPV> x_dot = Thyra::defaultProductVector(prod_space);
648 RCP<const VectorBase<Scalar> > adjoint_x_dot =
649 rcp_dynamic_cast<
const DPV>(adjoint_state->getXDot())->getVectorBlock(0);
650 assign(x_dot->getNonconstVectorBlock(0).ptr(), *(forward_state->getXDot()));
651 assign(x_dot->getNonconstVectorBlock(1).ptr(), *(adjoint_x_dot));
652 RCP<VectorBase<Scalar> > x_dot_b = x_dot;
656 if (forward_state->getXDotDot() != Teuchos::null) {
657 x_dot_dot = Thyra::defaultProductVector(prod_space);
658 RCP<const VectorBase<Scalar> > adjoint_x_dot_dot =
659 rcp_dynamic_cast<
const DPV>(
660 adjoint_state->getXDotDot())->getVectorBlock(0);
661 assign(x_dot_dot->getNonconstVectorBlock(0).ptr(),
662 *(forward_state->getXDotDot()));
663 assign(x_dot_dot->getNonconstVectorBlock(1).ptr(),
664 *(adjoint_x_dot_dot));
666 RCP<VectorBase<Scalar> > x_dot_dot_b = x_dot_dot;
668 RCP<SolutionState<Scalar> > prod_state = forward_state->clone();
669 prod_state->setX(x_b);
670 prod_state->setXDot(x_dot_b);
671 prod_state->setXDotDot(x_dot_dot_b);
672 prod_state->setPhysicsState(Teuchos::null);
673 solutionHistory_->addState(prod_state);
678 template<
class Scalar>
687 if (inputPL != Teuchos::null)
688 *spl = inputPL->
sublist(
"Sensitivities");
690 int p_index = spl->
get<
int>(
"Sensitivity Parameter Index", 0);
691 int g_index = spl->
get<
int>(
"Response Function Index", 0);
692 bool g_depends_on_p = spl->
get<
bool>(
"Response Depends on Parameters",
true);
693 bool f_depends_on_p = spl->
get<
bool>(
"Residual Depends on Parameters",
true);
694 bool ic_depends_on_p = spl->
get<
bool>(
"IC Depends on Parameters",
true);
695 bool mass_matrix_is_identity = spl->
get<
bool>(
"Mass Matrix Is Identity",
false);
697 auto state_integrator = createIntegratorBasic<Scalar>(inputPL, model);
700 if (spl->
isParameter(
"Response Depends on Parameters"))
701 spl->
remove(
"Response Depends on Parameters");
702 if (spl->
isParameter(
"Residual Depends on Parameters"))
703 spl->
remove(
"Residual Depends on Parameters");
705 spl->
remove(
"IC Depends on Parameters");
707 const Scalar tinit = state_integrator->getTimeStepControl()->getInitTime();
708 const Scalar tfinal = state_integrator->getTimeStepControl()->getFinalTime();
713 if (adjoint_model == Teuchos::null)
720 auto integrator_name = inputPL->get<std::string>(
"Integrator Name");
721 auto integratorPL = Teuchos::sublist(inputPL, integrator_name,
true);
722 auto shPL = Teuchos::sublist(integratorPL,
"Solution History",
true);
723 auto combined_solution_History = createSolutionHistoryPL<Scalar>(shPL);
725 auto adjoint_integrator = createIntegratorBasic<Scalar>(inputPL, adjoint_aux_model);
728 model, state_integrator, adjt_model, adjoint_aux_model, adjoint_integrator, combined_solution_History, p_index, g_index, g_depends_on_p,
729 f_depends_on_p, ic_depends_on_p, mass_matrix_is_identity));
735 template<
class Scalar>
745 #endif // Tempus_IntegratorAdjointSensitivity_impl_hpp
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const override
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getYDot() const
Get the current time derivative of the adjoint solution, ydot.
SensitivityStepMode getStepMode() const
What mode the current time integration step is in.
virtual Teuchos::RCP< TimeStepControl< Scalar > > getNonConstTimeStepControl() override
T & get(const std::string &name, T def_value)
Teuchos::RCP< IntegratorAdjointSensitivity< Scalar > > createIntegratorAdjointSensitivity(Teuchos::RCP< Teuchos::ParameterList > pList, const Teuchos::RCP< Thyra::ModelEvaluator< Scalar >> &model, const Teuchos::RCP< Thyra::ModelEvaluator< Scalar > > &adjoint_model=Teuchos::null)
Nonmember constructor.
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getXDotDot() const
Get the current second time derivative of the solution, xdotdot.
virtual void setStatus(const Status st) override
Set Status.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
virtual Teuchos::RCP< SolutionHistory< Scalar > > getNonConstSolutionHistory() override
Get the SolutionHistory.
ModelEvaluator for forming adjoint sensitivity equations.
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getXDot() const
Get the current time derivative of the solution, xdot.
void buildSolutionHistory(const Teuchos::RCP< const SolutionHistory< Scalar > > &state_solution_history, const Teuchos::RCP< const SolutionHistory< Scalar > > &adjoint_solution_history)
virtual Teuchos::RCP< IntegratorObserver< Scalar > > getObserver()
Get the Observer.
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getX() const
Get the current solution, x.
bool isParameter(const std::string &name) const
bool remove(std::string const &name, bool throwIfNotExists=true)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
std::string description() const override
virtual Teuchos::RCP< const SolutionHistory< Scalar > > getSolutionHistory() const override
Get the SolutionHistory.
Status
Status for the Integrator, the Stepper and the SolutionState.
IntegratorObserver class for time integrators.
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getYDotDot() const
Get the current second time derivative of the adjoint solution, ydotdot.
virtual void initializeSolutionHistory(Scalar t0, Teuchos::RCP< const Thyra::VectorBase< Scalar > > x0, Teuchos::RCP< const Thyra::VectorBase< Scalar > > xdot0=Teuchos::null, Teuchos::RCP< const Thyra::VectorBase< Scalar > > xdotdot0=Teuchos::null, Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > DxDp0=Teuchos::null, Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > DxdotDp0=Teuchos::null, Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > DxdotdotDp0=Teuchos::null)
Set the initial state from Thyra::VectorBase(s)
virtual Teuchos::RCP< const TimeStepControl< Scalar > > getTimeStepControl() const override
Get the TimeStepControl.
RCP< const MultiVectorBase< Scalar > > getMultiVector() const
SolutionHistory is basically a container of SolutionStates. SolutionHistory maintains a collection of...
virtual bool advanceTime()
Advance the solution to timeMax, and return true if successful.
Teuchos::RCP< TimeStepControl< Scalar > > getSensNonConstTimeStepControl()
Teuchos::RCP< const SolutionHistory< Scalar > > getSensSolutionHistory() const
virtual void setObserver(Teuchos::RCP< IntegratorObserver< Scalar > > obs=Teuchos::null)
Set the Observer.
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getY() const
Get the current adjoint solution, y.
Teuchos::RCP< const SolutionHistory< Scalar > > getStateSolutionHistory() const
RCP< ImplicitAdjointModelEvaluator< Scalar > > implicitAdjointModelEvaluator(const RCP< const ModelEvaluator< Scalar > > &model)
RCP< std::basic_ostream< char_type, traits_type > > getOStream()
Time integrator suitable for adjoint sensitivity analysis.
virtual int getIndex() const override
Get current index.
Teuchos::RCP< TimeStepControl< Scalar > > getStateNonConstTimeStepControl()
Teuchos::RCP< AdjointAuxSensitivityModelEvaluator< Scalar > > createAdjointModel(const Teuchos::RCP< Thyra::ModelEvaluator< Scalar > > &model, const Teuchos::RCP< Thyra::ModelEvaluator< Scalar > > &adjoint_model, const Teuchos::RCP< Teuchos::ParameterList > &inputPL)
ParameterList & sublist(const std::string &name, bool mustAlreadyExist=false, const std::string &docString="")
virtual Scalar getTime() const override
Get current time.
virtual Status getStatus() const override
Get Status.
virtual void initialize()
Initializes the Integrator after set* function calls.
virtual Teuchos::RCP< Stepper< Scalar > > getStepper() const override
Get the Stepper.
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getDgDp() const
Return adjoint sensitivity stored in gradient format.
IntegratorAdjointSensitivity()
Constructor that requires a subsequent setParameterList, setStepper, and initialize calls...