Tempus
Version of the Day
Time Integration
|
Time integrator suitable for pseudotransient forward sensitivity analysis. More...
#include <Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp>
Public Member Functions | |
IntegratorPseudoTransientForwardSensitivity (const Teuchos::RCP< Thyra::ModelEvaluator< Scalar >> &model, const Teuchos::RCP< SensitivityModelEvaluatorBase< Scalar >> &sens_model, const Teuchos::RCP< IntegratorBasic< Scalar >> &fwd_integrator, const Teuchos::RCP< IntegratorBasic< Scalar >> &sens_integrator, const bool reuse_solver, const bool force_W_update) | |
Constructor with ParameterList and model, and will be fully initialized. More... | |
IntegratorPseudoTransientForwardSensitivity () | |
Destructor. More... | |
virtual | ~IntegratorPseudoTransientForwardSensitivity () |
Destructor. More... | |
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) More... | |
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > | getX () const |
Get current the solution, x. More... | |
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > | getDxDp () const |
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > | getXDot () const |
Get current the time derivative of the solution, xdot. More... | |
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > | getDXDotDp () const |
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > | getXDotDot () const |
Get current the second time derivative of the solution, xdotdot. More... | |
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > | getDXDotDotDp () const |
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > | getG () const |
Return response function g. More... | |
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > | getDgDp () const |
Return forward sensitivity stored in Jacobian format. More... | |
SensitivityStepMode | getStepMode () const |
What mode the current time integration step is in. More... | |
Public Member Functions inherited from Tempus::Integrator< Scalar > | |
Public Member Functions inherited from Teuchos::Describable | |
void | describe (std::ostream &out, const EVerbosityLevel verbLevel=verbLevel_default) const |
virtual | ~Describable () |
LabeledObject () | |
virtual | ~LabeledObject () |
virtual void | setObjectLabel (const std::string &objectLabel) |
virtual std::string | getObjectLabel () const |
DescribableStreamManipulatorState | describe (const Describable &describable, const EVerbosityLevel verbLevel=Describable::verbLevel_default) |
std::ostream & | operator<< (std::ostream &os, const DescribableStreamManipulatorState &d) |
Public Member Functions inherited from Teuchos::VerboseObject< Tempus::Integrator< Scalar > > | |
VerboseObject (const EVerbosityLevel verbLevel=VERB_DEFAULT, const RCP< FancyOStream > &oStream=Teuchos::null) | |
virtual const VerboseObject & | setVerbLevel (const EVerbosityLevel verbLevel) const |
virtual const VerboseObject & | setOverridingVerbLevel (const EVerbosityLevel verbLevel) const |
virtual EVerbosityLevel | getVerbLevel () const |
TEUCHOSPARAMETERLIST_LIB_DLL_EXPORT RCP< const ParameterList > | getValidVerboseObjectSublist () |
TEUCHOSPARAMETERLIST_LIB_DLL_EXPORT void | setupVerboseObjectSublist (ParameterList *paramList) |
TEUCHOSPARAMETERLIST_LIB_DLL_EXPORT void | readVerboseObjectSublist (ParameterList *paramList, RCP< FancyOStream > *oStream, EVerbosityLevel *verbLevel) |
void | readVerboseObjectSublist (ParameterList *paramList, VerboseObject< ObjectType > *verboseObject) |
Public Member Functions inherited from Teuchos::VerboseObjectBase | |
virtual | ~VerboseObjectBase () |
VerboseObjectBase (const RCP< FancyOStream > &oStream=Teuchos::null) | |
virtual const VerboseObjectBase & | setOStream (const RCP< FancyOStream > &oStream) const |
virtual const VerboseObjectBase & | setOverridingOStream (const RCP< FancyOStream > &oStream) const |
virtual VerboseObjectBase & | setLinePrefix (const std::string &linePrefix) |
virtual RCP< FancyOStream > | getOStream () const |
virtual RCP< FancyOStream > | getOverridingOStream () const |
virtual std::string | getLinePrefix () const |
virtual OSTab | getOSTab (const int tabs=1, const std::string &linePrefix="") const |
Protected Member Functions | |
void | buildSolutionHistory () |
Protected Member Functions inherited from Teuchos::VerboseObject< Tempus::Integrator< Scalar > > | |
void | initializeVerboseObject (const EVerbosityLevel verbLevel=VERB_DEFAULT, const RCP< FancyOStream > &oStream=Teuchos::null) |
Protected Member Functions inherited from Teuchos::VerboseObjectBase | |
void | initializeVerboseObjectBase (const RCP< FancyOStream > &oStream=Teuchos::null) |
virtual void | informUpdatedVerbosityState () const |
Protected Attributes | |
Teuchos::RCP < Thyra::ModelEvaluator < Scalar > > | model_ |
Teuchos::RCP < SensitivityModelEvaluatorBase < Scalar > > | sens_model_ |
Teuchos::RCP< IntegratorBasic < Scalar > > | state_integrator_ |
Teuchos::RCP< IntegratorBasic < Scalar > > | sens_integrator_ |
Teuchos::RCP< SolutionHistory < Scalar > > | solutionHistory_ |
bool | reuse_solver_ |
bool | force_W_update_ |
SensitivityStepMode | stepMode_ |
Basic integrator methods | |
virtual bool | advanceTime () |
Advance the solution to timeMax, and return true if successful. More... | |
virtual bool | advanceTime (const Scalar timeFinal) override |
Advance the solution to timeFinal, and return true if successful. More... | |
virtual Scalar | getTime () const override |
Get current time. More... | |
virtual int | getIndex () const override |
Get current index. More... | |
virtual Status | getStatus () const override |
Get Status. More... | |
virtual void | setStatus (const Status st) override |
Set Status. More... | |
virtual Teuchos::RCP< Stepper < Scalar > > | getStepper () const override |
Get the Stepper. More... | |
Teuchos::RCP< Stepper< Scalar > > | getStateStepper () const |
Teuchos::RCP< Stepper< Scalar > > | getSensStepper () const |
virtual Teuchos::RCP< const SolutionHistory< Scalar > > | getSolutionHistory () const override |
Get the SolutionHistory. More... | |
Teuchos::RCP< const SolutionHistory< Scalar > > | getStateSolutionHistory () const |
Teuchos::RCP< const SolutionHistory< Scalar > > | getSensSolutionHistory () const |
virtual Teuchos::RCP < SolutionHistory< Scalar > > | getNonConstSolutionHistory () override |
Get the SolutionHistory. More... | |
virtual Teuchos::RCP< const TimeStepControl< Scalar > > | getTimeStepControl () const override |
Get the TimeStepControl. More... | |
virtual Teuchos::RCP < TimeStepControl< Scalar > > | getNonConstTimeStepControl () override |
Teuchos::RCP< TimeStepControl < Scalar > > | getStateNonConstTimeStepControl () |
Teuchos::RCP< TimeStepControl < Scalar > > | getSensNonConstTimeStepControl () |
virtual Teuchos::RCP < IntegratorObserver< Scalar > > | getObserver () |
Get the Observer. More... | |
virtual void | setObserver (Teuchos::RCP< IntegratorObserver< Scalar >> obs=Teuchos::null) |
Set the Observer. More... | |
virtual Teuchos::RCP < Teuchos::Time > | getIntegratorTimer () const override |
Returns the IntegratorTimer_ for this Integrator. More... | |
virtual Teuchos::RCP < Teuchos::Time > | getStepperTimer () const override |
Overridden from Teuchos::Describable | |
std::string | description () const override |
void | describe (Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const override |
Additional Inherited Members | |
Static Public Member Functions inherited from Teuchos::VerboseObject< Tempus::Integrator< Scalar > > | |
static void | setDefaultVerbLevel (const EVerbosityLevel defaultVerbLevel) |
static EVerbosityLevel | getDefaultVerbLevel () |
Static Public Member Functions inherited from Teuchos::VerboseObjectBase | |
static void | setDefaultOStream (const RCP< FancyOStream > &defaultOStream) |
static RCP< FancyOStream > | getDefaultOStream () |
Static Public Attributes inherited from Teuchos::Describable | |
static const EVerbosityLevel | verbLevel_default |
Time integrator suitable for pseudotransient forward sensitivity analysis.
For some problems, time integrators are used to compute steady-state solutions (also known as pseudo-transient solvers). When computing sensitivities, it is not necessary in these cases to propagate sensitivities all the way through the forward time integration. Instead the steady-state is first computed as usual, and then the sensitivities are computed using a similar pseudo-transient time integration applied to the sensitivity equations with the state frozen to the computed steady-state. This integrator specializes the transient sensitivity methods implemented by Tempus::IntegratorForwardSensitivity to this case.
Consider an implicit ODE f(x_dot,x,p) = 0 with a stable steady-state solution x = x^s, x_dot = 0 where f(0,x^s,p) = 0 and all of the eigenvalues of df/dx(0,x^s,p) are in the right half-plane (for an explicit ODE, the eigenvalues must be in the left half-plane). In the pseudo-transient method a time-integrator is applied to f(x_dot,x,p) = 0 until x_dot is sufficiently small. Now consider the forward sensitivity equations: df/dx_dot*z_dot + df/dx*z + df/dp = 0 where z = dx/dp. For pseudo-transient forward sensitivities, the above is integrated from z(0) = 0 until z_dot is sufficiently small, in which case z^s = -(df/dx)^{-1}*(df/dp). Then the final sensitivity of g is dg/dp + dg/dx*z^s. One can see that z^s is the only steady-state solution of the sensitivity equations, since df/dx and df/dp are constant, and must be linearly stable since it has the same Jacobian matrix as the forward equations.
One should use the getX() and getDxDp() methods for extracting the final sultion and its parameter sensitivity as a multi-vector. This data can also be extracted from the solution history, but is stored as a Thyra product vector which requires knowledge of the internal implementation.
Definition at line 58 of file Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp.
Tempus::IntegratorPseudoTransientForwardSensitivity< Scalar >::IntegratorPseudoTransientForwardSensitivity | ( | const Teuchos::RCP< Thyra::ModelEvaluator< Scalar >> & | model, |
const Teuchos::RCP< SensitivityModelEvaluatorBase< Scalar >> & | sens_model, | ||
const Teuchos::RCP< IntegratorBasic< Scalar >> & | fwd_integrator, | ||
const Teuchos::RCP< IntegratorBasic< Scalar >> & | sens_integrator, | ||
const bool | reuse_solver, | ||
const bool | force_W_update | ||
) |
Constructor with ParameterList and model, and will be fully initialized.
In addition to all of the regular integrator options, the supplied parameter list supports the following options contained within a sublist "Sensitivities" from the top-level parameter list:
Definition at line 23 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
Tempus::IntegratorPseudoTransientForwardSensitivity< Scalar >::IntegratorPseudoTransientForwardSensitivity | ( | ) |
Destructor.
Constructor that requires a subsequent setStepper, and initialize calls.
Definition at line 41 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
inlinevirtual |
Destructor.
Definition at line 108 of file Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp.
|
virtual |
Advance the solution to timeMax, and return true if successful.
Definition at line 51 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
overridevirtual |
Advance the solution to timeFinal, and return true if successful.
Implements Tempus::Integrator< Scalar >.
Definition at line 80 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
overridevirtual |
Get current time.
Implements Tempus::Integrator< Scalar >.
Definition at line 126 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
overridevirtual |
Get current index.
Implements Tempus::Integrator< Scalar >.
Definition at line 132 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
overridevirtual |
Get Status.
Implements Tempus::Integrator< Scalar >.
Definition at line 138 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
overridevirtual |
Set Status.
Implements Tempus::Integrator< Scalar >.
Definition at line 148 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
overridevirtual |
Get the Stepper.
Implements Tempus::Integrator< Scalar >.
Definition at line 157 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
Teuchos::RCP< Stepper< Scalar > > Tempus::IntegratorPseudoTransientForwardSensitivity< Scalar >::getStateStepper | ( | ) | const |
Definition at line 164 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
Teuchos::RCP< Stepper< Scalar > > Tempus::IntegratorPseudoTransientForwardSensitivity< Scalar >::getSensStepper | ( | ) | const |
Definition at line 171 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
overridevirtual |
Get the SolutionHistory.
Implements Tempus::Integrator< Scalar >.
Definition at line 178 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
Teuchos::RCP< const SolutionHistory< Scalar > > Tempus::IntegratorPseudoTransientForwardSensitivity< Scalar >::getStateSolutionHistory | ( | ) | const |
Definition at line 185 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
Teuchos::RCP< const SolutionHistory< Scalar > > Tempus::IntegratorPseudoTransientForwardSensitivity< Scalar >::getSensSolutionHistory | ( | ) | const |
Definition at line 193 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
overridevirtual |
Get the SolutionHistory.
Implements Tempus::Integrator< Scalar >.
Definition at line 202 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
overridevirtual |
Get the TimeStepControl.
Implements Tempus::Integrator< Scalar >.
Definition at line 209 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
overridevirtual |
Implements Tempus::Integrator< Scalar >.
Definition at line 217 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
Teuchos::RCP< TimeStepControl< Scalar > > Tempus::IntegratorPseudoTransientForwardSensitivity< Scalar >::getStateNonConstTimeStepControl | ( | ) |
Definition at line 225 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
Teuchos::RCP< TimeStepControl< Scalar > > Tempus::IntegratorPseudoTransientForwardSensitivity< Scalar >::getSensNonConstTimeStepControl | ( | ) |
Definition at line 233 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
virtual |
Get the Observer.
Definition at line 240 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
virtual |
Set the Observer.
Definition at line 246 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
inlineoverridevirtual |
Returns the IntegratorTimer_ for this Integrator.
Implements Tempus::Integrator< Scalar >.
Definition at line 149 of file Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp.
|
inlineoverridevirtual |
Implements Tempus::Integrator< Scalar >.
Definition at line 153 of file Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp.
|
virtual |
Set the initial state from Thyra::VectorBase(s)
Definition at line 255 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
virtual |
Get current the solution, x.
Definition at line 303 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
virtual |
Definition at line 310 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
virtual |
Get current the time derivative of the solution, xdot.
Definition at line 322 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
virtual |
Definition at line 329 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
virtual |
Get current the second time derivative of the solution, xdotdot.
Definition at line 342 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
virtual |
Definition at line 349 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
virtual |
Return response function g.
Definition at line 362 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
virtual |
Return forward sensitivity stored in Jacobian format.
Definition at line 387 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
overridevirtual |
Reimplemented from Teuchos::Describable.
Definition at line 413 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
overridevirtual |
Reimplemented from Teuchos::Describable.
Definition at line 421 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
SensitivityStepMode Tempus::IntegratorPseudoTransientForwardSensitivity< Scalar >::getStepMode | ( | ) | const |
What mode the current time integration step is in.
Definition at line 435 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
protected |
Definition at line 441 of file Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp.
|
protected |
Definition at line 200 of file Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp.
|
protected |
Definition at line 201 of file Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp.
|
protected |
Definition at line 202 of file Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp.
|
protected |
Definition at line 203 of file Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp.
|
protected |
Definition at line 204 of file Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp.
|
protected |
Definition at line 205 of file Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp.
|
protected |
Definition at line 206 of file Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp.
|
protected |
Definition at line 207 of file Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp.