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>
33 const int p_index,
const int g_index,
const bool g_depends_on_p,
34 const bool f_depends_on_p,
const bool ic_depends_on_p,
35 const bool mass_matrix_is_identity)
37 state_integrator_(state_integrator),
38 adjoint_model_(adjoint_model),
39 adjoint_aux_model_(adjoint_aux_model),
40 adjoint_integrator_(adjoint_integrator),
41 solutionHistory_(solutionHistory),
44 g_depends_on_p_(g_depends_on_p),
45 f_depends_on_p_(f_depends_on_p),
46 ic_depends_on_p_(ic_depends_on_p),
47 mass_matrix_is_identity_(mass_matrix_is_identity),
52 "Error - IntegratorAdjointSensitivity(): Cannot use FSAL with\n"
53 " IntegratorAdjointSensitivity, because the state and adjoint\n"
54 " integrators require ModelEvaluator evaluation in the\n"
55 " constructor to make the initial conditions consistent.\n"
56 " For the adjoint integrator, this requires special construction\n"
57 " which has not been implemented yet.\n");
60 template <
class Scalar>
63 state_integrator_ = createIntegratorBasic<Scalar>();
64 adjoint_integrator_ = createIntegratorBasic<Scalar>();
68 template <
class Scalar>
71 const Scalar tfinal = state_integrator_->getTimeStepControl()->getFinalTime();
72 return advanceTime(tfinal);
75 template <
class Scalar>
78 TEMPUS_FUNC_TIME_MONITOR_DIFF(
79 "Tempus::IntegratorAdjointSensitivity::advanceTime()", TEMPUS_AS_AT);
82 using Teuchos::rcp_dynamic_cast;
84 using Thyra::createMember;
85 using Thyra::createMembers;
99 RCP<const SolutionHistory<Scalar>> state_solution_history =
100 state_integrator_->getSolutionHistory();
101 RCP<const SolutionState<Scalar>> initial_state = (*state_solution_history)[0];
104 bool state_status =
true;
106 TEMPUS_FUNC_TIME_MONITOR_DIFF(
107 "Tempus::IntegratorAdjointSensitivity::advanceTime::state",
110 state_status = state_integrator_->advanceTime(timeFinal);
117 adjoint_aux_model_->setFinalTime(state_integrator_->getTime());
120 adjoint_aux_model_->setForwardSolutionHistory(state_solution_history);
123 RCP<const VectorSpaceBase<Scalar>> g_space = model_->get_g_space(g_index_);
124 RCP<const VectorSpaceBase<Scalar>> x_space = model_->get_x_space();
125 const int num_g = g_space->dim();
126 RCP<MultiVectorBase<Scalar>> dgdx = createMembers(x_space, num_g);
127 MEB::InArgs<Scalar> inargs = model_->getNominalValues();
128 RCP<const SolutionState<Scalar>> state =
129 state_solution_history->getCurrentState();
130 inargs.set_t(state->getTime());
131 inargs.set_x(state->getX());
132 inargs.set_x_dot(state->getXDot());
133 MEB::OutArgs<Scalar> outargs = model_->createOutArgs();
134 MEB::OutArgs<Scalar> adj_outargs = adjoint_model_->createOutArgs();
135 outargs.set_DgDx(g_index_,
136 MEB::Derivative<Scalar>(dgdx, MEB::DERIV_MV_GRADIENT_FORM));
137 model_->evalModel(inargs, outargs);
138 outargs.set_DgDx(g_index_, MEB::Derivative<Scalar>());
144 RCP<DPV> adjoint_init = rcp_dynamic_cast<DPV>(
145 Thyra::createMember(adjoint_aux_model_->get_x_space()));
146 RCP<MultiVectorBase<Scalar>> adjoint_init_mv =
147 rcp_dynamic_cast<DMVPV>(adjoint_init->getNonconstVectorBlock(0))
148 ->getNonconstMultiVector();
149 assign(adjoint_init->getNonconstVectorBlock(1).ptr(),
151 if (mass_matrix_is_identity_)
152 assign(adjoint_init_mv.ptr(), *dgdx);
154 inargs.set_alpha(1.0);
155 inargs.set_beta(0.0);
156 RCP<LinearOpWithSolveBase<Scalar>> W;
157 if (adj_outargs.supports(MEB::OUT_ARG_W)) {
159 W = adjoint_model_->create_W();
160 adj_outargs.set_W(W);
161 adjoint_model_->evalModel(inargs, adj_outargs);
162 adj_outargs.set_W(Teuchos::null);
166 RCP<const LinearOpWithSolveFactoryBase<Scalar>> lowsfb =
167 adjoint_model_->get_W_factory();
169 "Adjoint ME must support W out-arg or provide "
170 "a W_factory for non-identity mass matrix");
173 RCP<LinearOpBase<Scalar>> W_op = adjoint_model_->create_W_op();
174 adj_outargs.set_W_op(W_op);
175 RCP<PreconditionerFactoryBase<Scalar>> prec_factory =
176 lowsfb->getPreconditionerFactory();
177 RCP<PreconditionerBase<Scalar>> W_prec;
178 if (prec_factory != Teuchos::null)
179 W_prec = prec_factory->createPrec();
180 else if (adj_outargs.supports(MEB::OUT_ARG_W_prec)) {
181 W_prec = adjoint_model_->create_W_prec();
182 adj_outargs.set_W_prec(W_prec);
184 adjoint_model_->evalModel(inargs, adj_outargs);
185 adj_outargs.set_W_op(Teuchos::null);
186 if (adj_outargs.supports(MEB::OUT_ARG_W_prec))
187 adj_outargs.set_W_prec(Teuchos::null);
190 W = lowsfb->createOp();
191 if (W_prec != Teuchos::null) {
192 if (prec_factory != Teuchos::null)
193 prec_factory->initializePrec(
194 Thyra::defaultLinearOpSource<Scalar>(W_op), W_prec.get());
195 Thyra::initializePreconditionedOp<Scalar>(*lowsfb, W_op, W_prec,
199 Thyra::initializeOp<Scalar>(*lowsfb, W_op, W.ptr());
202 W == Teuchos::null, std::logic_error,
203 "A null W has been encountered in "
204 "Tempus::IntegratorAdjointSensitivity::advanceTime!\n");
212 bool sens_status =
true;
214 TEMPUS_FUNC_TIME_MONITOR_DIFF(
215 "Tempus::IntegratorAdjointSensitivity::advanceTime::adjoint",
219 adjoint_integrator_->getTimeStepControl()->getInitTime();
220 adjoint_integrator_->initializeSolutionHistory(tinit, adjoint_init);
221 sens_status = adjoint_integrator_->advanceTime(timeFinal);
223 RCP<const SolutionHistory<Scalar>> adjoint_solution_history =
224 adjoint_integrator_->getSolutionHistory();
227 RCP<const VectorSpaceBase<Scalar>> p_space = model_->get_p_space(p_index_);
228 dgdp_ = createMembers(p_space, num_g);
229 if (g_depends_on_p_) {
230 MEB::DerivativeSupport dgdp_support =
231 outargs.supports(MEB::OUT_ARG_DgDp, g_index_, p_index_);
232 if (dgdp_support.supports(MEB::DERIV_MV_GRADIENT_FORM)) {
235 MEB::Derivative<Scalar>(dgdp_, MEB::DERIV_MV_GRADIENT_FORM));
236 model_->evalModel(inargs, outargs);
238 else if (dgdp_support.supports(MEB::DERIV_MV_JACOBIAN_FORM)) {
239 const int num_p = p_space->dim();
240 RCP<MultiVectorBase<Scalar>> dgdp_trans = createMembers(g_space, num_p);
243 MEB::Derivative<Scalar>(dgdp_trans, MEB::DERIV_MV_JACOBIAN_FORM));
244 model_->evalModel(inargs, outargs);
247 for (
int i = 0; i < num_p; ++i)
248 for (
int j = 0; j < num_g; ++j) dgdp_view(i, j) = dgdp_trans_view(j, i);
252 "Invalid dg/dp support");
253 outargs.set_DgDp(g_index_, p_index_, MEB::Derivative<Scalar>());
256 assign(dgdp_.ptr(), Scalar(0.0));
260 if (ic_depends_on_p_ && dxdp_init_ != Teuchos::null) {
261 RCP<const SolutionState<Scalar>> adjoint_state =
262 adjoint_solution_history->getCurrentState();
263 RCP<const VectorBase<Scalar>> adjoint_x =
264 rcp_dynamic_cast<
const DPV>(adjoint_state->getX())->getVectorBlock(0);
265 RCP<const MultiVectorBase<Scalar>> adjoint_mv =
266 rcp_dynamic_cast<
const DMVPV>(adjoint_x)->getMultiVector();
267 if (mass_matrix_is_identity_)
271 inargs.set_t(initial_state->getTime());
272 inargs.set_x(initial_state->getX());
273 inargs.set_x_dot(initial_state->getXDot());
274 inargs.set_alpha(1.0);
275 inargs.set_beta(0.0);
276 RCP<LinearOpBase<Scalar>> W_op = adjoint_model_->create_W_op();
277 adj_outargs.set_W_op(W_op);
278 adjoint_model_->evalModel(inargs, adj_outargs);
279 adj_outargs.set_W_op(Teuchos::null);
280 RCP<MultiVectorBase<Scalar>> tmp = createMembers(x_space, num_g);
291 if (f_depends_on_p_) {
292 RCP<const SolutionState<Scalar>> adjoint_state =
293 adjoint_solution_history->getCurrentState();
294 RCP<const VectorBase<Scalar>> z =
295 rcp_dynamic_cast<
const DPV>(adjoint_state->getX())->getVectorBlock(1);
296 RCP<const MultiVectorBase<Scalar>> z_mv =
297 rcp_dynamic_cast<
const DMVPV>(z)->getMultiVector();
298 Thyra::V_VmV(dgdp_.ptr(), *dgdp_, *z_mv);
301 buildSolutionHistory(state_solution_history, adjoint_solution_history);
303 return state_status && sens_status;
306 template <
class Scalar>
309 return solutionHistory_->getCurrentTime();
312 template <
class Scalar>
315 return solutionHistory_->getCurrentIndex();
318 template <
class Scalar>
321 Status state_status = state_integrator_->getStatus();
322 Status sens_status = adjoint_integrator_->getStatus();
328 template <
class Scalar>
331 state_integrator_->setStatus(st);
332 adjoint_integrator_->setStatus(st);
335 template <
class Scalar>
339 return state_integrator_->getStepper();
342 template <
class Scalar>
346 return solutionHistory_;
349 template <
class Scalar>
353 return state_integrator_->getSolutionHistory();
356 template <
class Scalar>
360 return adjoint_integrator_->getSolutionHistory();
363 template <
class Scalar>
367 return solutionHistory_;
370 template <
class Scalar>
374 return state_integrator_->getTimeStepControl();
377 template <
class Scalar>
381 return state_integrator_->getNonConstTimeStepControl();
384 template <
class Scalar>
388 return state_integrator_->getNonConstTimeStepControl();
391 template <
class Scalar>
395 return adjoint_integrator_->getNonConstTimeStepControl();
398 template <
class Scalar>
407 state_integrator_->initializeSolutionHistory(t0, x0, xdot0, xdotdot0);
411 template <
class Scalar>
415 return state_integrator_->getObserver();
418 template <
class Scalar>
422 state_integrator_->setObserver(obs);
428 template <
class Scalar>
431 state_integrator_->initialize();
432 adjoint_integrator_->initialize();
435 template <
class Scalar>
439 return state_integrator_->getX();
442 template <
class Scalar>
446 return state_integrator_->getXDot();
449 template <
class Scalar>
453 return state_integrator_->getXDotDot();
456 template <
class Scalar>
461 using Teuchos::rcp_dynamic_cast;
464 RCP<const DPV> pv = rcp_dynamic_cast<
const DPV>(adjoint_integrator_->getX());
465 RCP<const DMVPV> mvpv = rcp_dynamic_cast<
const DMVPV>(pv->getVectorBlock(0));
469 template <
class Scalar>
474 using Teuchos::rcp_dynamic_cast;
478 rcp_dynamic_cast<
const DPV>(adjoint_integrator_->getXDot());
479 RCP<const DMVPV> mvpv = rcp_dynamic_cast<
const DMVPV>(pv->getVectorBlock(0));
483 template <
class Scalar>
488 using Teuchos::rcp_dynamic_cast;
492 rcp_dynamic_cast<
const DPV>(adjoint_integrator_->getXDotDot());
493 RCP<const DMVPV> mvpv = rcp_dynamic_cast<
const DMVPV>(pv->getVectorBlock(0));
497 template <
class Scalar>
504 template <
class Scalar>
507 std::string name =
"Tempus::IntegratorAdjointSensitivity";
511 template <
class Scalar>
515 auto l_out = Teuchos::fancyOStream(out.
getOStream());
517 l_out->setOutputToRootOnly(0);
519 *l_out << description() <<
"::describe" << std::endl;
520 state_integrator_->describe(*l_out, verbLevel);
521 adjoint_integrator_->describe(*l_out, verbLevel);
524 template <
class Scalar>
530 template <
class Scalar>
541 RCP<ParameterList> spl = Teuchos::parameterList();
542 if (inputPL != Teuchos::null) {
543 *spl = inputPL->
sublist(
"Sensitivities");
545 if (spl->isParameter(
"Response Depends on Parameters"))
546 spl->
remove(
"Response Depends on Parameters");
547 if (spl->isParameter(
"Residual Depends on Parameters"))
548 spl->remove(
"Residual Depends on Parameters");
549 if (spl->isParameter(
"IC Depends on Parameters"))
550 spl->remove(
"IC Depends on Parameters");
552 const Scalar tinit = state_integrator_->getTimeStepControl()->getInitTime();
553 const Scalar tfinal = state_integrator_->getTimeStepControl()->getFinalTime();
555 model, adjoint_model, tinit, tfinal, spl));
558 template <
class Scalar>
566 using Teuchos::rcp_dynamic_cast;
568 using Thyra::createMembers;
570 using Thyra::multiVectorProductVector;
576 RCP<const VectorSpaceBase<Scalar>> x_space = model_->get_x_space();
577 RCP<const VectorSpaceBase<Scalar>> adjoint_space =
578 rcp_dynamic_cast<
const DPVS>(adjoint_aux_model_->get_x_space())
582 spaces[1] = adjoint_space;
583 RCP<const DPVS> prod_space = Thyra::productVectorSpace(spaces());
585 int num_states = state_solution_history->getNumStates();
586 const Scalar t_init = state_integrator_->getTimeStepControl()->getInitTime();
587 const Scalar t_final = state_integrator_->getTime();
588 for (
int i = 0; i < num_states; ++i) {
589 RCP<const SolutionState<Scalar>> forward_state =
590 (*state_solution_history)[i];
591 RCP<const SolutionState<Scalar>> adjoint_state =
592 adjoint_solution_history->findState(t_final + t_init -
593 forward_state->getTime());
596 RCP<DPV> x = Thyra::defaultProductVector(prod_space);
597 RCP<const VectorBase<Scalar>> adjoint_x =
598 rcp_dynamic_cast<
const DPV>(adjoint_state->getX())->getVectorBlock(0);
599 assign(x->getNonconstVectorBlock(0).ptr(), *(forward_state->getX()));
600 assign(x->getNonconstVectorBlock(1).ptr(), *(adjoint_x));
601 RCP<VectorBase<Scalar>> x_b = x;
604 RCP<DPV> x_dot = Thyra::defaultProductVector(prod_space);
605 RCP<const VectorBase<Scalar>> adjoint_x_dot =
606 rcp_dynamic_cast<
const DPV>(adjoint_state->getXDot())
608 assign(x_dot->getNonconstVectorBlock(0).ptr(), *(forward_state->getXDot()));
609 assign(x_dot->getNonconstVectorBlock(1).ptr(), *(adjoint_x_dot));
610 RCP<VectorBase<Scalar>> x_dot_b = x_dot;
614 if (forward_state->getXDotDot() != Teuchos::null) {
615 x_dot_dot = Thyra::defaultProductVector(prod_space);
616 RCP<const VectorBase<Scalar>> adjoint_x_dot_dot =
617 rcp_dynamic_cast<
const DPV>(adjoint_state->getXDotDot())
619 assign(x_dot_dot->getNonconstVectorBlock(0).ptr(),
620 *(forward_state->getXDotDot()));
621 assign(x_dot_dot->getNonconstVectorBlock(1).ptr(), *(adjoint_x_dot_dot));
623 RCP<VectorBase<Scalar>> x_dot_dot_b = x_dot_dot;
625 RCP<SolutionState<Scalar>> prod_state = forward_state->clone();
626 prod_state->setX(x_b);
627 prod_state->setXDot(x_dot_b);
628 prod_state->setXDotDot(x_dot_dot_b);
629 prod_state->setPhysicsState(Teuchos::null);
630 solutionHistory_->addState(prod_state);
635 template <
class Scalar>
644 if (inputPL != Teuchos::null) *spl = inputPL->
sublist(
"Sensitivities");
646 int p_index = spl->
get<
int>(
"Sensitivity Parameter Index", 0);
647 int g_index = spl->
get<
int>(
"Response Function Index", 0);
648 bool g_depends_on_p = spl->
get<
bool>(
"Response Depends on Parameters",
true);
649 bool f_depends_on_p = spl->
get<
bool>(
"Residual Depends on Parameters",
true);
650 bool ic_depends_on_p = spl->
get<
bool>(
"IC Depends on Parameters",
true);
651 bool mass_matrix_is_identity =
652 spl->
get<
bool>(
"Mass Matrix Is Identity",
false);
654 auto state_integrator = createIntegratorBasic<Scalar>(inputPL, model);
657 if (spl->
isParameter(
"Response Depends on Parameters"))
658 spl->
remove(
"Response Depends on Parameters");
659 if (spl->
isParameter(
"Residual Depends on Parameters"))
660 spl->
remove(
"Residual Depends on Parameters");
662 spl->
remove(
"IC Depends on Parameters");
664 const Scalar tinit = state_integrator->getTimeStepControl()->getInitTime();
665 const Scalar tfinal = state_integrator->getTimeStepControl()->getFinalTime();
671 if (adjoint_model == Teuchos::null)
674 auto adjoint_aux_model =
676 model, adjt_model, tinit, tfinal, spl));
680 auto integrator_name = inputPL->get<std::string>(
"Integrator Name");
681 auto integratorPL = Teuchos::sublist(inputPL, integrator_name,
true);
682 auto shPL = Teuchos::sublist(integratorPL,
"Solution History",
true);
683 auto combined_solution_History = createSolutionHistoryPL<Scalar>(shPL);
685 auto adjoint_integrator =
686 createIntegratorBasic<Scalar>(inputPL, adjoint_aux_model);
690 model, state_integrator, adjt_model, adjoint_aux_model,
691 adjoint_integrator, combined_solution_History, p_index, g_index,
692 g_depends_on_p, f_depends_on_p, ic_depends_on_p,
693 mass_matrix_is_identity));
699 template <
class Scalar>
709 #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
virtual void setObserver(Teuchos::RCP< IntegratorObserver< Scalar >> obs=Teuchos::null)
Set the Observer.
T & get(const std::string &name, T def_value)
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)
void buildSolutionHistory(const Teuchos::RCP< const SolutionHistory< Scalar >> &state_solution_history, const Teuchos::RCP< const SolutionHistory< Scalar >> &adjoint_solution_history)
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.
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::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.
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.
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)
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 Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getY() const
Get the current adjoint solution, y.
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)
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()
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...