Tempus Version of the Day
Time Integration
Loading...
Searching...
No Matches
Tempus_IntegratorPseudoTransientForwardSensitivity_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_IntegratorPseudoTransientForwardSensitivity_impl_hpp
10#define Tempus_IntegratorPseudoTransientForwardSensitivity_impl_hpp
11
12#include "Thyra_DefaultMultiVectorProductVector.hpp"
13#include "Thyra_VectorStdOps.hpp"
14#include "Thyra_MultiVectorStdOps.hpp"
15
17
18
19namespace Tempus {
20
21template <class Scalar>
23 const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> > &model,
24 const Teuchos::RCP<SensitivityModelEvaluatorBase<Scalar> > &sens_model,
25 const Teuchos::RCP<IntegratorBasic<Scalar> > &fwd_integrator,
26 const Teuchos::RCP<IntegratorBasic<Scalar> > &sens_integrator,
27 const bool reuse_solver,
28 const bool force_W_update)
29 : model_(model)
30 , sens_model_(sens_model)
31 , state_integrator_(fwd_integrator)
32 , sens_integrator_(sens_integrator)
33 , reuse_solver_(reuse_solver)
34 , force_W_update_(force_W_update)
35 , stepMode_(SensitivityStepMode::Forward)
36{
37}
38
39template<class Scalar>
42 reuse_solver_(false),
43 force_W_update_(false),
45{
46 state_integrator_ = createIntegratorBasic<Scalar>();
47 sens_integrator_ = createIntegratorBasic<Scalar>();
48}
49
50template<class Scalar>
51bool
54{
55 using Teuchos::RCP;
57
58 // Run state integrator and get solution
60 bool state_status = state_integrator_->advanceTime();
61
62 // Set solution in sensitivity ME
63 sens_model_->setForwardSolutionState(state_integrator_->getCurrentState());
64
65 // Reuse state solver if requested
66 if (reuse_solver_ &&
67 state_integrator_->getStepper()->getSolver() != Teuchos::null) {
68 sens_model_->setSolver(state_integrator_->getStepper()->getSolver(),
69 force_W_update_);
70 }
71
72 // Run sensitivity integrator
74 bool sens_status = sens_integrator_->advanceTime();
75
76 buildSolutionHistory();
77
78 return state_status && sens_status;
79}
80
81template<class Scalar>
82bool
84advanceTime(const Scalar timeFinal)
85{
86 TEMPUS_FUNC_TIME_MONITOR_DIFF("Tempus::IntegratorPseudoTransientForwardSensitivity::advanceTime()", TEMPUS_PTFS_AT);
87
88 using Teuchos::RCP;
90
91 // Run state integrator and get solution
92 bool state_status = true;
93 {
94 TEMPUS_FUNC_TIME_MONITOR_DIFF("Tempus::IntegratorPseudoTransientForwardSensitivity::advanceTime::state", TEMPUS_PTFS_AT_FWD);
95 state_status = state_integrator_->advanceTime(timeFinal);
96 }
97
98 // Set solution in sensitivity ME
99 sens_model_->setForwardSolutionState(state_integrator_->getCurrentState());
100
101 // Reuse state solver if requested
102 if (reuse_solver_ &&
103 state_integrator_->getStepper()->getSolver() != Teuchos::null) {
104 sens_model_->setSolver(state_integrator_->getStepper()->getSolver(),
105 force_W_update_);
106 }
107
108 // Run sensitivity integrator
109 bool sens_status = true;
110 {
111 TEMPUS_FUNC_TIME_MONITOR_DIFF("Tempus::IntegratorPseudoTransientForwardSensitivity::advanceTime::sensitivity", TEMPUS_PTFS_AT_SEN);
112 sens_status = sens_integrator_->advanceTime(timeFinal);
113 }
114
115 buildSolutionHistory();
116
117 return state_status && sens_status;
118}
119
120template<class Scalar>
121Scalar
123getTime() const
124{
125 return solutionHistory_->getCurrentTime();
126}
127
128template<class Scalar>
129int
131getIndex() const
132{
133 return solutionHistory_->getCurrentIndex();
134}
135
136template<class Scalar>
137Status
139getStatus() const
140{
141 Status state_status = state_integrator_->getStatus();
142 Status sens_status = sens_integrator_->getStatus();
143 if (state_status == FAILED || sens_status == FAILED)
144 return FAILED;
145 if (state_status == WORKING || sens_status == WORKING)
146 return WORKING;
147 return PASSED;
148}
149
150template <class Scalar>
152 const Status st) {
153 state_integrator_->setStatus(st);
154 sens_integrator_->setStatus(st);
155}
156
157template<class Scalar>
158Teuchos::RCP<Stepper<Scalar> >
160getStepper() const
161{
162 return state_integrator_->getStepper();
163}
164
165template<class Scalar>
166Teuchos::RCP<Stepper<Scalar> >
168getStateStepper() const
169{
170 return state_integrator_->getStepper();
171}
172
173template<class Scalar>
174Teuchos::RCP<Stepper<Scalar> >
176getSensStepper() const
177{
178 return sens_integrator_->getStepper();
179}
180
181template<class Scalar>
182Teuchos::RCP<const SolutionHistory<Scalar> >
184getSolutionHistory() const
185{
186 return solutionHistory_;
187}
188
189template<class Scalar>
190Teuchos::RCP<const SolutionHistory<Scalar> >
193{
194 return state_integrator_->getSolutionHistory();
195}
196
197template<class Scalar>
198Teuchos::RCP<const SolutionHistory<Scalar> >
201{
202 return sens_integrator_->getSolutionHistory();
203}
204
205template<class Scalar>
206Teuchos::RCP<SolutionHistory<Scalar> >
209{
210 return solutionHistory_;
211}
212
213template<class Scalar>
214Teuchos::RCP<const TimeStepControl<Scalar> >
216getTimeStepControl() const
217{
218 return state_integrator_->getTimeStepControl();
219}
220
221template<class Scalar>
222Teuchos::RCP<TimeStepControl<Scalar> >
225{
226 return state_integrator_->getNonConstTimeStepControl();
227}
228
229template<class Scalar>
230Teuchos::RCP<TimeStepControl<Scalar> >
233{
234 return state_integrator_->getNonConstTimeStepControl();
235}
236
237template<class Scalar>
238Teuchos::RCP<TimeStepControl<Scalar> >
241{
242 return sens_integrator_->getNonConstTimeStepControl();
243}
244
245template<class Scalar>
246Teuchos::RCP<IntegratorObserver<Scalar> >
249{
250 return state_integrator_->getObserver();
251}
252
253template<class Scalar>
254void
256setObserver(Teuchos::RCP<IntegratorObserver<Scalar> > obs)
257{
258 state_integrator_->setObserver(obs);
259 sens_integrator_->setObserver(obs);
260}
261
262template<class Scalar>
265 Teuchos::RCP<const Thyra::VectorBase<Scalar> > x0,
266 Teuchos::RCP<const Thyra::VectorBase<Scalar> > xdot0,
267 Teuchos::RCP<const Thyra::VectorBase<Scalar> > xdotdot0,
268 Teuchos::RCP<const Thyra::MultiVectorBase<Scalar> > DxDp0,
269 Teuchos::RCP<const Thyra::MultiVectorBase<Scalar> > DxdotDp0,
270 Teuchos::RCP<const Thyra::MultiVectorBase<Scalar> > DxdotdotDp0)
271{
272 using Teuchos::RCP;
273 using Teuchos::rcp_dynamic_cast;
274 using Thyra::VectorSpaceBase;
275 using Thyra::assign;
276 using Thyra::createMember;
277 typedef Thyra::DefaultMultiVectorProductVector<Scalar> DMVPV;
278
279 //
280 // Create and initialize product X, Xdot, Xdotdot
281
282 RCP< const VectorSpaceBase<Scalar> > space = sens_model_->get_x_space();
283 RCP<DMVPV> X = rcp_dynamic_cast<DMVPV>(createMember(space));
284 RCP<DMVPV> Xdot = rcp_dynamic_cast<DMVPV>(createMember(space));
285 RCP<DMVPV> Xdotdot = rcp_dynamic_cast<DMVPV>(createMember(space));
286 const Scalar zero = Teuchos::ScalarTraits<Scalar>::zero();
287
288 // x
289 if (DxDp0 == Teuchos::null)
290 assign(X->getNonconstMultiVector().ptr(), zero);
291 else
292 assign(X->getNonconstMultiVector().ptr(), *DxDp0);
293
294 // xdot
295 if (DxdotDp0 == Teuchos::null)
296 assign(Xdot->getNonconstMultiVector().ptr(), zero);
297 else
298 assign(Xdot->getNonconstMultiVector().ptr(), *DxdotDp0);
299
300 // xdotdot
301 if (DxdotDp0 == Teuchos::null)
302 assign(Xdotdot->getNonconstMultiVector().ptr(), zero);
303 else
304 assign(Xdotdot->getNonconstMultiVector().ptr(), *DxdotdotDp0);
305
306 state_integrator_->initializeSolutionHistory(t0, x0, xdot0, xdotdot0);
307 sens_integrator_->initializeSolutionHistory(t0, X, Xdot, Xdotdot);
308}
309
310template<class Scalar>
311Teuchos::RCP<const Thyra::VectorBase<Scalar> >
313getX() const
314{
315 return state_integrator_->getX();
316}
317
318template<class Scalar>
319Teuchos::RCP<const Thyra::MultiVectorBase<Scalar> >
321getDxDp() const
322{
323 using Teuchos::RCP;
324 using Teuchos::rcp_dynamic_cast;
325 typedef Thyra::DefaultMultiVectorProductVector<Scalar> DMVPV;
326
327 RCP<const DMVPV> X =
328 rcp_dynamic_cast<const DMVPV>(sens_integrator_->getX());
329 return X->getMultiVector();
330}
331
332template<class Scalar>
333Teuchos::RCP<const Thyra::VectorBase<Scalar> >
335getXDot() const
336{
337 return state_integrator_->getXDot();
338}
339
340template<class Scalar>
341Teuchos::RCP<const Thyra::MultiVectorBase<Scalar> >
343getDXDotDp() const
344{
345 using Teuchos::RCP;
346 using Teuchos::rcp_dynamic_cast;
347 typedef Thyra::DefaultMultiVectorProductVector<Scalar> DMVPV;
348
349 RCP<const DMVPV> Xdot =
350 rcp_dynamic_cast<const DMVPV>(sens_integrator_->getXDot());
351 return Xdot->getMultiVector();
352}
353
354template<class Scalar>
355Teuchos::RCP<const Thyra::VectorBase<Scalar> >
357getXDotDot() const
358{
359 return state_integrator_->getXDotDot();
360}
361
362template<class Scalar>
363Teuchos::RCP<const Thyra::MultiVectorBase<Scalar> >
365getDXDotDotDp() const
366{
367 using Teuchos::RCP;
368 using Teuchos::rcp_dynamic_cast;
369 typedef Thyra::DefaultMultiVectorProductVector<Scalar> DMVPV;
370
371 RCP<const DMVPV> Xdotdot =
372 rcp_dynamic_cast<const DMVPV>(sens_integrator_->getXDotDot());
373 return Xdotdot->getMultiVector();
374}
375
376template<class Scalar>
377Teuchos::RCP<const Thyra::VectorBase<Scalar> >
379getG() const
380{
381 typedef Thyra::ModelEvaluatorBase MEB;
382
383 // Compute g which is computed by response 1 of the
384 // sensitivity model evaluator
385 MEB::InArgs<Scalar> inargs = sens_model_->getNominalValues();
386 MEB::OutArgs<Scalar> outargs = sens_model_->createOutArgs();
387 inargs.set_t(sens_integrator_->getTime());
388 inargs.set_x(sens_integrator_->getX());
389 if (inargs.supports(MEB::IN_ARG_x_dot))
390 inargs.set_x_dot(sens_integrator_->getXDot());
391 if (inargs.supports(MEB::IN_ARG_x_dot_dot))
392 inargs.set_x_dot_dot(sens_integrator_->getXDotDot());
393
394 Teuchos::RCP<Thyra::VectorBase<Scalar> > g =
395 Thyra::createMember(sens_model_->get_g_space(1));
396 outargs.set_g(1, g);
397
398 sens_model_->evalModel(inargs, outargs);
399 return g;
400}
401
402template<class Scalar>
403Teuchos::RCP<const Thyra::MultiVectorBase<Scalar> >
405getDgDp() const
406{
407 typedef Thyra::ModelEvaluatorBase MEB;
408 typedef Thyra::DefaultMultiVectorProductVector<Scalar> DMVPV;
409
410 // Compute final dg/dp which is computed by response 0 of the
411 // sensitivity model evaluator
412 MEB::InArgs<Scalar> inargs = sens_model_->getNominalValues();
413 MEB::OutArgs<Scalar> outargs = sens_model_->createOutArgs();
414 inargs.set_t(sens_integrator_->getTime());
415 inargs.set_x(sens_integrator_->getX());
416 if (inargs.supports(MEB::IN_ARG_x_dot))
417 inargs.set_x_dot(sens_integrator_->getXDot());
418 if (inargs.supports(MEB::IN_ARG_x_dot_dot))
419 inargs.set_x_dot_dot(sens_integrator_->getXDotDot());
420
421 Teuchos::RCP<Thyra::VectorBase<Scalar> > G =
422 Thyra::createMember(sens_model_->get_g_space(0));
423 Teuchos::RCP<DMVPV> dgdp = Teuchos::rcp_dynamic_cast<DMVPV>(G);
424 outargs.set_g(0, G);
425
426 sens_model_->evalModel(inargs, outargs);
427 return dgdp->getMultiVector();
428}
429
430template<class Scalar>
431std::string
433description() const
434{
435 std::string name = "Tempus::IntegratorPseudoTransientForwardSensitivity";
436 return(name);
437}
438
439template<class Scalar>
440void
443 Teuchos::FancyOStream &out,
444 const Teuchos::EVerbosityLevel verbLevel) const
445{
446 auto l_out = Teuchos::fancyOStream( out.getOStream() );
447 Teuchos::OSTab ostab(*l_out, 2, this->description());
448 l_out->setOutputToRootOnly(0);
449
450 *l_out << description() << "::describe" << std::endl;
451 state_integrator_->describe(*l_out, verbLevel);
452 sens_integrator_->describe(*l_out, verbLevel);
453}
454
455template<class Scalar>
458getStepMode() const
459{
460 return stepMode_;
461}
462
463template<class Scalar>
464void
467{
468 using Teuchos::RCP;
469 using Teuchos::rcp;
470 using Teuchos::rcp_dynamic_cast;
471 using Teuchos::ParameterList;
472 using Thyra::VectorBase;
474 using Thyra::VectorSpaceBase;
475 using Thyra::createMembers;
476 using Thyra::multiVectorProductVector;
477 using Thyra::assign;
478 typedef Thyra::DefaultMultiVectorProductVector<Scalar> DMVPV;
479
480 //TODO: get the solution history PL or create it
481
482 // Create combined solution histories, first for the states with zero
483 // sensitivities and then for the sensitivities with frozen states
484 RCP<ParameterList> shPL;
485 //Teuchos::sublist(state_integrator_->getIntegratorParameterList(), "Solution History", true);
486 solutionHistory_ = createSolutionHistoryPL<Scalar>(shPL);
487
488 const int num_param =
489 rcp_dynamic_cast<const DMVPV>(sens_integrator_->getX())->getMultiVector()->domain()->dim();
490 RCP<const VectorSpaceBase<Scalar> > x_space = model_->get_x_space();
491 RCP<const Thyra::DefaultMultiVectorProductVectorSpace<Scalar> > prod_space =
492 Thyra::multiVectorProductVectorSpace(x_space, num_param+1);
493 const Teuchos::Range1D rng(1,num_param);
494 const Scalar zero = Teuchos::ScalarTraits<Scalar>::zero();
495
496 RCP<const SolutionHistory<Scalar> > state_solution_history =
497 state_integrator_->getSolutionHistory();
498 int num_states = state_solution_history->getNumStates();
499 for (int i=0; i<num_states; ++i) {
500 RCP<const SolutionState<Scalar> > state = (*state_solution_history)[i];
501
502 // X
503 RCP< MultiVectorBase<Scalar> > x_mv =
504 createMembers(x_space, num_param+1);
505 assign(x_mv->col(0).ptr(), *(state->getX()));
506 assign(x_mv->subView(rng).ptr(), zero);
507 RCP<VectorBase<Scalar> > x = multiVectorProductVector(prod_space, x_mv);
508
509 // X-Dot
510 RCP<VectorBase<Scalar> > x_dot;
511 if (state->getXDot() != Teuchos::null) {
512 RCP< MultiVectorBase<Scalar> > x_dot_mv =
513 createMembers(x_space, num_param+1);
514 assign(x_dot_mv->col(0).ptr(), *(state->getXDot()));
515 assign(x_dot_mv->subView(rng).ptr(), zero);
516 x_dot = multiVectorProductVector(prod_space, x_dot_mv);
517 }
518
519 // X-Dot-Dot
520 RCP<VectorBase<Scalar> > x_dot_dot;
521 if (state->getXDotDot() != Teuchos::null) {
522 RCP< MultiVectorBase<Scalar> > x_dot_dot_mv =
523 createMembers(x_space, num_param+1);
524 assign(x_dot_dot_mv->col(0).ptr(), *(state->getXDotDot()));
525 assign(x_dot_dot_mv->subView(rng).ptr(), zero);
526 x_dot_dot = multiVectorProductVector(prod_space, x_dot_dot_mv);
527 }
528
529 RCP<SolutionState<Scalar> > prod_state = state->clone();
530 prod_state->setX(x);
531 prod_state->setXDot(x_dot);
532 prod_state->setXDotDot(x_dot_dot);
533 solutionHistory_->addState(prod_state);
534 }
535
536 RCP<const VectorBase<Scalar> > frozen_x =
537 state_solution_history->getCurrentState()->getX();
538 RCP<const VectorBase<Scalar> > frozen_x_dot =
539 state_solution_history->getCurrentState()->getXDot();
540 RCP<const VectorBase<Scalar> > frozen_x_dot_dot =
541 state_solution_history->getCurrentState()->getXDotDot();
542 RCP<const SolutionHistory<Scalar> > sens_solution_history =
543 sens_integrator_->getSolutionHistory();
544 num_states = sens_solution_history->getNumStates();
545 for (int i=0; i<num_states; ++i) {
546 RCP<const SolutionState<Scalar> > state = (*sens_solution_history)[i];
547
548 // X
549 RCP< MultiVectorBase<Scalar> > x_mv =
550 createMembers(x_space, num_param+1);
551 RCP<const MultiVectorBase<Scalar> > dxdp =
552 rcp_dynamic_cast<const DMVPV>(state->getX())->getMultiVector();
553 assign(x_mv->col(0).ptr(), *(frozen_x));
554 assign(x_mv->subView(rng).ptr(), *dxdp);
555 RCP<VectorBase<Scalar> > x = multiVectorProductVector(prod_space, x_mv);
556
557 // X-Dot
558 RCP<VectorBase<Scalar> > x_dot;
559 if (state->getXDot() != Teuchos::null) {
560 RCP< MultiVectorBase<Scalar> > x_dot_mv =
561 createMembers(x_space, num_param+1);
562 RCP<const MultiVectorBase<Scalar> > dxdotdp =
563 rcp_dynamic_cast<const DMVPV>(state->getXDot())->getMultiVector();
564 assign(x_dot_mv->col(0).ptr(), *(frozen_x_dot));
565 assign(x_dot_mv->subView(rng).ptr(), *dxdotdp);
566 x_dot = multiVectorProductVector(prod_space, x_dot_mv);
567 }
568
569 // X-Dot-Dot
570 RCP<VectorBase<Scalar> > x_dot_dot;
571 if (state->getXDotDot() != Teuchos::null) {
572 RCP< MultiVectorBase<Scalar> > x_dot_dot_mv =
573 createMembers(x_space, num_param+1);
574 RCP<const MultiVectorBase<Scalar> > dxdotdotdp =
575 rcp_dynamic_cast<const DMVPV>(state->getXDotDot())->getMultiVector();
576 assign(x_dot_dot_mv->col(0).ptr(), *(frozen_x_dot_dot));
577 assign(x_dot_dot_mv->subView(rng).ptr(), *dxdotdotdp);
578 x_dot_dot = multiVectorProductVector(prod_space, x_dot_dot_mv);
579 }
580
581 RCP<SolutionState<Scalar> > prod_state = state->clone();
582 prod_state->setX(x);
583 prod_state->setXDot(x_dot);
584 prod_state->setXDotDot(x_dot_dot);
585 solutionHistory_->addState(prod_state, false);
586 }
587}
588
590template<class Scalar>
591Teuchos::RCP<Tempus::IntegratorPseudoTransientForwardSensitivity<Scalar> >
593 Teuchos::RCP<Teuchos::ParameterList> pList,
594 const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& model,
595 const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& sens_residual_model,
596 const Teuchos::RCP<Thyra::ModelEvaluator<Scalar> >& sens_solve_model)
597{
598
599 auto fwd_integrator = createIntegratorBasic<Scalar>(pList, model);
600 Teuchos::RCP<SensitivityModelEvaluatorBase<Scalar> > sens_model;
601 Teuchos::RCP<IntegratorBasic<Scalar> > sens_integrator;
602
603 {
604 Teuchos::RCP<Teuchos::ParameterList> pl = Teuchos::rcp(new Teuchos::ParameterList);
605 Teuchos::RCP<const Teuchos::ParameterList> integrator_pl = fwd_integrator->getValidParameters();
606 Teuchos::RCP<const Teuchos::ParameterList> sensitivity_pl =
608 pl->setParameters(*integrator_pl);
609 pl->sublist("Sensitivities").setParameters(*sensitivity_pl);
610 pl->sublist("Sensitivities").set("Reuse State Linear Solver", false);
611 pl->sublist("Sensitivities").set("Force W Update", false);
612 pl->sublist("Sensitivities").set("Cache Matrices", false);
613 pList->setParametersNotAlreadySet(*pl);
614 }
615
616 bool reuse_solver = pList->sublist("Sensitivities").get("Reuse State Linear Solver", false);
617 bool force_W_update = pList->sublist("Sensitivities").get("Force W Update", false);
618 bool cache_matrices = pList->sublist("Sensitivities").get("Cache Matrices", false);
619
620 {
621 Teuchos::RCP<Teuchos::ParameterList> pl = Teuchos::parameterList();
622 if (pList!= Teuchos::null)
623 {
624 *pl = pList->sublist("Sensitivities");
625 }
626 pl->remove("Reuse State Linear Solver");
627 pl->remove("Force W Update");
628 pl->remove("Cache Matrices");
630 model, sens_residual_model, sens_solve_model, cache_matrices, pl);
631 sens_integrator = createIntegratorBasic<Scalar>(pList, sens_model);
632 }
633
634 Teuchos::RCP<Tempus::IntegratorPseudoTransientForwardSensitivity<Scalar>> integrator =
635 Teuchos::rcp(new Tempus::IntegratorPseudoTransientForwardSensitivity<Scalar>(model, sens_model, fwd_integrator, sens_integrator, reuse_solver, force_W_update));
636
637 return(integrator);
638}
639
640
642template<class Scalar>
643Teuchos::RCP<Tempus::IntegratorPseudoTransientForwardSensitivity<Scalar> >
645{
646 Teuchos::RCP<Tempus::IntegratorPseudoTransientForwardSensitivity<Scalar> > integrator =
648 return(integrator);
649}
650
651} // namespace Tempus
652#endif // Tempus_IntegratorPseudoTransientForwardSensitivity_impl_hpp
IntegratorObserver class for time integrators.
Time integrator suitable for pseudotransient forward sensitivity analysis.
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getXDotDot() const
Get current the second time derivative of the solution, xdotdot.
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getG() const
Return response function g.
virtual Teuchos::RCP< TimeStepControl< Scalar > > getNonConstTimeStepControl() override
SensitivityStepMode getStepMode() const
What mode the current time integration step is in.
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 Thyra::VectorBase< Scalar > > getX() const
Get current the solution, x.
virtual Teuchos::RCP< SolutionHistory< Scalar > > getNonConstSolutionHistory() override
Get the SolutionHistory.
virtual Teuchos::RCP< IntegratorObserver< Scalar > > getObserver()
Get the Observer.
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getDgDp() const
Return forward sensitivity stored in Jacobian format.
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > getXDot() const
Get current the time derivative of the solution, xdot.
virtual Teuchos::RCP< const TimeStepControl< Scalar > > getTimeStepControl() const override
Get the TimeStepControl.
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getDXDotDp() const
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getDXDotDotDp() const
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const override
virtual Teuchos::RCP< Stepper< Scalar > > getStepper() const override
Get the Stepper.
virtual Teuchos::RCP< const SolutionHistory< Scalar > > getSolutionHistory() const override
Get the SolutionHistory.
virtual void setObserver(Teuchos::RCP< IntegratorObserver< Scalar > > obs=Teuchos::null)
Set the Observer.
virtual bool advanceTime()
Advance the solution to timeMax, and return true if successful.
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > getDxDp() const
A ModelEvaluator decorator for sensitivity analysis.
Status
Status for the Integrator, the Stepper and the SolutionState.
Teuchos::RCP< SensitivityModelEvaluatorBase< Scalar > > wrapStaggeredFSAModelEvaluator(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 bool is_pseudotransient, const Teuchos::RCP< const Teuchos::ParameterList > &pList=Teuchos::null)
Teuchos::RCP< Tempus::IntegratorPseudoTransientForwardSensitivity< Scalar > > createIntegratorPseudoTransientForwardSensitivity()
Nonmember constructor.