9#ifndef Tempus_IntegratorPseudoTransientAdjointSensitivity_impl_hpp
10#define Tempus_IntegratorPseudoTransientAdjointSensitivity_impl_hpp
12#include "Thyra_DefaultMultiVectorProductVector.hpp"
13#include "Thyra_VectorStdOps.hpp"
14#include "Thyra_MultiVectorStdOps.hpp"
23 Teuchos::RCP<Teuchos::ParameterList> inputPL,
29 adjoint_residual_model_ = adjoint_residual_model;
30 adjoint_solve_model_ = adjoint_solve_model;
31 state_integrator_ = createIntegratorBasic<Scalar>(inputPL, model_);
32 sens_model_ = createSensitivityModel(model_, adjoint_residual_model_,
33 adjoint_solve_model_, inputPL);
34 sens_integrator_ = createIntegratorBasic<Scalar>(inputPL, sens_model_);
36 do_forward_integration_ =
true;
37 do_adjoint_integration_ =
true;
46 std::string stepperType)
49 adjoint_residual_model_ = adjoint_residual_model;
50 adjoint_solve_model_ = adjoint_solve_model;
51 state_integrator_ = createIntegratorBasic<Scalar>(model_, stepperType);
52 sens_model_ = createSensitivityModel(model_, adjoint_residual_model_,
53 adjoint_solve_model_, Teuchos::null);
54 sens_integrator_ = createIntegratorBasic<Scalar>(sens_model_, stepperType);
56 do_forward_integration_ =
true;
57 do_adjoint_integration_ =
true;
63 Teuchos::RCP<Teuchos::ParameterList> inputPL,
67 inputPL, model, adjoint_model, adjoint_model)
76 std::string stepperType) :
78 model, adjoint_model, adjoint_model, stepperType)
85 Teuchos::RCP<Teuchos::ParameterList> inputPL,
88 inputPL, model,
Thyra::implicitAdjointModelEvaluator(model))
96 std::string stepperType) :
98 model,
Thyra::implicitAdjointModelEvaluator(model), stepperType)
102template<
class Scalar>
106 state_integrator_ = createIntegratorBasic<Scalar>();
107 sens_integrator_ = createIntegratorBasic<Scalar>();
111template<
class Scalar>
116 const Scalar tfinal =
117 state_integrator_->getTimeStepControl()->getFinalTime();
118 return advanceTime(tfinal);
121template<
class Scalar>
126 TEMPUS_FUNC_TIME_MONITOR_DIFF(
"Tempus::IntegratorPseudoTransientAdjointSensitivity::advanceTime()", TEMPUS_PTAS_AT);
130 typedef Thyra::ModelEvaluatorBase MEB;
132 bool state_status =
true;
133 if (do_forward_integration_) {
134 TEMPUS_FUNC_TIME_MONITOR_DIFF(
"Tempus::IntegratorPseudoTransientAdjointSensitivity::advanceTime::state", TEMPUS_PTAS_AT_FWD);
138 state_status = state_integrator_->advanceTime(timeFinal);
141 bool sens_status =
true;
142 if (do_adjoint_integration_) {
143 TEMPUS_FUNC_TIME_MONITOR_DIFF(
"Tempus::IntegratorPseudoTransientAdjointSensitivity::advanceTime::adjoint", TEMPUS_PTAS_AT_ADJ);
149 sens_model_->setFinalTime(state_integrator_->getTime());
152 sens_model_->setForwardSolutionHistory(
153 state_integrator_->getSolutionHistory());
157 sens_status = sens_integrator_->advanceTime(timeFinal);
161 MEB::InArgs<Scalar> inargs = sens_model_->getNominalValues();
162 MEB::OutArgs<Scalar> outargs = sens_model_->createOutArgs();
163 inargs.set_t(sens_integrator_->getTime());
164 inargs.set_x(sens_integrator_->getX());
165 if (inargs.supports(MEB::IN_ARG_x_dot))
166 inargs.set_x_dot(sens_integrator_->getXDot());
167 if (inargs.supports(MEB::IN_ARG_x_dot_dot))
168 inargs.set_x_dot_dot(sens_integrator_->getXDotDot());
169 RCP<VectorBase<Scalar> > G = dgdp_;
170 if (G == Teuchos::null) {
171 G = Thyra::createMember(sens_model_->get_g_space(0));
172 dgdp_ = Teuchos::rcp_dynamic_cast<DMVPV>(G);
174 if (g_ == Teuchos::null)
175 g_ = Thyra::createMember(sens_model_->get_g_space(1));
177 outargs.set_g(1, g_);
178 sens_model_->evalModel(inargs, outargs);
180 buildSolutionHistory();
183 return state_status && sens_status;
186template<
class Scalar>
191 return solutionHistory_->getCurrentTime();
194template<
class Scalar>
199 return solutionHistory_->getCurrentIndex();
202template<
class Scalar>
207 Status state_status = state_integrator_->getStatus();
208 Status sens_status = sens_integrator_->getStatus();
216template <
class Scalar>
219 state_integrator_->setStatus(st);
220 sens_integrator_->setStatus(st);
223template<
class Scalar>
224Teuchos::RCP<Stepper<Scalar> >
228 return state_integrator_->getStepper();
231template<
class Scalar>
232Teuchos::RCP<Stepper<Scalar> >
236 return state_integrator_->getStepper();
239template<
class Scalar>
240Teuchos::RCP<Stepper<Scalar> >
244 return sens_integrator_->getStepper();
247template<
class Scalar>
248Teuchos::RCP<const SolutionHistory<Scalar> >
252 return solutionHistory_;
255template<
class Scalar>
256Teuchos::RCP<const SolutionHistory<Scalar> >
260 return state_integrator_->getSolutionHistory();
263template<
class Scalar>
264Teuchos::RCP<const SolutionHistory<Scalar> >
268 return sens_integrator_->getSolutionHistory();
271template<
class Scalar>
272Teuchos::RCP<SolutionHistory<Scalar> >
276 return solutionHistory_;
279template<
class Scalar>
280Teuchos::RCP<const TimeStepControl<Scalar> >
284 return state_integrator_->getTimeStepControl();
287template<
class Scalar>
288Teuchos::RCP<TimeStepControl<Scalar> >
292 return state_integrator_->getNonConstTimeStepControl();
295template<
class Scalar>
296Teuchos::RCP<TimeStepControl<Scalar> >
300 return state_integrator_->getNonConstTimeStepControl();
303template<
class Scalar>
304Teuchos::RCP<TimeStepControl<Scalar> >
308 return sens_integrator_->getNonConstTimeStepControl();
311template<
class Scalar>
312Teuchos::RCP<IntegratorObserver<Scalar> >
316 return state_integrator_->getObserver();
319template<
class Scalar>
324 state_integrator_->setObserver(obs);
325 sens_integrator_->setObserver(obs);
328template<
class Scalar>
339 using Teuchos::rcp_dynamic_cast;
340 using Thyra::VectorSpaceBase;
342 using Thyra::createMember;
347 RCP< const VectorSpaceBase<Scalar> > space = sens_model_->get_x_space();
348 RCP<DMVPV> Y = rcp_dynamic_cast<DMVPV>(createMember(space));
349 RCP<DMVPV> Ydot = rcp_dynamic_cast<DMVPV>(createMember(space));
350 RCP<DMVPV> Ydotdot = rcp_dynamic_cast<DMVPV>(createMember(space));
351 const Scalar zero = Teuchos::ScalarTraits<Scalar>::zero();
354 if (y0 == Teuchos::null)
355 assign(Y->getNonconstMultiVector().ptr(), zero);
357 assign(Y->getNonconstMultiVector().ptr(), *y0);
360 if (ydot0 == Teuchos::null)
361 assign(Ydot->getNonconstMultiVector().ptr(), zero);
363 assign(Ydot->getNonconstMultiVector().ptr(), *ydot0);
366 if (ydotdot0 == Teuchos::null)
367 assign(Ydotdot->getNonconstMultiVector().ptr(), zero);
369 assign(Ydotdot->getNonconstMultiVector().ptr(), *ydotdot0);
371 state_integrator_->initializeSolutionHistory(t0, x0, xdot0, xdotdot0);
372 sens_integrator_->initializeSolutionHistory(t0, Y, Ydot, Ydotdot);
375template<
class Scalar>
376Teuchos::RCP<const Thyra::VectorBase<Scalar> >
380 return state_integrator_->getX();
383template<
class Scalar>
384Teuchos::RCP<const Thyra::VectorBase<Scalar> >
388 return state_integrator_->getXDot();
391template<
class Scalar>
392Teuchos::RCP<const Thyra::VectorBase<Scalar> >
396 return state_integrator_->getXDotDot();
399template<
class Scalar>
400Teuchos::RCP<const Thyra::MultiVectorBase<Scalar> >
405 using Teuchos::rcp_dynamic_cast;
406 RCP<const DMVPV> mvpv =
407 rcp_dynamic_cast<const DMVPV>(sens_integrator_->getX());
408 return mvpv->getMultiVector();
411template<
class Scalar>
412Teuchos::RCP<const Thyra::MultiVectorBase<Scalar> >
417 using Teuchos::rcp_dynamic_cast;
418 RCP<const DMVPV> mvpv =
419 rcp_dynamic_cast<const DMVPV>(sens_integrator_->getXDot());
420 return mvpv->getMultiVector();
423template<
class Scalar>
424Teuchos::RCP<const Thyra::MultiVectorBase<Scalar> >
429 using Teuchos::rcp_dynamic_cast;
430 RCP<const DMVPV> mvpv =
431 rcp_dynamic_cast<const DMVPV>(sens_integrator_->getXDotDot());
432 return mvpv->getMultiVector();
435template<
class Scalar>
436Teuchos::RCP<const Thyra::VectorBase<Scalar> >
443template<
class Scalar>
444Teuchos::RCP<const Thyra::MultiVectorBase<Scalar> >
448 return dgdp_->getMultiVector();
451template<
class Scalar>
456 std::string name =
"Tempus::IntegratorPseudoTransientAdjointSensitivity";
460template<
class Scalar>
464 Teuchos::FancyOStream &out,
465 const Teuchos::EVerbosityLevel verbLevel)
const
467 auto l_out = Teuchos::fancyOStream( out.getOStream() );
468 Teuchos::OSTab ostab(*l_out, 2, this->description());
469 l_out->setOutputToRootOnly(0);
471 *l_out << description() <<
"::describe" << std::endl;
472 state_integrator_->describe(*l_out, verbLevel);
473 sens_integrator_->describe(*l_out, verbLevel);
476template<
class Scalar>
484 auto model = Teuchos::rcp_const_cast<Thyra::ModelEvaluator<Scalar>> (
485 state_integrator_->getStepper()->getModel());
486 auto tmp_state_integrator = createIntegratorBasic<Scalar>(inputPL, model);
487 state_integrator_->copy(tmp_state_integrator);
489 model = Teuchos::rcp_const_cast<Thyra::ModelEvaluator<Scalar>> (
490 sens_integrator_->getStepper()->getModel());
491 auto tmp_sens_integrator = createIntegratorBasic<Scalar>(inputPL, model);
492 sens_integrator_->copy(tmp_sens_integrator);
495template<
class Scalar>
496Teuchos::RCP<Teuchos::ParameterList>
503 auto tmp_state_integrator = createIntegratorBasic<Scalar>();
504 auto model = state_integrator_->getStepper()->getModel();
505 tmp_state_integrator->setModel(model);
506 state_integrator_->copy(tmp_state_integrator);
508 auto tmp_sens_integrator = createIntegratorBasic<Scalar>();
509 model = sens_integrator_->getStepper()->getModel();
510 tmp_sens_integrator->setModel(model);
511 sens_integrator_->copy(tmp_sens_integrator);
513 auto pl = Teuchos::rcp_const_cast<Teuchos::ParameterList> (
514 sens_integrator_->getValidParameters());
518template<
class Scalar>
519Teuchos::RCP<const Teuchos::ParameterList>
523 Teuchos::RCP<Teuchos::ParameterList> pl =
524 Teuchos::rcp(
new Teuchos::ParameterList);
525 Teuchos::RCP<const Teuchos::ParameterList> integrator_pl =
526 state_integrator_->getValidParameters();
527 Teuchos::RCP<const Teuchos::ParameterList> sensitivity_pl =
529 pl->setParameters(*integrator_pl);
530 pl->sublist(
"Sensitivities").setParameters(*sensitivity_pl);
535template<
class Scalar>
536Teuchos::RCP<Teuchos::ParameterList>
540 auto pl = Teuchos::rcp_const_cast<Teuchos::ParameterList> (
541 state_integrator_->getValidParameters());
545template<
class Scalar>
553template <
class Scalar>
554Teuchos::RCP<AdjointSensitivityModelEvaluator<Scalar> >
560 const Teuchos::RCP<Teuchos::ParameterList>& inputPL)
564 Teuchos::RCP<Teuchos::ParameterList> pl = Teuchos::parameterList();
565 if (inputPL != Teuchos::null) {
566 *pl = inputPL->sublist(
"Sensitivities");
568 const Scalar tinit = state_integrator_->getTimeStepControl()->getInitTime();
569 const Scalar tfinal = state_integrator_->getTimeStepControl()->getFinalTime();
571 model, adjoint_residual_model, adjoint_solve_model_,
572 tinit, tfinal,
true, pl));
575template<
class Scalar>
582 using Teuchos::rcp_dynamic_cast;
583 using Teuchos::ParameterList;
585 using Thyra::VectorSpaceBase;
587 using Thyra::defaultProductVector;
588 typedef Thyra::DefaultProductVector<Scalar> DPV;
589 typedef Thyra::DefaultProductVectorSpace<Scalar> DPVS;
593 auto shPL = Teuchos::rcp_const_cast<Teuchos::ParameterList> (
594 state_integrator_->getSolutionHistory()->getValidParameters());
595 solutionHistory_ = createSolutionHistoryPL<Scalar>(shPL);
597 RCP<const VectorSpaceBase<Scalar> > x_space =
598 model_->get_x_space();
599 RCP<const VectorSpaceBase<Scalar> > adjoint_space =
600 sens_model_->get_x_space();
601 Teuchos::Array< RCP<const VectorSpaceBase<Scalar> > > spaces(2);
603 spaces[1] = adjoint_space;
604 RCP<const DPVS > prod_space = Thyra::productVectorSpace(spaces());
605 const Scalar zero = Teuchos::ScalarTraits<Scalar>::zero();
607 RCP<const SolutionHistory<Scalar> > state_solution_history =
608 state_integrator_->getSolutionHistory();
609 int num_states = state_solution_history->getNumStates();
610 for (
int i=0; i<num_states; ++i) {
611 RCP<const SolutionState<Scalar> > state = (*state_solution_history)[i];
614 RCP<DPV> x = defaultProductVector(prod_space);
615 assign(x->getNonconstVectorBlock(0).ptr(), *(state->getX()));
616 assign(x->getNonconstVectorBlock(1).ptr(), zero);
617 RCP<VectorBase<Scalar> > x_b = x;
620 RCP<VectorBase<Scalar> > x_dot_b;
621 if (state->getXDot() != Teuchos::null) {
622 RCP<DPV> x_dot = defaultProductVector(prod_space);
623 assign(x_dot->getNonconstVectorBlock(0).ptr(), *(state->getXDot()));
624 assign(x_dot->getNonconstVectorBlock(1).ptr(), zero);
629 RCP<VectorBase<Scalar> > x_dot_dot_b;
630 if (state->getXDotDot() != Teuchos::null) {
631 RCP<DPV> x_dot_dot = defaultProductVector(prod_space);
632 assign(x_dot_dot->getNonconstVectorBlock(0).ptr(),*(state->getXDotDot()));
633 assign(x_dot_dot->getNonconstVectorBlock(1).ptr(), zero);
634 x_dot_dot_b = x_dot_dot;
637 RCP<SolutionState<Scalar> > prod_state = state->clone();
638 prod_state->setX(x_b);
639 prod_state->setXDot(x_dot_b);
640 prod_state->setXDotDot(x_dot_dot_b);
641 solutionHistory_->addState(prod_state);
644 RCP<const VectorBase<Scalar> > frozen_x =
645 state_solution_history->getCurrentState()->getX();
646 RCP<const VectorBase<Scalar> > frozen_x_dot =
647 state_solution_history->getCurrentState()->getXDot();
648 RCP<const VectorBase<Scalar> > frozen_x_dot_dot =
649 state_solution_history->getCurrentState()->getXDotDot();
650 RCP<const SolutionHistory<Scalar> > sens_solution_history =
651 sens_integrator_->getSolutionHistory();
652 num_states = sens_solution_history->getNumStates();
653 for (
int i=0; i<num_states; ++i) {
654 RCP<const SolutionState<Scalar> > state = (*sens_solution_history)[i];
657 RCP<DPV> x = defaultProductVector(prod_space);
658 assign(x->getNonconstVectorBlock(0).ptr(), *frozen_x);
659 assign(x->getNonconstVectorBlock(1).ptr(), *(state->getX()));
660 RCP<VectorBase<Scalar> > x_b = x;
663 RCP<VectorBase<Scalar> > x_dot_b;
664 if (state->getXDot() != Teuchos::null) {
665 RCP<DPV> x_dot = defaultProductVector(prod_space);
666 assign(x_dot->getNonconstVectorBlock(0).ptr(), *frozen_x_dot);
667 assign(x_dot->getNonconstVectorBlock(1).ptr(), *(state->getXDot()));
672 RCP<VectorBase<Scalar> > x_dot_dot_b;
673 if (state->getXDotDot() != Teuchos::null) {
674 RCP<DPV> x_dot_dot = defaultProductVector(prod_space);
675 assign(x_dot_dot->getNonconstVectorBlock(0).ptr(), *frozen_x_dot_dot);
676 assign(x_dot_dot->getNonconstVectorBlock(1).ptr(),*(state->getXDotDot()));
677 x_dot_dot_b = x_dot_dot;
680 RCP<SolutionState<Scalar> > prod_state = state->clone();
681 prod_state->setX(x_b);
682 prod_state->setXDot(x_dot_b);
683 prod_state->setXDotDot(x_dot_dot_b);
684 solutionHistory_->addState(prod_state,
false);
689template<
class Scalar>
690Teuchos::RCP<IntegratorPseudoTransientAdjointSensitivity<Scalar> >
692 Teuchos::RCP<Teuchos::ParameterList> pList,
695 Teuchos::RCP<IntegratorPseudoTransientAdjointSensitivity<Scalar> > integrator =
701template<
class Scalar>
702Teuchos::RCP<IntegratorPseudoTransientAdjointSensitivity<Scalar> >
705 std::string stepperType)
707 Teuchos::RCP<IntegratorPseudoTransientAdjointSensitivity<Scalar> > integrator =
713template<
class Scalar>
714Teuchos::RCP<IntegratorPseudoTransientAdjointSensitivity<Scalar> >
716 Teuchos::RCP<Teuchos::ParameterList> pList,
720 Teuchos::RCP<IntegratorPseudoTransientAdjointSensitivity<Scalar> > integrator =
726template<
class Scalar>
727Teuchos::RCP<IntegratorPseudoTransientAdjointSensitivity<Scalar> >
731 std::string stepperType)
733 Teuchos::RCP<IntegratorPseudoTransientAdjointSensitivity<Scalar> > integrator =
739template<
class Scalar>
740Teuchos::RCP<IntegratorPseudoTransientAdjointSensitivity<Scalar> >
742 Teuchos::RCP<Teuchos::ParameterList> pList,
747 Teuchos::RCP<IntegratorPseudoTransientAdjointSensitivity<Scalar> > integrator =
753template<
class Scalar>
754Teuchos::RCP<IntegratorPseudoTransientAdjointSensitivity<Scalar> >
759 std::string stepperType)
761 Teuchos::RCP<IntegratorPseudoTransientAdjointSensitivity<Scalar> > integrator =
767template<
class Scalar>
768Teuchos::RCP<IntegratorPseudoTransientAdjointSensitivity<Scalar> >
771 Teuchos::RCP<IntegratorPseudoTransientAdjointSensitivity<Scalar> > integrator =
ModelEvaluator for forming adjoint sensitivity equations.
static Teuchos::RCP< const Teuchos::ParameterList > getValidParameters()
IntegratorObserver class for time integrators.
Time integrator suitable for pseudotransient adjoint sensitivity analysis.
virtual Teuchos::RCP< SolutionHistory< Scalar > > getNonConstSolutionHistory() override
Get the SolutionHistory.
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getY() const
Get the current adjoint solution, y.
Teuchos::RCP< Stepper< Scalar > > getSensStepper() const
virtual Scalar getTime() const override
Get current time.
SensitivityStepMode getStepMode() const
What mode the current time integration step is in.
IntegratorPseudoTransientAdjointSensitivity()
Destructor.
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getYDot() const
Get the current time derivative of the adjoint solution, ydot.
Teuchos::RCP< Teuchos::ParameterList > getNonconstParameterList() override
Teuchos::RCP< AdjointSensitivityModelEvaluator< Scalar > > createSensitivityModel(const Teuchos::RCP< Thyra::ModelEvaluator< Scalar > > &model, const Teuchos::RCP< Thyra::ModelEvaluator< Scalar > > &adjoint_residual_model, const Teuchos::RCP< Thyra::ModelEvaluator< Scalar > > &adjoint_solve_model, const Teuchos::RCP< Teuchos::ParameterList > &inputPL)
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getXDotDot() const
Get the current second time derivative of the solution, xdotdot.
Teuchos::RCP< const SolutionHistory< Scalar > > getSensSolutionHistory() const
virtual Teuchos::RCP< Stepper< Scalar > > getStepper() const override
Get the Stepper.
virtual Teuchos::RCP< const TimeStepControl< Scalar > > getTimeStepControl() const override
Get the TimeStepControl.
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getG() const
Return response function g.
virtual void setObserver(Teuchos::RCP< IntegratorObserver< Scalar > > obs=Teuchos::null)
Set the Observer.
virtual Teuchos::RCP< const SolutionHistory< Scalar > > getSolutionHistory() const override
Get the SolutionHistory.
Teuchos::RCP< Teuchos::ParameterList > unsetParameterList() override
virtual Teuchos::RCP< IntegratorObserver< Scalar > > getObserver()
Get the Observer.
Teuchos::RCP< Stepper< Scalar > > getStateStepper() const
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 > > y0=Teuchos::null, Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > ydot0=Teuchos::null, Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > ydotdot0=Teuchos::null)
Set the initial state from Thyra::VectorBase(s)
Teuchos::RCP< TimeStepControl< Scalar > > getSensNonConstTimeStepControl()
Teuchos::RCP< TimeStepControl< Scalar > > getStateNonConstTimeStepControl()
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getXDot() const
Get the current time derivative of the solution, xdot.
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > &pl) override
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getDgDp() const
Return adjoint sensitivity stored in gradient format.
virtual Teuchos::RCP< TimeStepControl< Scalar > > getNonConstTimeStepControl() override
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getX() const
Get the current solution, x.
void buildSolutionHistory()
virtual Status getStatus() const override
Get Status.
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getYDotDot() const
Get the current second time derivative of the adjoint solution, ydotdot.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const override
virtual bool advanceTime()
Advance the solution to timeMax, and return true if successful.
virtual void setStatus(const Status st) override
Set Status.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
Teuchos::RCP< const SolutionHistory< Scalar > > getStateSolutionHistory() const
virtual int getIndex() const override
Get current index.
std::string description() const override
Status
Status for the Integrator, the Stepper and the SolutionState.
Teuchos::RCP< Tempus::IntegratorPseudoTransientAdjointSensitivity< Scalar > > integratorPseudoTransientAdjointSensitivity()
Nonmember constructor.