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... | |
SensitivityStepMode | getStepMode () const |
What mode the current time integration step is in. More... | |
![]() | |
![]() | |
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) |
![]() | |
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) |
![]() | |
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) |
![]() | |
void | initializeVerboseObject (const EVerbosityLevel verbLevel=VERB_DEFAULT, const RCP< FancyOStream > &oStream=Teuchos::null) |
![]() | |
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 void | setDefaultVerbLevel (const EVerbosityLevel defaultVerbLevel) |
static EVerbosityLevel | getDefaultVerbLevel () |
![]() | |
static void | setDefaultOStream (const RCP< FancyOStream > &defaultOStream) |
static RCP< FancyOStream > | getDefaultOStream () |
![]() | |
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 25 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
Tempus::IntegratorAdjointSensitivity< Scalar >::IntegratorAdjointSensitivity | ( | ) |
Constructor that requires a subsequent setParameterList, setStepper, and initialize calls.
Definition at line 65 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
inlinevirtual |
Destructor.
Definition at line 132 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
virtual |
Advance the solution to timeMax, and return true if successful.
Definition at line 75 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Advance the solution to timeFinal, and return true if successful.
Implements Tempus::Integrator< Scalar >.
Definition at line 85 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Get current time.
Implements Tempus::Integrator< Scalar >.
Definition at line 313 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Get current index.
Implements Tempus::Integrator< Scalar >.
Definition at line 321 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Get Status.
Implements Tempus::Integrator< Scalar >.
Definition at line 329 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Set Status.
Implements Tempus::Integrator< Scalar >.
Definition at line 341 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Get the Stepper.
Implements Tempus::Integrator< Scalar >.
Definition at line 349 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Get the SolutionHistory.
Implements Tempus::Integrator< Scalar >.
Definition at line 357 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
Teuchos::RCP< const SolutionHistory< Scalar > > Tempus::IntegratorAdjointSensitivity< Scalar >::getStateSolutionHistory | ( | ) | const |
Definition at line 365 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
Teuchos::RCP< const SolutionHistory< Scalar > > Tempus::IntegratorAdjointSensitivity< Scalar >::getSensSolutionHistory | ( | ) | const |
Definition at line 373 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Get the SolutionHistory.
Implements Tempus::Integrator< Scalar >.
Definition at line 381 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Get the TimeStepControl.
Implements Tempus::Integrator< Scalar >.
Definition at line 389 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Implements Tempus::Integrator< Scalar >.
Definition at line 397 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
Teuchos::RCP< TimeStepControl< Scalar > > Tempus::IntegratorAdjointSensitivity< Scalar >::getStateNonConstTimeStepControl | ( | ) |
Definition at line 405 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
Teuchos::RCP< TimeStepControl< Scalar > > Tempus::IntegratorAdjointSensitivity< Scalar >::getSensNonConstTimeStepControl | ( | ) |
Definition at line 413 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
inlineoverridevirtual |
Returns the IntegratorTimer_ for this Integrator.
Implements Tempus::Integrator< Scalar >.
Definition at line 163 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
inlineoverridevirtual |
Implements Tempus::Integrator< Scalar >.
Definition at line 165 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
virtual |
Set the initial state from Thyra::VectorBase(s)
Definition at line 420 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Get the Observer.
Definition at line 435 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Set the Observer.
Definition at line 442 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Initializes the Integrator after set* function calls.
Definition at line 452 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Get the current solution, x.
Definition at line 461 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Get the current time derivative of the solution, xdot.
Definition at line 469 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Get the current second time derivative of the solution, xdotdot.
Definition at line 477 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Get the current adjoint solution, y.
Definition at line 485 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Get the current time derivative of the adjoint solution, ydot.
Definition at line 501 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Get the current second time derivative of the adjoint solution, ydotdot.
Definition at line 517 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
virtual |
Return adjoint sensitivity stored in gradient format.
Definition at line 533 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Reimplemented from Teuchos::Describable.
Definition at line 541 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
overridevirtual |
Reimplemented from Teuchos::Describable.
Definition at line 550 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
SensitivityStepMode Tempus::IntegratorAdjointSensitivity< Scalar >::getStepMode | ( | ) | const |
What mode the current time integration step is in.
Definition at line 566 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
protected |
Definition at line 574 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
protected |
Definition at line 604 of file Tempus_IntegratorAdjointSensitivity_impl.hpp.
|
protected |
Definition at line 228 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 229 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 230 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 231 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 232 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 233 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 234 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 235 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 236 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 237 of file Tempus_IntegratorAdjointSensitivity_decl.hpp.
|
protected |
Definition at line 238 of file Tempus_IntegratorAdjointSensitivity_decl.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.