9 #ifndef Tempus_StepperStaggeredForwardSensitivity_impl_hpp
10 #define Tempus_StepperStaggeredForwardSensitivity_impl_hpp
12 #include "Tempus_config.hpp"
16 #include "Teuchos_VerboseObjectParameterListHelpers.hpp"
18 #include "Thyra_DefaultMultiVectorProductVectorSpace.hpp"
19 #include "Thyra_DefaultMultiVectorProductVector.hpp"
24 template<
class Scalar>
class StepperFactory;
27 template<
class Scalar>
31 this->setParams(Teuchos::null, Teuchos::null);
36 template<
class Scalar>
39 const Teuchos::RCP<
const Thyra::ModelEvaluator<Scalar> >& appModel,
40 const Teuchos::RCP<Teuchos::ParameterList>& pList,
41 const Teuchos::RCP<Teuchos::ParameterList>& sens_pList)
44 this->setParams(pList, sens_pList);
45 this->setModel(appModel);
50 template<
class Scalar>
53 const Teuchos::RCP<
const Thyra::ModelEvaluator<Scalar> >& appModel)
57 using Teuchos::ParameterList;
60 Teuchos::RCP<Teuchos::ParameterList> spl = Teuchos::parameterList();
62 spl->remove(
"Reuse State Linear Solver");
63 spl->remove(
"Force W Update");
74 if (stateStepper_ == Teuchos::null)
75 stateStepper_ = sf->createStepper(stepperPL_, appModel);
77 stateStepper_->setModel(appModel);
78 if (sensitivityStepper_ == Teuchos::null)
79 sensitivityStepper_ = sf->createStepper(stepperPL_, fsa_model_);
81 sensitivityStepper_->setModel(fsa_model_);
85 template<
class Scalar>
88 const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& appModel)
90 this->setModel(appModel);
94 template<
class Scalar>
98 stateStepper_->setSolver(solverName);
99 sensitivityStepper_->setSolver(solverName);
102 template<
class Scalar>
103 Teuchos::RCP<const Thyra::ModelEvaluator<Scalar> >
107 return combined_fsa_model_;
111 template<
class Scalar>
114 Teuchos::RCP<Teuchos::ParameterList> solverPL)
116 stateStepper_->setSolver(solverPL);
117 sensitivityStepper_->setSolver(solverPL);
121 template<
class Scalar>
124 Teuchos::RCP<Thyra::NonlinearSolverBase<Scalar> > solver)
126 stateStepper_->setSolver(solver);
127 sensitivityStepper_->setSolver(solver);
131 template<
class Scalar>
139 template<
class Scalar>
146 using Teuchos::rcp_dynamic_cast;
147 using Thyra::VectorBase;
148 using Thyra::MultiVectorBase;
150 using Thyra::createMember;
151 using Thyra::multiVectorProductVector;
152 using Thyra::multiVectorProductVectorSpace;
153 typedef Thyra::DefaultMultiVectorProductVector<Scalar> DMVPV;
154 typedef Thyra::DefaultMultiVectorProductVectorSpace<Scalar> DMVPVS;
160 if (stateSolutionHistory_ == Teuchos::null) {
161 RCP<Teuchos::ParameterList> shPL =
165 RCP<SolutionState<Scalar> > state =
solutionHistory->getCurrentState();
166 RCP<DMVPV> X, XDot, XDotDot;
167 X = rcp_dynamic_cast<DMVPV>(state->getX(),
true);
168 XDot = rcp_dynamic_cast<DMVPV>(state->getXDot(),
true);
169 if (state->getXDotDot() != Teuchos::null)
170 XDotDot = rcp_dynamic_cast<DMVPV>(state->getXDotDot(),
true);
174 RCP<VectorBase<Scalar> > x, xdot, xdotdot;
175 x = X->getNonconstMultiVector()->col(0);
176 xdot = XDot->getNonconstMultiVector()->col(0);
177 if (XDotDot != Teuchos::null)
178 xdotdot = XDotDot->getNonconstMultiVector()->col(0);
181 RCP<SolutionState<Scalar> > state_state =
184 state->getStepperState()->clone()));
186 stateSolutionHistory_->addState(state_state);
188 const int num_param = X->getMultiVector()->domain()->dim()-1;
189 TEUCHOS_ASSERT(num_param > 0);
190 const Teuchos::Range1D rng(1,num_param);
193 RCP<MultiVectorBase<Scalar> > dxdp, dxdotdp, dxdotdotdp;
194 dxdp = X->getNonconstMultiVector()->subView(rng);
195 dxdotdp = XDot->getNonconstMultiVector()->subView(rng);
196 if (XDotDot != Teuchos::null)
197 dxdotdotdp = XDotDot->getNonconstMultiVector()->subView(rng);
200 RCP<VectorBase<Scalar> > dxdp_vec, dxdotdp_vec, dxdotdotdp_vec;
201 RCP<const DMVPVS> prod_space =
202 multiVectorProductVectorSpace(X->getMultiVector()->range(), num_param);
203 dxdp_vec = multiVectorProductVector(prod_space, dxdp);
204 dxdotdp_vec = multiVectorProductVector(prod_space, dxdotdp);
205 if (XDotDot != Teuchos::null)
206 dxdotdotdp_vec = multiVectorProductVector(prod_space, dxdotdotdp);
209 RCP<SolutionState<Scalar> > sens_state =
211 dxdp_vec, dxdotdp_vec, dxdotdotdp_vec,
212 state->getStepperState()->clone()));
214 sensSolutionHistory_->addState(sens_state);
218 RCP<SolutionState<Scalar> > prod_state =
solutionHistory->getWorkingState();
219 RCP<DMVPV> X, XDot, XDotDot;
220 X = rcp_dynamic_cast<DMVPV>(prod_state->getX(),
true);
221 XDot = rcp_dynamic_cast<DMVPV>(prod_state->getXDot(),
true);
222 if (prod_state->getXDotDot() != Teuchos::null)
223 XDotDot = rcp_dynamic_cast<DMVPV>(prod_state->getXDotDot(),
true);
226 stateSolutionHistory_->initWorkingState();
227 RCP<SolutionState<Scalar> > state = stateSolutionHistory_->getWorkingState();
228 state->getMetaData()->copy(prod_state->getMetaData());
229 stateStepper_->takeStep(stateSolutionHistory_);
232 assign(X->getNonconstMultiVector()->col(0).ptr(), *(state->getX()));
233 assign(XDot->getNonconstMultiVector()->col(0).ptr(), *(state->getXDot()));
234 if (XDotDot != Teuchos::null)
235 assign(XDotDot->getNonconstMultiVector()->col(0).ptr(),
236 *(state->getXDotDot()));
237 prod_state->setOrder(state->getOrder());
244 stateSolutionHistory_->promoteWorkingState();
247 fsa_model_->setForwardSolutionHistory(stateSolutionHistory_);
248 if (reuse_solver_ && stateStepper_->getSolver() != Teuchos::null)
249 fsa_model_->setSolver(stateStepper_->getSolver(), force_W_update_);
252 sensSolutionHistory_->initWorkingState();
253 RCP<SolutionState<Scalar> > sens_state =
254 sensSolutionHistory_->getWorkingState();
255 sens_state->getMetaData()->copy(prod_state->getMetaData());
256 sensitivityStepper_->takeStep(sensSolutionHistory_);
259 RCP<const MultiVectorBase<Scalar> > dxdp =
260 rcp_dynamic_cast<
const DMVPV>(sens_state->getX(),
true)->getMultiVector();
261 const int num_param = dxdp->domain()->dim();
262 const Teuchos::Range1D rng(1,num_param);
263 assign(X->getNonconstMultiVector()->subView(rng).ptr(), *dxdp);
264 RCP<const MultiVectorBase<Scalar> > dxdotdp =
265 rcp_dynamic_cast<
const DMVPV>(sens_state->getXDot(),
true)->getMultiVector();
266 assign(XDot->getNonconstMultiVector()->subView(rng).ptr(), *dxdotdp);
267 if (sens_state->getXDotDot() != Teuchos::null) {
268 RCP<const MultiVectorBase<Scalar> > dxdotdotdp =
269 rcp_dynamic_cast<
const DMVPV>(sens_state->getXDotDot(),
true)->getMultiVector();
270 assign(XDotDot->getNonconstMultiVector()->subView(rng).ptr(), *dxdotdotdp);
272 prod_state->setOrder(std::min(state->getOrder(), sens_state->getOrder()));
279 sensSolutionHistory_->promoteWorkingState();
285 template<
class Scalar>
286 Teuchos::RCP<Tempus::StepperState<Scalar> >
292 Teuchos::RCP<Tempus::StepperState<Scalar> > stepperState =
298 template<
class Scalar>
302 std::string name =
"StepperStaggeredForwardSensitivity";
307 template<
class Scalar>
310 Teuchos::FancyOStream &out,
311 const Teuchos::EVerbosityLevel verbLevel)
const
313 out << description() <<
"::describe:" << std::endl;
314 stateStepper_->describe(out, verbLevel);
316 sensitivityStepper_->describe(out, verbLevel);
320 template <
class Scalar>
323 Teuchos::RCP<Teuchos::ParameterList>
const& pList)
325 if (pList == Teuchos::null) {
327 if (stepperPL_ == Teuchos::null) stepperPL_ = this->getDefaultParameters();
336 template<
class Scalar>
337 Teuchos::RCP<const Teuchos::ParameterList>
341 return stateStepper_->getValidParameters();
345 template<
class Scalar>
346 Teuchos::RCP<Teuchos::ParameterList>
350 return stateStepper_->getDefaultParameters();
354 template <
class Scalar>
355 Teuchos::RCP<Teuchos::ParameterList>
363 template <
class Scalar>
364 Teuchos::RCP<Teuchos::ParameterList>
368 Teuchos::RCP<Teuchos::ParameterList> temp_plist = stepperPL_;
369 stepperPL_ = Teuchos::null;
374 template <
class Scalar>
377 Teuchos::RCP<Teuchos::ParameterList>
const& pList,
378 Teuchos::RCP<Teuchos::ParameterList>
const& spList)
380 if (pList == Teuchos::null)
381 stepperPL_ = this->getDefaultParameters();
385 if (spList == Teuchos::null)
386 sensPL_ = Teuchos::parameterList();
390 reuse_solver_ = sensPL_->get(
"Reuse State Linear Solver",
false);
391 force_W_update_ = sensPL_->get(
"Force W Update",
true);
398 template <
class Scalar>
399 Teuchos::RCP<const Thyra::VectorSpaceBase<Scalar> >
404 using Teuchos::rcp_dynamic_cast;
405 typedef Thyra::DefaultMultiVectorProductVectorSpace<Scalar> DMVPVS;
407 RCP<const Thyra::VectorSpaceBase<Scalar> > x_space =
408 stateStepper_->getModel()->get_x_space();
409 RCP<const DMVPVS> dxdp_space =
410 rcp_dynamic_cast<
const DMVPVS>(sensitivityStepper_->getModel()->get_x_space(),
true);
411 const int num_param = dxdp_space->numBlocks();
412 RCP<const Thyra::VectorSpaceBase<Scalar> > prod_space =
413 multiVectorProductVectorSpace(x_space, num_param+1);
419 #endif // Tempus_StepperStaggeredForwardSensitivity_impl_hpp
virtual void initialize()
Initialize during construction and after changing input parameters.
virtual void setSolver(std::string solverName)
Set solver via ParameterList solver name.
Teuchos::RCP< Teuchos::ParameterList > getDefaultParameters() const
Teuchos::RCP< Teuchos::ParameterList > getNonconstParameterList()
virtual void setNonConstModel(const Teuchos::RCP< Thyra::ModelEvaluator< Scalar > > &appModel)
virtual Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > getModel()
virtual void takeStep(const Teuchos::RCP< SolutionHistory< Scalar > > &solutionHistory)
Take the specified timestep, dt, and return true if successful.
StepperStaggeredForwardSensitivity()
Default constructor.
Teuchos::RCP< SensitivityModelEvaluatorBase< Scalar > > wrapStaggeredFSAModelEvaluator(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, const Teuchos::RCP< const Teuchos::ParameterList > &pList=Teuchos::null)
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > get_x_space() const
virtual Teuchos::RCP< Tempus::StepperState< Scalar > > getDefaultStepperState()
Get a default (initial) StepperState.
virtual void setModel(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &appModel)
StepperState is a simple class to hold state information about the stepper.
Teuchos::RCP< Teuchos::ParameterList > unsetParameterList()
Teuchos::RCP< SolutionHistory< Scalar > > solutionHistory(Teuchos::RCP< Teuchos::ParameterList > pList=Teuchos::null)
Nonmember constructor.
SolutionHistory is basically a container of SolutionStates. SolutionHistory maintains a collection of...
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
Teuchos::RCP< SensitivityModelEvaluatorBase< Scalar > > wrapCombinedFSAModelEvaluator(const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &model, const Teuchos::RCP< const Teuchos::ParameterList > &pList=Teuchos::null)
void setParams(const Teuchos::RCP< Teuchos::ParameterList > &pl, const Teuchos::RCP< Teuchos::ParameterList > &spl)
Solution state for integrators and steppers. SolutionState contains the metadata for solutions and th...
virtual std::string description() const
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > &pl)
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const