Tempus  Version of the Day
Time Integration
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Tempus_IntegratorAdjointSensitivity_impl.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ****************************************************************************
3 // Tempus: Copyright (2017) Sandia Corporation
4 //
5 // Distributed under BSD 3-clause license (See accompanying file Copyright.txt)
6 // ****************************************************************************
7 // @HEADER
8 
9 #ifndef Tempus_IntegratorAdjointSensitivity_impl_hpp
10 #define Tempus_IntegratorAdjointSensitivity_impl_hpp
11 
12 #include "Teuchos_VerboseObjectParameterListHelpers.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"
19 #include "NOX_Thyra.H"
20 
21 namespace Tempus {
22 
23 template<class Scalar>
26  Teuchos::RCP<Teuchos::ParameterList> inputPL,
27  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& model)
28 {
29  model_ = model;
30  setParameterList(inputPL);
31  state_integrator_ = integratorBasic<Scalar>(inputPL, model_);
32 
33  TEUCHOS_TEST_FOR_EXCEPTION( getStepper()->getUseFSAL(), std::logic_error,
34  "Error - IntegratorAdjointSensitivity(): Cannot use FSAL with\n"
35  " IntegratorAdjointSensitivity, because the state and adjoint\n"
36  " integrators require ModelEvaluator evaluation in the\n"
37  " constructor to make the initial conditions consistent.\n"
38  " For the adjoint integrator, this requires special construction\n"
39  " which has not been implemented yet.\n");
40 
41  adjoint_model_ = createAdjointModel(model_, inputPL);
42  adjoint_integrator_ = integratorBasic<Scalar>(inputPL, adjoint_model_);
43 }
44 
45 template<class Scalar>
48 {
49  state_integrator_ = integratorBasic<Scalar>();
50  adjoint_integrator_ = integratorBasic<Scalar>();
51 }
52 
53 template<class Scalar>
54 bool
57 {
58  const Scalar tfinal =
59  state_integrator_->getTimeStepControl()->getFinalTime();
60  return advanceTime(tfinal);
61 }
62 
63 template<class Scalar>
64 bool
66 advanceTime(const Scalar timeFinal)
67 {
68  using Teuchos::RCP;
69  using Teuchos::rcp_dynamic_cast;
70  using Thyra::VectorBase;
71  using Thyra::VectorSpaceBase;
72  using Thyra::MultiVectorBase;
73  using Thyra::LinearOpBase;
74  using Thyra::LinearOpWithSolveBase;
75  using Thyra::createMember;
76  using Thyra::createMembers;
77  using Thyra::assign;
78  typedef Thyra::ModelEvaluatorBase MEB;
79  typedef Thyra::DefaultMultiVectorProductVector<Scalar> DMVPV;
80  typedef Thyra::DefaultProductVector<Scalar> DPV;
81 
82  // Get initial state for later
83  RCP<const SolutionHistory<Scalar> > state_solution_history =
84  state_integrator_->getSolutionHistory();
85  RCP<const SolutionState<Scalar> > initial_state =
86  (*state_solution_history)[0];
87 
88  // Run state integrator and get solution
89  bool state_status = state_integrator_->advanceTime(timeFinal);
90 
91  // Set solution history in adjoint stepper
92  adjoint_model_->setForwardSolutionHistory(state_solution_history);
93 
94  // Compute dg/dx
95  RCP<const VectorSpaceBase<Scalar> > g_space = model_->get_g_space(g_index_);
96  RCP<const VectorSpaceBase<Scalar> > x_space = model_->get_x_space();
97  const int num_g = g_space->dim();
98  RCP<MultiVectorBase<Scalar> > dgdx = createMembers(x_space, num_g);
99  MEB::InArgs<Scalar> inargs = model_->getNominalValues();
100  RCP<const SolutionState<Scalar> > state =
101  state_solution_history->getCurrentState();
102  inargs.set_t(state->getTime());
103  inargs.set_x(state->getX());
104  inargs.set_x_dot(state->getXDot());
105  MEB::OutArgs<Scalar> outargs = model_->createOutArgs();
106  outargs.set_DgDx(g_index_,
107  MEB::Derivative<Scalar>(dgdx, MEB::DERIV_MV_GRADIENT_FORM));
108  model_->evalModel(inargs, outargs);
109  outargs.set_DgDx(g_index_, MEB::Derivative<Scalar>());
110 
111  // Compute ICs == [ (df/dx_dot)^{-T} (dg/dx)^T; 0 ]
112  // For explicit form, we are relying on the user to inform us the
113  // the mass matrix is the identity. It would be nice to be able to determine
114  // somehow automatically that we are using an explicit stepper.
115  RCP<DPV> adjoint_init =
116  rcp_dynamic_cast<DPV>(Thyra::createMember(adjoint_model_->get_x_space()));
117  RCP<MultiVectorBase<Scalar> > adjoint_init_mv =
118  rcp_dynamic_cast<DMVPV>(adjoint_init->getNonconstVectorBlock(0))->getNonconstMultiVector();
119  assign(adjoint_init->getNonconstVectorBlock(1).ptr(),
120  Teuchos::ScalarTraits<Scalar>::zero());
121  if (mass_matrix_is_identity_)
122  assign(adjoint_init_mv.ptr(), *dgdx);
123  else {
124  inargs.set_alpha(1.0);
125  inargs.set_beta(0.0);
126  RCP<LinearOpWithSolveBase<Scalar> > W = model_->create_W();
127  outargs.set_W(W);
128  model_->evalModel(inargs, outargs);
129  W->solve(Thyra::CONJTRANS, *dgdx, adjoint_init_mv.ptr());
130  outargs.set_W(Teuchos::null);
131  }
132 
133  // Run sensitivity integrator and get solution
134  adjoint_integrator_->initializeSolutionHistory(Scalar(0.0), adjoint_init);
135  bool sens_status = adjoint_integrator_->advanceTime(timeFinal);
136  RCP<const SolutionHistory<Scalar> > adjoint_solution_history =
137  adjoint_integrator_->getSolutionHistory();
138 
139  // Compute dg/dp at final time T
140  RCP<const VectorSpaceBase<Scalar> > p_space = model_->get_p_space(p_index_);
141  dgdp_ = createMembers(p_space, num_g);
142  if (g_depends_on_p_) {
143  MEB::DerivativeSupport dgdp_support =
144  outargs.supports(MEB::OUT_ARG_DgDp, g_index_, p_index_);
145  if (dgdp_support.supports(MEB::DERIV_MV_GRADIENT_FORM)) {
146  outargs.set_DgDp(g_index_, p_index_,
147  MEB::Derivative<Scalar>(dgdp_,
148  MEB::DERIV_MV_GRADIENT_FORM));
149  model_->evalModel(inargs, outargs);
150  }
151  else if (dgdp_support.supports(MEB::DERIV_MV_JACOBIAN_FORM)) {
152  const int num_p = p_space->dim();
153  RCP<MultiVectorBase<Scalar> > dgdp_trans =
154  createMembers(g_space, num_p);
155  outargs.set_DgDp(g_index_, p_index_,
156  MEB::Derivative<Scalar>(dgdp_trans,
157  MEB::DERIV_MV_JACOBIAN_FORM));
158  model_->evalModel(inargs, outargs);
159  Thyra::DetachedMultiVectorView<Scalar> dgdp_view(*dgdp_);
160  Thyra::DetachedMultiVectorView<Scalar> dgdp_trans_view(*dgdp_trans);
161  for (int i=0; i<num_p; ++i)
162  for (int j=0; j<num_g; ++j)
163  dgdp_view(i,j) = dgdp_trans_view(j,i);
164  }
165  else
166  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
167  "Invalid dg/dp support");
168  outargs.set_DgDp(g_index_, p_index_, MEB::Derivative<Scalar>());
169  }
170  else
171  assign(dgdp_.ptr(), Scalar(0.0));
172 
173  // Add in initial condition term = (dx/dp^T(0))*(df/dx_dot^T(0))*y(0)
174  // If dxdp_init_ is null, assume it is zero
175  if (ic_depends_on_p_ && dxdp_init_ != Teuchos::null) {
176  RCP<const SolutionState<Scalar> > adjoint_state =
177  adjoint_solution_history->getCurrentState();
178  RCP<const VectorBase<Scalar> > adjoint_x =
179  rcp_dynamic_cast<const DPV>(adjoint_state->getX())->getVectorBlock(0);
180  RCP<const MultiVectorBase<Scalar> > adjoint_mv =
181  rcp_dynamic_cast<const DMVPV>(adjoint_x)->getMultiVector();
182  if (mass_matrix_is_identity_)
183  dxdp_init_->apply(Thyra::CONJTRANS, *adjoint_mv, dgdp_.ptr(), Scalar(1.0),
184  Scalar(1.0));
185  else {
186  inargs.set_t(initial_state->getTime());
187  inargs.set_x(initial_state->getX());
188  inargs.set_x_dot(initial_state->getXDot());
189  inargs.set_alpha(1.0);
190  inargs.set_beta(0.0);
191  RCP<LinearOpBase<Scalar> > W_op = model_->create_W_op();
192  outargs.set_W_op(W_op);
193  model_->evalModel(inargs, outargs);
194  outargs.set_W_op(Teuchos::null);
195  RCP<MultiVectorBase<Scalar> > tmp = createMembers(x_space, num_g);
196  W_op->apply(Thyra::CONJTRANS, *adjoint_mv, tmp.ptr(), Scalar(1.0),
197  Scalar(0.0));
198  dxdp_init_->apply(Thyra::CONJTRANS, *tmp, dgdp_.ptr(), Scalar(1.0),
199  Scalar(1.0));
200  }
201  }
202 
203  // Add in model parameter term = \int_0^T( (df/dp^T(t)*y(t) )dt which
204  // is computed during the adjoint integration as an auxiliary integral
205  // (2nd block of the solution vector)
206  if (f_depends_on_p_) {
207  RCP<const SolutionState<Scalar> > adjoint_state =
208  adjoint_solution_history->getCurrentState();
209  RCP<const VectorBase<Scalar> > z =
210  rcp_dynamic_cast<const DPV>(adjoint_state->getX())->getVectorBlock(1);
211  RCP<const MultiVectorBase<Scalar> > z_mv =
212  rcp_dynamic_cast<const DMVPV>(z)->getMultiVector();
213  Thyra::V_VmV(dgdp_.ptr(), *dgdp_, *z_mv);
214  }
215 
216  buildSolutionHistory(state_solution_history, adjoint_solution_history);
217 
218  return state_status && sens_status;
219 }
220 
221 template<class Scalar>
222 Scalar
224 getTime() const
225 {
226  return solutionHistory_->getCurrentTime();
227 }
228 
229 template<class Scalar>
230 Scalar
232 getIndex() const
233 {
234  return solutionHistory_->getCurrentIndex();
235 }
236 
237 template<class Scalar>
238 Status
240 getStatus() const
241 {
242  Status state_status = state_integrator_->getStatus();
243  Status sens_status = adjoint_integrator_->getStatus();
244  if (state_status == FAILED || sens_status == FAILED)
245  return FAILED;
246  if (state_status == WORKING || sens_status == WORKING)
247  return WORKING;
248  return PASSED;
249 }
250 
251 template<class Scalar>
252 Teuchos::RCP<Stepper<Scalar> >
254 getStepper() const
255 {
256  return state_integrator_->getStepper();
257 }
258 
259 template<class Scalar>
260 Teuchos::RCP<Teuchos::ParameterList>
263 {
264  return state_integrator_->getTempusParameterList();
265 }
266 
267 template<class Scalar>
268 void
270 setTempusParameterList(Teuchos::RCP<Teuchos::ParameterList> pl)
271 {
272  state_integrator_->setTempusParameterList(pl);
273  adjoint_integrator_->setTempusParameterList(pl);
274 }
275 
276 template<class Scalar>
277 Teuchos::RCP<const SolutionHistory<Scalar> >
280 {
281  return solutionHistory_;
282 }
283 
284 template<class Scalar>
285 Teuchos::RCP<const TimeStepControl<Scalar> >
288 {
289  return state_integrator_->getTimeStepControl();
290 }
291 
292 template<class Scalar>
295  Teuchos::RCP<const Thyra::VectorBase<Scalar> > x0,
296  Teuchos::RCP<const Thyra::VectorBase<Scalar> > xdot0,
297  Teuchos::RCP<const Thyra::VectorBase<Scalar> > xdotdot0,
298  Teuchos::RCP<const Thyra::MultiVectorBase<Scalar> > DxDp0,
299  Teuchos::RCP<const Thyra::MultiVectorBase<Scalar> > /* DxdotDp0 */,
300  Teuchos::RCP<const Thyra::MultiVectorBase<Scalar> > /* DxdotdotDp0 */)
301 {
302  state_integrator_->initializeSolutionHistory(t0, x0, xdot0, xdotdot0);
303  dxdp_init_ = DxDp0;
304 }
305 
306 template<class Scalar>
307 Teuchos::RCP<const Thyra::VectorBase<Scalar> >
309 getX() const
310 {
311  return state_integrator_->getX();
312 }
313 
314 template<class Scalar>
315 Teuchos::RCP<const Thyra::VectorBase<Scalar> >
317 getXdot() const
318 {
319  return state_integrator_->getXdot();
320 }
321 
322 template<class Scalar>
323 Teuchos::RCP<const Thyra::VectorBase<Scalar> >
325 getXdotdot() const
326 {
327  return state_integrator_->getXdotdot();
328 }
329 
330 template<class Scalar>
331 Teuchos::RCP<const Thyra::MultiVectorBase<Scalar> >
333 getDgDp() const
334 {
335  return dgdp_;
336 }
337 
338 template<class Scalar>
339 std::string
341 description() const
342 {
343  std::string name = "Tempus::IntegratorAdjointSensitivity";
344  return(name);
345 }
346 
347 template<class Scalar>
348 void
351  Teuchos::FancyOStream &out,
352  const Teuchos::EVerbosityLevel verbLevel) const
353 {
354  out << description() << "::describe" << std::endl;
355  state_integrator_->describe(out, verbLevel);
356  adjoint_integrator_->describe(out, verbLevel);
357 }
358 
359 template<class Scalar>
360 void
362 setParameterList(const Teuchos::RCP<Teuchos::ParameterList> & inputPL)
363 {
364  if (state_integrator_ != Teuchos::null)
365  state_integrator_->setParameterList(inputPL);
366  if (adjoint_integrator_ != Teuchos::null)
367  adjoint_integrator_->setParameterList(inputPL);
368  Teuchos::ParameterList& spl = inputPL->sublist("Sensitivities");
369  p_index_ = spl.get<int>("Sensitivity Parameter Index", 0);
370  g_index_ = spl.get<int>("Response Function Index", 0);
371  g_depends_on_p_ = spl.get<bool>("Response Depends on Parameters", true);
372  f_depends_on_p_ = spl.get<bool>("Residual Depends on Parameters", true);
373  ic_depends_on_p_ = spl.get<bool>("IC Depends on Parameters", true);
374  mass_matrix_is_identity_ = spl.get<bool>("Mass Matrix Is Identity", false);
375 }
376 
377 template<class Scalar>
378 Teuchos::RCP<Teuchos::ParameterList>
381 {
382  state_integrator_->unsetParameterList();
383  return adjoint_integrator_->unsetParameterList();
384 }
385 
386 template<class Scalar>
387 Teuchos::RCP<const Teuchos::ParameterList>
390 {
391  using Teuchos::RCP;
392  using Teuchos::ParameterList;
393  using Teuchos::parameterList;
394 
395  RCP<ParameterList> pl = parameterList();
396  RCP<const ParameterList> integrator_pl =
397  state_integrator_->getValidParameters();
398  RCP<const ParameterList> sensitivity_pl =
400  pl->setParameters(*integrator_pl);
401  ParameterList& spl = pl->sublist("Sensitivities");
402  spl.setParameters(*sensitivity_pl);
403  spl.set<bool>("Response Depends on Parameters", true);
404  spl.set<bool>("Residual Depends on Parameters", true);
405  spl.set<bool>("IC Depends on Parameters", true);
406 
407  return pl;
408 }
409 
410 template<class Scalar>
411 Teuchos::RCP<Teuchos::ParameterList>
414 {
415  return state_integrator_->getNonconstParameterList();
416 }
417 
418 template <class Scalar>
419 Teuchos::RCP<Tempus::AdjointAuxSensitivityModelEvaluator<Scalar> >
422  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& model,
423  const Teuchos::RCP<Teuchos::ParameterList>& inputPL)
424 {
425  using Teuchos::RCP;
426  using Teuchos::rcp;
427  using Teuchos::ParameterList;
428 
429  RCP<ParameterList> spl = Teuchos::parameterList();
430  if (inputPL != Teuchos::null) {
431  *spl = inputPL->sublist("Sensitivities");
432  }
433  spl->remove("Response Depends on Parameters");
434  spl->remove("Residual Depends on Parameters");
435  spl->remove("IC Depends on Parameters");
436  const Scalar tfinal = state_integrator_->getTimeStepControl()->getFinalTime();
438  model, tfinal, spl));
439 }
440 
441 template<class Scalar>
442 void
445  const Teuchos::RCP<const SolutionHistory<Scalar> >& state_solution_history,
446  const Teuchos::RCP<const SolutionHistory<Scalar> >& adjoint_solution_history)
447 {
448  using Teuchos::RCP;
449  using Teuchos::rcp;
450  using Teuchos::rcp_dynamic_cast;
451  using Teuchos::ParameterList;
452  using Thyra::VectorBase;
453  using Thyra::MultiVectorBase;
454  using Thyra::VectorSpaceBase;
455  using Thyra::createMembers;
456  using Thyra::multiVectorProductVector;
457  using Thyra::assign;
458  typedef Thyra::DefaultProductVectorSpace<Scalar> DPVS;
459  typedef Thyra::DefaultProductVector<Scalar> DPV;
460 
461  // Create combined solution histories combining the forward and adjoint
462  // solutions. We do not include the auxiliary part from the adjoint solution.
463  RCP<ParameterList> shPL =
464  Teuchos::sublist(state_integrator_->getIntegratorParameterList(),
465  "Solution History", true);
466  solutionHistory_ = rcp(new SolutionHistory<Scalar>(shPL));
467 
468  RCP<const VectorSpaceBase<Scalar> > x_space = model_->get_x_space();
469  RCP<const VectorSpaceBase<Scalar> > adjoint_space =
470  rcp_dynamic_cast<const DPVS>(adjoint_model_->get_x_space())->getBlock(0);
471  Teuchos::Array< RCP<const VectorSpaceBase<Scalar> > > spaces(2);
472  spaces[0] = x_space;
473  spaces[1] = adjoint_space;
474  RCP<const DPVS > prod_space = Thyra::productVectorSpace(spaces());
475 
476  int num_states = state_solution_history->getNumStates();
477  const Scalar t_final = state_integrator_->getTime();
478  for (int i=0; i<num_states; ++i) {
479  RCP<const SolutionState<Scalar> > forward_state =
480  (*state_solution_history)[i];
481  RCP<const SolutionState<Scalar> > adjoint_state =
482  adjoint_solution_history->findState(t_final-forward_state->getTime());
483 
484  // X
485  RCP<DPV> x = Thyra::defaultProductVector(prod_space);
486  RCP<const VectorBase<Scalar> > adjoint_x =
487  rcp_dynamic_cast<const DPV>(adjoint_state->getX())->getVectorBlock(0);
488  assign(x->getNonconstVectorBlock(0).ptr(), *(forward_state->getX()));
489  assign(x->getNonconstVectorBlock(1).ptr(), *(adjoint_x));
490  RCP<VectorBase<Scalar> > x_b = x;
491 
492  // X-Dot
493  RCP<DPV> x_dot = Thyra::defaultProductVector(prod_space);
494  RCP<const VectorBase<Scalar> > adjoint_x_dot =
495  rcp_dynamic_cast<const DPV>(adjoint_state->getXDot())->getVectorBlock(0);
496  assign(x_dot->getNonconstVectorBlock(0).ptr(), *(forward_state->getXDot()));
497  assign(x_dot->getNonconstVectorBlock(1).ptr(), *(adjoint_x_dot));
498  RCP<VectorBase<Scalar> > x_dot_b = x_dot;
499 
500  // X-Dot-Dot
501  RCP<DPV> x_dot_dot;
502  if (forward_state->getXDotDot() != Teuchos::null) {
503  x_dot_dot = Thyra::defaultProductVector(prod_space);
504  RCP<const VectorBase<Scalar> > adjoint_x_dot_dot =
505  rcp_dynamic_cast<const DPV>(
506  adjoint_state->getXDotDot())->getVectorBlock(0);
507  assign(x_dot_dot->getNonconstVectorBlock(0).ptr(),
508  *(forward_state->getXDotDot()));
509  assign(x_dot_dot->getNonconstVectorBlock(1).ptr(),
510  *(adjoint_x_dot_dot));
511  }
512  RCP<VectorBase<Scalar> > x_dot_dot_b = x_dot_dot;
513 
514  RCP<SolutionState<Scalar> > prod_state =
515  rcp(new SolutionState<Scalar>(forward_state->getMetaData()->clone(),
516  x_b, x_dot_b, x_dot_dot_b,
517  forward_state->getStepperState()->clone(),
518  Teuchos::null));
519  solutionHistory_->addState(prod_state);
520  }
521 }
522 
523 /// Non-member constructor
524 template<class Scalar>
525 Teuchos::RCP<Tempus::IntegratorAdjointSensitivity<Scalar> >
527  Teuchos::RCP<Teuchos::ParameterList> pList,
528  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& model)
529 {
530  Teuchos::RCP<Tempus::IntegratorAdjointSensitivity<Scalar> > integrator =
531  Teuchos::rcp(new Tempus::IntegratorAdjointSensitivity<Scalar>(pList, model));
532  return(integrator);
533 }
534 
535 /// Non-member constructor
536 template<class Scalar>
537 Teuchos::RCP<Tempus::IntegratorAdjointSensitivity<Scalar> >
539 {
540  Teuchos::RCP<Tempus::IntegratorAdjointSensitivity<Scalar> > integrator =
542  return(integrator);
543 }
544 
545 } // namespace Tempus
546 #endif // Tempus_IntegratorAdjointSensitivity_impl_hpp
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const override
Teuchos::RCP< Tempus::IntegratorAdjointSensitivity< Scalar > > integratorAdjointSensitivity(Teuchos::RCP< Teuchos::ParameterList > pList, const Teuchos::RCP< Thyra::ModelEvaluator< Scalar > > &model)
Non-member constructor.
Teuchos::RCP< Tempus::AdjointAuxSensitivityModelEvaluator< Scalar > > createAdjointModel(const Teuchos::RCP< Thyra::ModelEvaluator< Scalar > > &model, const Teuchos::RCP< Teuchos::ParameterList > &inputPL)
Teuchos::RCP< Teuchos::ParameterList > getNonconstParameterList() override
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
static Teuchos::RCP< const Teuchos::ParameterList > getValidParameters()
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getXdot() const
Get current the time derivative of the solution, xdot.
ModelEvaluator for forming adjoint sensitivity equations.
Teuchos::RCP< Teuchos::ParameterList > unsetParameterList() override
virtual Scalar getIndex() const override
Get current index.
void buildSolutionHistory(const Teuchos::RCP< const SolutionHistory< Scalar > > &state_solution_history, const Teuchos::RCP< const SolutionHistory< Scalar > > &adjoint_solution_history)
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getX() const
Get current the solution, x.
virtual Teuchos::RCP< const SolutionHistory< Scalar > > getSolutionHistory() const override
Get the SolutionHistory.
Status
Status for the Integrator, the Stepper and the SolutionState.
virtual Teuchos::RCP< Teuchos::ParameterList > getTempusParameterList() override
Return a copy of the Tempus ParameterList.
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.
virtual void setTempusParameterList(Teuchos::RCP< Teuchos::ParameterList > pl) override
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.
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > &pl) override
Time integrator suitable for adjoint sensitivity analysis.
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getXdotdot() const
Get current the second time derivative of the solution, xdotdot.
virtual Scalar getTime() const override
Get current time.
virtual Status getStatus() const override
Get Status.
Solution state for integrators and steppers. SolutionState contains the metadata for solutions and th...
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.