9 #ifndef Tempus_IntegratorForwardSensitivity_impl_hpp
10 #define Tempus_IntegratorForwardSensitivity_impl_hpp
13 #include "Thyra_DefaultMultiVectorProductVector.hpp"
14 #include "Thyra_VectorStdOps.hpp"
15 #include "Thyra_MultiVectorStdOps.hpp"
17 #include "Tempus_CombinedForwardSensitivityModelEvaluator.hpp"
22 template <
class Scalar>
28 bool use_combined_method)
30 , integrator_(integrator)
31 , sens_model_(sens_model)
32 , sens_stepper_(sens_stepper)
33 , use_combined_method_(use_combined_method)
38 template<
class Scalar>
42 integrator_ = createIntegratorBasic<Scalar>();
43 integrator_->initialize();
46 template<
class Scalar>
51 if (use_combined_method_)
52 integrator_->setModel(sens_model_);
54 integrator_->setStepper(sens_stepper_);
57 template<
class Scalar>
68 using Teuchos::rcp_dynamic_cast;
71 using Thyra::createMember;
77 RCP< const VectorSpaceBase<Scalar> > space;
78 if (use_combined_method_)
79 space = sens_model_->get_x_space();
81 space = sens_stepper_->get_x_space();
82 RCP<DMVPV> X = rcp_dynamic_cast<DMVPV>(createMember(space));
83 RCP<DMVPV> Xdot = rcp_dynamic_cast<DMVPV>(createMember(space));
84 RCP<DMVPV> Xdotdot = rcp_dynamic_cast<DMVPV>(createMember(space));
91 assign(X->getNonconstMultiVector()->col(0).ptr(), *x0);
92 if (DxDp0 == Teuchos::null)
93 assign(X->getNonconstMultiVector()->subView(rng).ptr(), zero);
95 assign(X->getNonconstMultiVector()->subView(rng).ptr(), *DxDp0);
98 if (xdot0 == Teuchos::null)
99 assign(Xdot->getNonconstMultiVector()->col(0).ptr(), zero);
101 assign(Xdot->getNonconstMultiVector()->col(0).ptr(), *xdot0);
102 if (DxdotDp0 == Teuchos::null)
103 assign(Xdot->getNonconstMultiVector()->subView(rng).ptr(), zero);
105 assign(Xdot->getNonconstMultiVector()->subView(rng).ptr(), *DxdotDp0);
108 if (xdotdot0 == Teuchos::null)
109 assign(Xdotdot->getNonconstMultiVector()->col(0).ptr(), zero);
111 assign(Xdotdot->getNonconstMultiVector()->col(0).ptr(), *xdotdot0);
112 if (DxdotDp0 == Teuchos::null)
113 assign(Xdotdot->getNonconstMultiVector()->subView(rng).ptr(), zero);
115 assign(Xdotdot->getNonconstMultiVector()->subView(rng).ptr(), *DxdotdotDp0);
117 integrator_->initializeSolutionHistory(t0, X, Xdot, Xdotdot);
120 template<
class Scalar>
126 using Teuchos::rcp_dynamic_cast;
129 RCP<const DMVPV> X = rcp_dynamic_cast<
const DMVPV>(integrator_->getX());
133 template<
class Scalar>
139 using Teuchos::rcp_dynamic_cast;
142 RCP<const DMVPV> X = rcp_dynamic_cast<
const DMVPV>(integrator_->getX());
145 return X->getMultiVector()->subView(rng);
148 template<
class Scalar>
154 using Teuchos::rcp_dynamic_cast;
157 RCP<const DMVPV> Xdot = rcp_dynamic_cast<
const DMVPV>(integrator_->getXDot());
161 template<
class Scalar>
167 using Teuchos::rcp_dynamic_cast;
170 RCP<const DMVPV> Xdot = rcp_dynamic_cast<
const DMVPV>(integrator_->getXDot());
173 return Xdot->getMultiVector()->subView(rng);
176 template<
class Scalar>
182 using Teuchos::rcp_dynamic_cast;
185 RCP<const DMVPV> Xdotdot =
186 rcp_dynamic_cast<
const DMVPV>(integrator_->getXDotDot());
190 template<
class Scalar>
196 using Teuchos::rcp_dynamic_cast;
199 RCP<const DMVPV> Xdotdot =
200 rcp_dynamic_cast<
const DMVPV>(integrator_->getXDotDot());
201 const int num_param = Xdotdot->
getMultiVector()->domain()->dim()-1;
203 return Xdotdot->getMultiVector()->subView(rng);
206 template<
class Scalar>
216 if (use_combined_method_)
217 smodel = sens_model_;
219 smodel = sens_stepper_->getModel();
220 MEB::InArgs<Scalar> inargs = smodel->getNominalValues();
221 MEB::OutArgs<Scalar> outargs = smodel->createOutArgs();
222 inargs.set_t(integrator_->getTime());
223 inargs.set_x(integrator_->getX());
224 if (inargs.supports(MEB::IN_ARG_x_dot))
225 inargs.set_x_dot(integrator_->getXDot());
226 if (inargs.supports(MEB::IN_ARG_x_dot_dot))
227 inargs.set_x_dot_dot(integrator_->getXDotDot());
230 Thyra::createMember(smodel->get_g_space(1));
233 smodel->evalModel(inargs, outargs);
237 template<
class Scalar>
248 if (use_combined_method_)
249 smodel = sens_model_;
251 smodel = sens_stepper_->getModel();
252 MEB::InArgs<Scalar> inargs = smodel->getNominalValues();
253 MEB::OutArgs<Scalar> outargs = smodel->createOutArgs();
254 inargs.set_t(integrator_->getTime());
255 inargs.set_x(integrator_->getX());
256 if (inargs.supports(MEB::IN_ARG_x_dot))
257 inargs.set_x_dot(integrator_->getXDot());
258 if (inargs.supports(MEB::IN_ARG_x_dot_dot))
259 inargs.set_x_dot_dot(integrator_->getXDotDot());
262 Thyra::createMember(smodel->get_g_space(0));
266 smodel->evalModel(inargs, outargs);
267 return dgdp->getMultiVector();
270 template<
class Scalar>
275 std::string name =
"Tempus::IntegratorForwardSensitivity";
279 template<
class Scalar>
286 auto l_out = Teuchos::fancyOStream( out.
getOStream() );
288 l_out->setOutputToRootOnly(0);
290 *l_out << description() <<
"::describe" << std::endl;
291 integrator_->describe(*l_out, verbLevel);
294 template <
class Scalar>
299 if (use_combined_method_)
301 return sens_stepper_->getStepMode();
305 template<
class Scalar>
319 auto fwd_integrator = createIntegratorBasic<Scalar>(pList, model);
327 pl->setParameters(*integrator_pl);
330 spl.
set(
"Sensitivity Method",
"Combined");
331 spl.
set(
"Reuse State Linear Solver",
false);
335 auto sens_pl = Teuchos::sublist(pList,
"Sensitivities",
false);
336 std::string integratorName = pList->get<std::string>(
"Integrator Name",
"Default Integrator");
337 std::string stepperName = pList->sublist(integratorName).get<std::string>(
"Stepper Name");
338 auto stepper_pl = Teuchos::sublist(pList, stepperName,
true);
339 std::string sensitivity_method = sens_pl->get<std::string>(
"Sensitivity Method");
340 bool use_combined_method = sensitivity_method ==
"Combined";
345 sens_pl->remove(
"Sensitivity Method");
347 if (use_combined_method)
349 sens_pl->remove(
"Reuse State Linear Solver");
351 model, sens_residual_model, sens_solve_model, sens_pl);
352 fwd_integrator = createIntegratorBasic<Scalar>(pList, sens_model);
358 fwd_integrator = createIntegratorBasic<Scalar>(pList, fsa_staggered_me);
359 fwd_integrator->setStepper(sens_stepper);
370 template<
class Scalar>
381 #endif // Tempus_IntegratorForwardSensitivity_impl_hpp
Teuchos::RCP< IntegratorBasic< Scalar > > integrator_
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getXDotDot() const
Get current the second time derivative of the solution, xdotdot, only. This is the first block only a...
std::string description() const override
virtual void initializeSolutionHistory(Teuchos::RCP< SolutionState< Scalar > > state=Teuchos::null)
Set the initial state which has the initial conditions.
SensitivityStepMode getStepMode() const
What mode the current time integration step is in.
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getDXDotDp() const
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getDXDotDotDp() const
static Teuchos::RCP< const Teuchos::ParameterList > getValidParameters()
A stepper implementing staggered forward sensitivity analysis.
virtual void setStepper(Teuchos::RCP< Thyra::ModelEvaluator< Scalar > > model)
Set the Stepper.
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getDgDp() const
Return forward sensitivity stored in Jacobian format.
IntegratorForwardSensitivity()
Destructor.
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getDxDp() const
Get the forward sensitivities .
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
RCP< const MultiVectorBase< Scalar > > getMultiVector() const
Teuchos::RCP< IntegratorForwardSensitivity< Scalar > > createIntegratorForwardSensitivity(Teuchos::RCP< Teuchos::ParameterList > pList, const Teuchos::RCP< Thyra::ModelEvaluator< Scalar >> &model, const Teuchos::RCP< Thyra::ModelEvaluator< Scalar >> &sens_residual_model, const Teuchos::RCP< Thyra::ModelEvaluator< Scalar >> &sens_solve_model)
Nonmember constructor.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const override
ParameterList & setParameters(const ParameterList &source)
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getG() const
Return response function g.
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getX() const
Get the current solution, x, only. If looking for the solution vector and the sensitivities, use SolutionState->getX() which will return a Block MultiVector with the first block containing the current solution, x, and the remaining blocks are the forward sensitivities .
RCP< std::basic_ostream< char_type, traits_type > > getOStream()
ParameterList & setParametersNotAlreadySet(const ParameterList &source)
RCP< MultiVectorBase< Scalar > > getNonconstMultiVector()
ParameterList & sublist(const std::string &name, bool mustAlreadyExist=false, const std::string &docString="")
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getXDot() const
Get current the time derivative of the solution, xdot, only. This is the first block only and not the...
A ModelEvaluator decorator for sensitivity analysis.
Time integrator implementing forward sensitivity analysis.
Teuchos::RCP< SensitivityModelEvaluatorBase< Scalar > > wrapCombinedFSAModelEvaluator(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)