Tempus
Version of the Day
Time Integration
|
Time integrator suitable for adjoint sensitivity analysis. More...
#include <Tempus_IntegratorAdjointSensitivity_decl.hpp>
Public Member Functions | |
IntegratorAdjointSensitivity (const Teuchos::RCP< Thyra::ModelEvaluator< Scalar >> &model, const Teuchos::RCP< IntegratorBasic< Scalar >> &state_integrator, const Teuchos::RCP< Thyra::ModelEvaluator< Scalar >> &adjoint_model, const Teuchos::RCP< AdjointAuxSensitivityModelEvaluator< Scalar >> &adjoint_aux_model, const Teuchos::RCP< IntegratorBasic< Scalar >> &adjoint_integrator, const Teuchos::RCP< SolutionHistory< Scalar >> &solution_history, const int p_index, const int g_index, const bool g_depends_on_p, const bool f_depends_on_p, const bool ic_depends_on_p, const bool mass_matrix_is_identity) | |
Full Constructor with model, and will be fully initialized. More... | |
IntegratorAdjointSensitivity () | |
Constructor that requires a subsequent setParameterList, setStepper, and initialize calls. More... | |
virtual | ~IntegratorAdjointSensitivity () |
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 < IntegratorObserver< Scalar > > | getObserver () |
Get the Observer. More... | |
virtual void | setObserver (Teuchos::RCP< IntegratorObserver< Scalar >> obs=Teuchos::null) |
Set the Observer. More... | |
virtual void | initialize () |
Initializes the Integrator after set* function calls. More... | |
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > | getX () const |
Get the current solution, x. More... | |
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > | getXDot () const |
Get the current time derivative of the solution, xdot. More... | |
virtual Teuchos::RCP< const Thyra::VectorBase< Scalar > > | getXDotDot () const |
Get the current second time derivative of the solution, xdotdot. More... | |
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > | getY () const |
Get the current adjoint solution, y. More... | |
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > | getYDot () const |
Get the current time derivative of the adjoint solution, ydot. More... | |
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > | getYDotDot () const |
Get the current second time derivative of the adjoint solution, ydotdot. More... | |
virtual Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > | getDgDp () const |
Return adjoint sensitivity stored in gradient format. More... | |
Teuchos::RCP < Thyra::ModelEvaluator < Scalar > > | getAdjointModel () const |
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 | |
Teuchos::RCP < AdjointAuxSensitivityModelEvaluator < Scalar > > | createAdjointModel (const Teuchos::RCP< Thyra::ModelEvaluator< Scalar >> &model, const Teuchos::RCP< Thyra::ModelEvaluator< Scalar >> &adjoint_model, const Teuchos::RCP< Teuchos::ParameterList > &inputPL) |
void | buildSolutionHistory (const Teuchos::RCP< const SolutionHistory< Scalar >> &state_solution_history, const Teuchos::RCP< const SolutionHistory< Scalar >> &adjoint_solution_history) |
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< IntegratorBasic < Scalar > > | state_integrator_ |
Teuchos::RCP < Thyra::ModelEvaluator < Scalar > > | adjoint_model_ |
Teuchos::RCP < AdjointAuxSensitivityModelEvaluator < Scalar > > | adjoint_aux_model_ |
Teuchos::RCP< IntegratorBasic < Scalar > > | adjoint_integrator_ |
Teuchos::RCP< SolutionHistory < Scalar > > | solutionHistory_ |
int | p_index_ |
int | g_index_ |
bool | g_depends_on_p_ |
bool | f_depends_on_p_ |
bool | ic_depends_on_p_ |
bool | mass_matrix_is_identity_ |
Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > | dxdp_init_ |
Teuchos::RCP < Thyra::MultiVectorBase < Scalar > > | dgdp_ |
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... | |
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 < 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 adjoint sensitivity analysis.
This integrator implements transient adjoint sensitivities. Given a model evaluator encapsulating the equations f(x_dot,x,p) = 0, and a response function g(x,p), these equations are integrated forward in time to some final time T. Then the adjoint equations F(y) = d/dt( df/dx_dot^T*y ) - df/dx^T*y = 0 are integrated backward in time starting at t = T where y is the adjoint variable. Nominally y is a multi-vector belonging to the vector space of f and the number of columns given by the number of entries in the response g. The initial conditions for y at t = T are given by y(T) = df/dx_dot(x_dot(T),x(T),p)^{-T} * dg/dx(x(T),p)^T. Then the final sensitivity of g is dg/dp(T) - int_{0}^T(df/dp^T*y)dt + dx/dp(0)^T*df/dx_dot(0)^T*y(0).
This integrator supports both implicit and explicit steppers, and provides a method to compute dg/dp as described above after the reverse integration. The solution history contains solutions stored as product vectors (x,y), with y further stored as a product vector for each component of g.
Because of limitations on the steppers, the implementation currently assumes df/dxdot is a constant matrix.
To extract the final solution x(T) and sensitivity dg/dp(T) one should use the getX() and getDgDp() methods, which return these quantities directly. One can also extract this data for all times from the solution history, however the data is stored in Thyra product vectors which requires knowledge of the internal implementation.
Definition at line 52 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
Tempus::IntegratorAdjointSensitivity< Scalar >::IntegratorAdjointSensitivity | ( | const Teuchos::RCP< Thyra::ModelEvaluator< Scalar >> & | model, |
const Teuchos::RCP< IntegratorBasic< Scalar >> & | state_integrator, | ||
const Teuchos::RCP< Thyra::ModelEvaluator< Scalar >> & | adjoint_model, | ||
const Teuchos::RCP< AdjointAuxSensitivityModelEvaluator< Scalar >> & | adjoint_aux_model, | ||
const Teuchos::RCP< IntegratorBasic< Scalar >> & | adjoint_integrator, | ||
const Teuchos::RCP< SolutionHistory< Scalar >> & | solution_history, | ||
const int | p_index, | ||
const int | g_index, | ||
const bool | g_depends_on_p, | ||
const bool | f_depends_on_p, | ||
const bool | ic_depends_on_p, | ||
const bool | mass_matrix_is_identity | ||
) |
Full Constructor with model, and will be fully initialized.
[in] | model | The forward physics ModelEvaluator |
[in] | state_integrator | Forward state Integrator for the forward problem |
[in] | adjoint_model | ModelEvaluator for the adjoint physics/problem |
[in] | adjoint_aux_model | ModelEvaluator for the auxiliary adjoint physics/problem |
[in] | adjoint_integrator | Time integrator for the adjoint problem |
[in] | solution_history | The forward state solution history |
[in] | p_index | Sensitivity parameter index |
[in] | g_index | Response function index |
[in] | g_depends_on_p | Does response depends on parameters? |
[in] | f_depends_on_p | Does residual depends on parameters? |
[in] | ic_depends_on_p | Does the initial condition depends on parameters? |
adjoint_model | ModelEvaluator for the adjoint problem. Optional. Default value is null. | |
[in] | mass_matrix_is_identity | Is the mass matrix an identity matrix? |
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:
To support use-cases with explicitly computed adjoint operators, the constructor takes an additional model evaluator for computing the adjoint W/W_op. It is assumed the operator returned by this model evaluator is the adjoint, and so will not be transposed. It is also assumed this model evaluator accepts the same inArgs as the forward model, however it only requires supporting the adjoint W/W_op outArgs.
Version of the constructor taking a single model evaluator.
This version takes a single model evaluator for the case when the adjoint is implicitly determined from the forward operator by the (conjugate) transpose.
Definition at line 26 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
Tempus::IntegratorAdjointSensitivity< Scalar >::IntegratorAdjointSensitivity | ( | ) |
Constructor that requires a subsequent setParameterList, setStepper, and initialize calls.
Definition at line 62 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
inlinevirtual |
Destructor.
Definition at line 130 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
virtual |
Advance the solution to timeMax, and return true if successful.
Definition at line 70 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Advance the solution to timeFinal, and return true if successful.
Implements Tempus::Integrator< Scalar >.
Definition at line 77 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Get current time.
Implements Tempus::Integrator< Scalar >.
Definition at line 308 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Get current index.
Implements Tempus::Integrator< Scalar >.
Definition at line 314 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Get Status.
Implements Tempus::Integrator< Scalar >.
Definition at line 320 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Set Status.
Implements Tempus::Integrator< Scalar >.
Definition at line 330 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Get the Stepper.
Implements Tempus::Integrator< Scalar >.
Definition at line 337 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Get the SolutionHistory.
Implements Tempus::Integrator< Scalar >.
Definition at line 345 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
Teuchos::RCP< const SolutionHistory< Scalar > > Tempus::IntegratorAdjointSensitivity< Scalar >::getStateSolutionHistory | ( | ) | const |
Definition at line 352 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
Teuchos::RCP< const SolutionHistory< Scalar > > Tempus::IntegratorAdjointSensitivity< Scalar >::getSensSolutionHistory | ( | ) | const |
Definition at line 359 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Get the SolutionHistory.
Implements Tempus::Integrator< Scalar >.
Definition at line 366 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Get the TimeStepControl.
Implements Tempus::Integrator< Scalar >.
Definition at line 373 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Implements Tempus::Integrator< Scalar >.
Definition at line 380 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
Teuchos::RCP< TimeStepControl< Scalar > > Tempus::IntegratorAdjointSensitivity< Scalar >::getStateNonConstTimeStepControl | ( | ) |
Definition at line 387 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
Teuchos::RCP< TimeStepControl< Scalar > > Tempus::IntegratorAdjointSensitivity< Scalar >::getSensNonConstTimeStepControl | ( | ) |
Definition at line 394 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
inlineoverridevirtual |
Returns the IntegratorTimer_ for this Integrator.
Implements Tempus::Integrator< Scalar >.
Definition at line 165 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
inlineoverridevirtual |
Implements Tempus::Integrator< Scalar >.
Definition at line 169 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
virtual |
Set the initial state from Thyra::VectorBase(s)
Definition at line 400 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Get the Observer.
Definition at line 414 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Set the Observer.
Definition at line 420 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Initializes the Integrator after set* function calls.
Definition at line 430 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Get the current solution, x.
Definition at line 438 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Get the current time derivative of the solution, xdot.
Definition at line 445 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Get the current second time derivative of the solution, xdotdot.
Definition at line 452 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Get the current adjoint solution, y.
Definition at line 459 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Get the current time derivative of the adjoint solution, ydot.
Definition at line 472 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Get the current second time derivative of the adjoint solution, ydotdot.
Definition at line 486 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Return adjoint sensitivity stored in gradient format.
Definition at line 500 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
inline |
Definition at line 212 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
overridevirtual |
Reimplemented from Teuchos::Describable.
Definition at line 506 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Reimplemented from Teuchos::Describable.
Definition at line 513 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
SensitivityStepMode Tempus::IntegratorAdjointSensitivity< Scalar >::getStepMode | ( | ) | const |
What mode the current time integration step is in.
Definition at line 526 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
protected |
Definition at line 533 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
protected |
Definition at line 560 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
protected |
Definition at line 239 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 240 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 241 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 242 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 243 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 244 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 245 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 246 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 247 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 248 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 249 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 250 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 251 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 252 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 253 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.