Rythmos - Transient Integration for Differential Equations  Version of the Day
 All Classes Functions Variables Typedefs Pages
List of all members
Rythmos::StateAndForwardSensitivityModelEvaluator< Scalar > Class Template Reference

Combined State and Forward sensitivity transient ModelEvaluator subclass. More...

#include <Rythmos_StateAndForwardSensitivityModelEvaluator.hpp>

Inherits StateFuncModelEvaluatorBase< Scalar >.

Constructors/Intializers/Accessors

 StateAndForwardSensitivityModelEvaluator ()
 
void initializeStructure (const Teuchos::RCP< const ForwardSensitivityModelEvaluatorBase< Scalar > > &sensModel)
 Set up the structure of the state and sensitivity equations. More...
 
Teuchos::RCP< const
Thyra::DefaultProductVector
< Scalar > > 
create_x_bar_vec (const Teuchos::RCP< const Thyra::VectorBase< Scalar > > &x_vec, const Teuchos::RCP< const Thyra::VectorBase< Scalar > > &s_bar_vec) const
 Create a wrapped product vector of the form x_bar = [ x; s_bar ]. More...
 

Public functions overridden from ModelEvaulator.

int Np () const
 Returns 0 . More...
 
Teuchos::RCP< const
Thyra::VectorSpaceBase< Scalar > > 
get_p_space (int l) const
 
Teuchos::RCP< const
Teuchos::Array< std::string > > 
get_p_names (int l) const
 
Teuchos::RCP< const
Thyra::VectorSpaceBase< Scalar > > 
get_x_space () const
 
Teuchos::RCP< const
Thyra::VectorSpaceBase< Scalar > > 
get_f_space () const
 
Thyra::ModelEvaluatorBase::InArgs
< Scalar > 
getNominalValues () const
 
Teuchos::RCP
< Thyra::LinearOpWithSolveBase
< Scalar > > 
create_W () const
 
Thyra::ModelEvaluatorBase::InArgs
< Scalar > 
createInArgs () const
 

Detailed Description

template<class Scalar>
class Rythmos::StateAndForwardSensitivityModelEvaluator< Scalar >

Combined State and Forward sensitivity transient ModelEvaluator subclass.

This class provides an implemenation of a combined state and forward sensitivity model evaluator for a DAE.

The form of the parameterized state equation is:

  f(x_dot(t),x(t),p) = 0, over t = [t0,tf]

  x(t0) = x_init(p)

The forward sensitivity equations, as written in multi-vector form, are:

  d(f)/d(x_dot)*S_dot + d(f)/d(x)*S + d(f)/d(p) = 0, over t = [t0,tf]

  S(t0) = d(x_init)/d(p)

where S is a multi-vector with np columns where each column S(:,j) = d(x)/d(p_j) is the sensitivity of x(t) with respect to the p_j parameter.

This model evaluator class represents the full state plus forward sensitivity system given as:

  f_bar(x_bar_dot(t),x_bar(t)) = 0, over t = [t0,tf]

  x_bar(t0) = x_bar_init

where

  x_bar = [ x; s_bar ] 

  s_bar = [ S(:,0); S(:,0); ...; S(:,np-1) ]

and f_bar(...) is the obvious concatenated state and sensitivity systems expresses in product vector form.

The vector x_bar is represented as a Thyra::ProductVectorBase object with two vector blocks x and s_bar. The flattened out sensitivity vector s_bar is then represented as a specialized product vector of type Thyra::DefaultMultiVectorProductVector.

If x_bar is an RCP<Thyra::VectorBase<Scalar> > object, then x can be access as

RCP<Thyra::VectorBase<Scalar> >
x = Thyra::nonconstProductVectorBase<Scalar>(x_bar)->getVectorBlock(0);

If x_bar is an RCP<const Thyra::VectorBase<Scalar> > object, then x can be access as

RCP<const Thyra::VectorBase<Scalar> >
x = Thyra::productVectorBase<Scalar>(x_bar)->getVectorBlock(0);

Likewise, s_bar can be access as

RCP<Thyra::VectorBase<Scalar> >
s_bar = Thyra::nonconstProductVectorBase<Scalar>(x_bar)->getVectorBlock(1);

when non-const and when const as:

RCP<const Thyra::VectorBase<Scalar> >
s_bar = Thyra::productVectorBase<Scalar>(x_bar)->getVectorBlock(1);

Given the flattened out vector form s_bar, one can get the underlying mulit-vector S as:

typedef Thyra::DefaultMultiVectorProductVector<Scalar> DMVPV;
RCP<Thyra::MultiVectorBase<Scalar> >
S = rcp_dynamic_cast<DMVPV>(s_bar)->getNonconstMultiVector();

for a nonconst vector/multi-vector and for a const vector/multi-vector as:

typedef Thyra::DefaultMultiVectorProductVector<Scalar> DMVPV;
RCP<const Thyra::MultiVectorBase<Scalar> >
S = rcp_dynamic_cast<const DMVPV>(s_bar)->getMultiVector();

ToDo: Replace the above documentation with the helper functions that will do all of this!

Currently, this class does not implement the full ModelEvaluator interface and it really just provides the spaces for x_bar and f_bar and the InArgs and OutArgs creation functions to allow for the full specification of the ForwardSensitivityStepper class. This is especially important in order to correctly set the full initial condition for the state and the forward sensitivities. Later this class can be completely finished in which case it would necessarly implement the evaluations and the linear solves which would automatically support the simultaneous corrector method.

ToDo: Finish documentation!

Definition at line 184 of file Rythmos_StateAndForwardSensitivityModelEvaluator.hpp.

Constructor & Destructor Documentation

Member Function Documentation

template<class Scalar >
void Rythmos::StateAndForwardSensitivityModelEvaluator< Scalar >::initializeStructure ( const Teuchos::RCP< const ForwardSensitivityModelEvaluatorBase< Scalar > > &  sensModel)

Set up the structure of the state and sensitivity equations.

Parameters
sensModel[in,persisting] The structure-initialized forward sensitivity model. From this object, the state model and other information will be extracted.

Definition at line 286 of file Rythmos_StateAndForwardSensitivityModelEvaluator.hpp.

template<class Scalar >
Teuchos::RCP< const Thyra::DefaultProductVector< Scalar > > Rythmos::StateAndForwardSensitivityModelEvaluator< Scalar >::create_x_bar_vec ( const Teuchos::RCP< const Thyra::VectorBase< Scalar > > &  x_vec,
const Teuchos::RCP< const Thyra::VectorBase< Scalar > > &  s_bar_vec 
) const

Create a wrapped product vector of the form x_bar = [ x; s_bar ].

Note: This does not copy any vector data, it only creates the wrapped product vector.

Definition at line 319 of file Rythmos_StateAndForwardSensitivityModelEvaluator.hpp.

template<class Scalar >
int Rythmos::StateAndForwardSensitivityModelEvaluator< Scalar >::Np ( ) const

Returns 0 .

Definition at line 340 of file Rythmos_StateAndForwardSensitivityModelEvaluator.hpp.

template<class Scalar >
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > Rythmos::StateAndForwardSensitivityModelEvaluator< Scalar >::get_p_space ( int  l) const
template<class Scalar >
Teuchos::RCP< const Teuchos::Array< std::string > > Rythmos::StateAndForwardSensitivityModelEvaluator< Scalar >::get_p_names ( int  l) const
template<class Scalar >
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > Rythmos::StateAndForwardSensitivityModelEvaluator< Scalar >::get_x_space ( ) const
template<class Scalar >
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > Rythmos::StateAndForwardSensitivityModelEvaluator< Scalar >::get_f_space ( ) const
template<class Scalar >
Thyra::ModelEvaluatorBase::InArgs< Scalar > Rythmos::StateAndForwardSensitivityModelEvaluator< Scalar >::getNominalValues ( ) const
template<class Scalar >
Teuchos::RCP< Thyra::LinearOpWithSolveBase< Scalar > > Rythmos::StateAndForwardSensitivityModelEvaluator< Scalar >::create_W ( ) const
template<class Scalar >
Thyra::ModelEvaluatorBase::InArgs< Scalar > Rythmos::StateAndForwardSensitivityModelEvaluator< Scalar >::createInArgs ( ) const

The documentation for this class was generated from the following file: