Tempus  Version of the Day
Time Integration
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Tempus_IntegratorPseudoTransientAdjointSensitivity_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_IntegratorPseudoTransientAdjointSensitivity_impl_hpp
10 #define Tempus_IntegratorPseudoTransientAdjointSensitivity_impl_hpp
11 
12 #include "Thyra_DefaultMultiVectorProductVector.hpp"
13 #include "Thyra_VectorStdOps.hpp"
14 #include "Thyra_MultiVectorStdOps.hpp"
16 
17 
18 namespace Tempus {
19 
20 template<class Scalar>
25  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& adjoint_residual_model,
26  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& adjoint_solve_model)
27 {
28  model_ = model;
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_);
35  stepMode_ = SensitivityStepMode::Forward;
36  do_forward_integration_ = true;
37  do_adjoint_integration_ = true;
38 }
39 
40 template<class Scalar>
44  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& adjoint_residual_model,
45  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& adjoint_solve_model,
46  std::string stepperType)
47 {
48  model_ = model;
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);
55  stepMode_ = SensitivityStepMode::Forward;
56  do_forward_integration_ = true;
57  do_adjoint_integration_ = true;
58 }
59 
60 template<class Scalar>
65  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& adjoint_model) :
67  inputPL, model, adjoint_model, adjoint_model)
68 {
69 }
70 
71 template<class Scalar>
75  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& adjoint_model,
76  std::string stepperType) :
78  model, adjoint_model, adjoint_model, stepperType)
79 {
80 }
81 
82 template<class Scalar>
88  inputPL, model, Thyra::implicitAdjointModelEvaluator(model))
89 {
90 }
91 
92 template<class Scalar>
96  std::string stepperType) :
98  model, Thyra::implicitAdjointModelEvaluator(model), stepperType)
99 {
100 }
101 
102 template<class Scalar>
105 {
106  state_integrator_ = createIntegratorBasic<Scalar>();
107  sens_integrator_ = createIntegratorBasic<Scalar>();
108  stepMode_ = SensitivityStepMode::Forward;
109 }
110 
111 template<class Scalar>
112 bool
115 {
116  const Scalar tfinal =
117  state_integrator_->getTimeStepControl()->getFinalTime();
118  return advanceTime(tfinal);
119 }
120 
121 template<class Scalar>
122 bool
124 advanceTime(const Scalar timeFinal)
125 {
126  TEMPUS_FUNC_TIME_MONITOR_DIFF("Tempus::IntegratorPseudoTransientAdjointSensitivity::advanceTime()", TEMPUS_PTAS_AT);
127 
128  using Teuchos::RCP;
129  using Thyra::VectorBase;
130  typedef Thyra::ModelEvaluatorBase MEB;
131 
132  bool state_status = true;
133  if (do_forward_integration_) {
134  TEMPUS_FUNC_TIME_MONITOR_DIFF("Tempus::IntegratorPseudoTransientAdjointSensitivity::advanceTime::state", TEMPUS_PTAS_AT_FWD);
135 
136  // Run state integrator and get solution
137  stepMode_ = SensitivityStepMode::Forward;
138  state_status = state_integrator_->advanceTime(timeFinal);
139  }
140 
141  bool sens_status = true;
142  if (do_adjoint_integration_) {
143  TEMPUS_FUNC_TIME_MONITOR_DIFF("Tempus::IntegratorPseudoTransientAdjointSensitivity::advanceTime::adjoint", TEMPUS_PTAS_AT_ADJ);
144 
145  // For at least some time-stepping methods, the time of the last time step
146  // may not be timeFinal (e.g., it may be greater by at most delta_t).
147  // But since the adjoint model requires timeFinal in its formulation, reset
148  // it to the achieved final time.
149  sens_model_->setFinalTime(state_integrator_->getTime());
150 
151  // Set solution in sensitivity ME
152  sens_model_->setForwardSolutionHistory(
153  state_integrator_->getSolutionHistory());
154 
155  // Run sensitivity integrator
156  stepMode_ = SensitivityStepMode::Adjoint;
157  sens_status = sens_integrator_->advanceTime(timeFinal);
158 
159  // Compute final dg/dp, g which is computed by response 0, 1 of the adjoint
160  // model evaluator
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);
173  }
174  if (g_ == Teuchos::null)
175  g_ = Thyra::createMember(sens_model_->get_g_space(1));
176  outargs.set_g(0, G);
177  outargs.set_g(1, g_);
178  sens_model_->evalModel(inargs, outargs);
179 
180  buildSolutionHistory();
181  }
182 
183  return state_status && sens_status;
184 }
185 
186 template<class Scalar>
187 Scalar
189 getTime() const
190 {
191  return solutionHistory_->getCurrentTime();
192 }
193 
194 template<class Scalar>
195 int
197 getIndex() const
198 {
199  return solutionHistory_->getCurrentIndex();
200 }
201 
202 template<class Scalar>
203 Status
205 getStatus() const
206 {
207  Status state_status = state_integrator_->getStatus();
208  Status sens_status = sens_integrator_->getStatus();
209  if (state_status == FAILED || sens_status == FAILED)
210  return FAILED;
211  if (state_status == WORKING || sens_status == WORKING)
212  return WORKING;
213  return PASSED;
214 }
215 
216 template <class Scalar>
218  const Status st) {
219  state_integrator_->setStatus(st);
220  sens_integrator_->setStatus(st);
221 }
222 
223 template<class Scalar>
226 getStepper() const
227 {
228  return state_integrator_->getStepper();
229 }
230 
231 template<class Scalar>
235 {
236  return state_integrator_->getStepper();
237 }
238 
239 template<class Scalar>
243 {
244  return sens_integrator_->getStepper();
245 }
246 
247 template<class Scalar>
251 {
252  return solutionHistory_;
253 }
254 
255 template<class Scalar>
259 {
260  return state_integrator_->getSolutionHistory();
261 }
262 
263 template<class Scalar>
267 {
268  return sens_integrator_->getSolutionHistory();
269 }
270 
271 template<class Scalar>
275 {
276  return solutionHistory_;
277 }
278 
279 template<class Scalar>
283 {
284  return state_integrator_->getTimeStepControl();
285 }
286 
287 template<class Scalar>
291 {
292  return state_integrator_->getNonConstTimeStepControl();
293 }
294 
295 template<class Scalar>
299 {
300  return state_integrator_->getNonConstTimeStepControl();
301 }
302 
303 template<class Scalar>
307 {
308  return sens_integrator_->getNonConstTimeStepControl();
309 }
310 
311 template<class Scalar>
315 {
316  return state_integrator_->getObserver();
317 }
318 
319 template<class Scalar>
320 void
323 {
324  state_integrator_->setObserver(obs);
325  sens_integrator_->setObserver(obs);
326 }
327 
328 template<class Scalar>
333  Teuchos::RCP<const Thyra::VectorBase<Scalar> > xdotdot0,
337 {
338  using Teuchos::RCP;
339  using Teuchos::rcp_dynamic_cast;
341  using Thyra::assign;
342  using Thyra::createMember;
343 
344  //
345  // Create and initialize product X, Xdot, Xdotdot
346 
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();
352 
353  // y
354  if (y0 == Teuchos::null)
355  assign(Y->getNonconstMultiVector().ptr(), zero);
356  else
357  assign(Y->getNonconstMultiVector().ptr(), *y0);
358 
359  // ydot
360  if (ydot0 == Teuchos::null)
361  assign(Ydot->getNonconstMultiVector().ptr(), zero);
362  else
363  assign(Ydot->getNonconstMultiVector().ptr(), *ydot0);
364 
365  // ydotdot
366  if (ydotdot0 == Teuchos::null)
367  assign(Ydotdot->getNonconstMultiVector().ptr(), zero);
368  else
369  assign(Ydotdot->getNonconstMultiVector().ptr(), *ydotdot0);
370 
371  state_integrator_->initializeSolutionHistory(t0, x0, xdot0, xdotdot0);
372  sens_integrator_->initializeSolutionHistory(t0, Y, Ydot, Ydotdot);
373 }
374 
375 template<class Scalar>
378 getX() const
379 {
380  return state_integrator_->getX();
381 }
382 
383 template<class Scalar>
386 getXDot() const
387 {
388  return state_integrator_->getXDot();
389 }
390 
391 template<class Scalar>
394 getXDotDot() const
395 {
396  return state_integrator_->getXDotDot();
397 }
398 
399 template<class Scalar>
402 getY() const
403 {
404  using Teuchos::RCP;
405  using Teuchos::rcp_dynamic_cast;
406  RCP<const DMVPV> mvpv =
407  rcp_dynamic_cast<const DMVPV>(sens_integrator_->getX());
408  return mvpv->getMultiVector();
409 }
410 
411 template<class Scalar>
414 getYDot() const
415 {
416  using Teuchos::RCP;
417  using Teuchos::rcp_dynamic_cast;
418  RCP<const DMVPV> mvpv =
419  rcp_dynamic_cast<const DMVPV>(sens_integrator_->getXDot());
420  return mvpv->getMultiVector();
421 }
422 
423 template<class Scalar>
426 getYDotDot() const
427 {
428  using Teuchos::RCP;
429  using Teuchos::rcp_dynamic_cast;
430  RCP<const DMVPV> mvpv =
431  rcp_dynamic_cast<const DMVPV>(sens_integrator_->getXDotDot());
432  return mvpv->getMultiVector();
433 }
434 
435 template<class Scalar>
438 getG() const
439 {
440  return g_;
441 }
442 
443 template<class Scalar>
446 getDgDp() const
447 {
448  return dgdp_->getMultiVector();
449 }
450 
451 template<class Scalar>
452 std::string
454 description() const
455 {
456  std::string name = "Tempus::IntegratorPseudoTransientAdjointSensitivity";
457  return(name);
458 }
459 
460 template<class Scalar>
461 void
465  const Teuchos::EVerbosityLevel verbLevel) const
466 {
467  auto l_out = Teuchos::fancyOStream( out.getOStream() );
468  Teuchos::OSTab ostab(*l_out, 2, this->description());
469  l_out->setOutputToRootOnly(0);
470 
471  *l_out << description() << "::describe" << std::endl;
472  state_integrator_->describe(*l_out, verbLevel);
473  sens_integrator_->describe(*l_out, verbLevel);
474 }
475 
476 template<class Scalar>
477 void
480 {
481  // IntegratorBasic is no longer a Teuchos::ParameterListAcceptor.
482  // Since setting the ParameterList is essentially a complete reset,
483  // we will rebuild from scratch and reuse the ModelEvaluator.
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);
488 
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);
493 }
494 
495 template<class Scalar>
499 {
500  // IntegratorBasic is no longer a Teuchos::ParameterListAcceptor.
501  // We will treat unsetting the ParameterList as a reset to default
502  // settings, and reuse the ModelEvaluator.
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);
507 
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);
512 
513  auto pl = Teuchos::rcp_const_cast<Teuchos::ParameterList> (
514  sens_integrator_->getValidParameters());
515  return pl;
516 }
517 
518 template<class Scalar>
522 {
526  state_integrator_->getValidParameters();
529  pl->setParameters(*integrator_pl);
530  pl->sublist("Sensitivities").setParameters(*sensitivity_pl);
531 
532  return pl;
533 }
534 
535 template<class Scalar>
539 {
540  auto pl = Teuchos::rcp_const_cast<Teuchos::ParameterList> (
541  state_integrator_->getValidParameters());
542  return pl;
543 }
544 
545 template<class Scalar>
548 getStepMode() const
549 {
550  return stepMode_;
551 }
552 
553 template <class Scalar>
558  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& adjoint_residual_model,
559  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& adjoint_solve_model,
561 {
562  using Teuchos::rcp;
563 
564  Teuchos::RCP<Teuchos::ParameterList> pl = Teuchos::parameterList();
565  if (inputPL != Teuchos::null) {
566  *pl = inputPL->sublist("Sensitivities");
567  }
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));
573 }
574 
575 template<class Scalar>
576 void
579 {
580  using Teuchos::RCP;
581  using Teuchos::rcp;
582  using Teuchos::rcp_dynamic_cast;
584  using Thyra::VectorBase;
586  using Thyra::assign;
587  using Thyra::defaultProductVector;
590 
591  // Create combined solution histories, first for the states with zero
592  // adjoint and then for the adjoint with frozen states
593  auto shPL = Teuchos::rcp_const_cast<Teuchos::ParameterList> (
594  state_integrator_->getSolutionHistory()->getValidParameters());
595  solutionHistory_ = createSolutionHistoryPL<Scalar>(shPL);
596 
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();
602  spaces[0] = x_space;
603  spaces[1] = adjoint_space;
604  RCP<const DPVS > prod_space = Thyra::productVectorSpace(spaces());
605  const Scalar zero = Teuchos::ScalarTraits<Scalar>::zero();
606 
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];
612 
613  // X
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;
618 
619  // X-Dot
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);
625  x_dot_b = x_dot;
626  }
627 
628  // X-Dot-Dot
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;
635  }
636 
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);
642  }
643 
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];
655 
656  // X
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;
661 
662  // X-Dot
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()));
668  x_dot_b = x_dot;
669  }
670 
671  // X-Dot-Dot
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;
678  }
679 
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);
685  }
686 }
687 
689 template<class Scalar>
694 {
697  return(integrator);
698 }
699 
701 template<class Scalar>
705  std::string stepperType)
706 {
709  return(integrator);
710 }
711 
713 template<class Scalar>
718  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& adjoint_model)
719 {
721  Teuchos::rcp(new IntegratorPseudoTransientAdjointSensitivity<Scalar>(pList, model, adjoint_model));
722  return(integrator);
723 }
724 
726 template<class Scalar>
730  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& adjoint_model,
731  std::string stepperType)
732 {
734  Teuchos::rcp(new IntegratorPseudoTransientAdjointSensitivity<Scalar>(model, adjoint_model, stepperType));
735  return(integrator);
736 }
737 
739 template<class Scalar>
744  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& adjoint_residual_model,
745  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& adjoint_solve_model)
746 {
748  Teuchos::rcp(new IntegratorPseudoTransientAdjointSensitivity<Scalar>(pList, model, adjoint_residual_model, adjoint_solve_model));
749  return(integrator);
750 }
751 
753 template<class Scalar>
757  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& adjoint_residual_model,
758  const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& adjoint_solve_model,
759  std::string stepperType)
760 {
762  Teuchos::rcp(new IntegratorPseudoTransientAdjointSensitivity<Scalar>(model, adjoint_residual_model, adjoint_solve_model, stepperType));
763  return(integrator);
764 }
765 
767 template<class Scalar>
770 {
773  return(integrator);
774 }
775 
776 } // namespace Tempus
777 #endif // Tempus_IntegratorPseudoTransientAdjointSensitivity_impl_hpp
static Teuchos::RCP< const Teuchos::ParameterList > getValidParameters()
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > &pl) override
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const override
virtual Teuchos::RCP< IntegratorObserver< Scalar > > getObserver()
Get the Observer.
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getY() const
Get the current adjoint solution, y.
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getXDot() const
Get the current time derivative of the solution, xdot.
Teuchos::RCP< Tempus::IntegratorPseudoTransientAdjointSensitivity< Scalar > > integratorPseudoTransientAdjointSensitivity(Teuchos::RCP< Teuchos::ParameterList > pList, const Teuchos::RCP< Thyra::ModelEvaluator< Scalar > > &model)
Nonmember constructor.
virtual void setObserver(Teuchos::RCP< IntegratorObserver< Scalar > > obs=Teuchos::null)
Set the Observer.
virtual Teuchos::RCP< SolutionHistory< Scalar > > getNonConstSolutionHistory() override
Get the SolutionHistory.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getDgDp() const
Return adjoint sensitivity stored in gradient format.
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getYDot() const
Get the current time derivative of the adjoint solution, ydot.
Time integrator suitable for pseudotransient adjoint sensitivity analysis.
SensitivityStepMode getStepMode() const
What mode the current time integration step is in.
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getG() const
Return response function g.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
virtual Teuchos::RCP< TimeStepControl< Scalar > > getNonConstTimeStepControl() override
Status
Status for the Integrator, the Stepper and the SolutionState.
IntegratorObserver class for time integrators.
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getXDotDot() const
Get the current second time derivative of the solution, xdotdot.
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)
RCP< const MultiVectorBase< Scalar > > getMultiVector() const
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getYDotDot() const
Get the current second time derivative of the adjoint solution, ydotdot.
ParameterList & setParameters(const ParameterList &source)
RCP< ImplicitAdjointModelEvaluator< Scalar > > implicitAdjointModelEvaluator(const RCP< const ModelEvaluator< Scalar > > &model)
RCP< std::basic_ostream< char_type, traits_type > > getOStream()
virtual Teuchos::RCP< Stepper< Scalar > > getStepper() const override
Get the Stepper.
virtual Teuchos::RCP< const SolutionHistory< Scalar > > getSolutionHistory() const override
Get the SolutionHistory.
virtual Teuchos::RCP< const TimeStepControl< Scalar > > getTimeStepControl() const override
Get the TimeStepControl.
ModelEvaluator for forming adjoint sensitivity equations.
ParameterList & sublist(const std::string &name, bool mustAlreadyExist=false, const std::string &docString="")
virtual bool advanceTime()
Advance the solution to timeMax, and return true if successful.
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 > > getX() const
Get the current solution, x.