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

Concrete Thyra::ModelEvaluator subclass that turns a forward ODE/DAE with an observation into a parameterized evaluation of p -> g with forward sensitivities DgDp. More...

#include <Rythmos_ForwardSensitivityIntegratorAsModelEvaluator.hpp>

Inherits ResponseOnlyModelEvaluatorBase< Scalar >, and ParameterListAcceptor.

Constructors, Initialization, Misc.

 ForwardSensitivityIntegratorAsModelEvaluator ()
 
void initialize (const RCP< StepperBase< Scalar > > &stateStepper, const RCP< IntegratorBase< Scalar > > &stateIntegrator, const RCP< ForwardSensitivityStepper< Scalar > > &stateAndSensStepper, const RCP< IntegratorBase< Scalar > > &stateAndSensIntegrator, const Thyra::ModelEvaluatorBase::InArgs< Scalar > &stateAndSensInitCond, const Array< Scalar > &responseTimes, const Array< RCP< const Thyra::ModelEvaluator< Scalar > > > &responseFuncs, const Array< Thyra::ModelEvaluatorBase::InArgs< Scalar > > &responseFuncBasePoints, const ForwardSensitivityIntegratorAsModelEvaluatorTypes::EResponseType responseType)
 Initalized with all objects needed to define evaluation. More...
 
const Array< Scalar > & getResponseTimes () const
 

Overridden from ParameterListAcceptor

void setParameterList (RCP< Teuchos::ParameterList > const &paramList)
 
RCP< Teuchos::ParameterList > getNonconstParameterList ()
 
RCP< Teuchos::ParameterList > unsetParameterList ()
 
RCP< const Teuchos::ParameterList > getParameterList () const
 
RCP< const Teuchos::ParameterList > getValidParameters () const
 

Public functions overridden from ModelEvaulator.

int Np () const
 
int Ng () const
 
RCP< const
Thyra::VectorSpaceBase< Scalar > > 
get_p_space (int l) const
 
RCP< const
Thyra::VectorSpaceBase< Scalar > > 
get_g_space (int j) const
 
Thyra::ModelEvaluatorBase::InArgs
< Scalar > 
createInArgs () const
 

Detailed Description

template<class Scalar>
class Rythmos::ForwardSensitivityIntegratorAsModelEvaluator< Scalar >

Concrete Thyra::ModelEvaluator subclass that turns a forward ODE/DAE with an observation into a parameterized evaluation of p -> g with forward sensitivities DgDp.

The form of the parameterized state equation is:

  f(x_dot(t),x(t),{p_l},t) = 0, over t = [t0,tf]

  x(t0) = x_init(p)

The full response function takes one of two possible forms:

The individual reduced response functions are:

  g_hat_k(p) =  g_k(x_dot(t_k),x(t_k),p,t_k), k=0...NO-1

Above, the reduced response function g_hat_k(p) is the evaluation of the full response function at the state solution.

The first derivative of the individual reduced response functions is:

  d(g_hat_k)/d(p) =
     beta_k
     *
     (
       d(g_k)/d(x_dot) * d(x_dot)/d(p)
       + d(g_k)/d(x) * d(x)/d(p)
       + d(g_k)/d(p)
       ),

   for k=0...NO-1

Above, S_dot = d(x_dot)/d(p) and S = d(x_dot)/d(p) are computed by solving the forward sensitivity equations (see ForwardSensitivityStepper).

This class implements the reduced response function p -> d_hat(p) by integrating the state+sens equations and assembling the reduced response functions.

ToDo: Finish Documentation!

Definition at line 139 of file Rythmos_ForwardSensitivityIntegratorAsModelEvaluator.hpp.

Constructor & Destructor Documentation

Member Function Documentation

template<class Scalar >
void Rythmos::ForwardSensitivityIntegratorAsModelEvaluator< Scalar >::initialize ( const RCP< StepperBase< Scalar > > &  stateStepper,
const RCP< IntegratorBase< Scalar > > &  stateIntegrator,
const RCP< ForwardSensitivityStepper< Scalar > > &  stateAndSensStepper,
const RCP< IntegratorBase< Scalar > > &  stateAndSensIntegrator,
const Thyra::ModelEvaluatorBase::InArgs< Scalar > &  stateAndSensInitCond,
const Array< Scalar > &  responseTimes,
const Array< RCP< const Thyra::ModelEvaluator< Scalar > > > &  responseFuncs,
const Array< Thyra::ModelEvaluatorBase::InArgs< Scalar > > &  responseFuncBasePoints,
const ForwardSensitivityIntegratorAsModelEvaluatorTypes::EResponseType  responseType 
)

Initalized with all objects needed to define evaluation.

Parameters
stateStepper[inout,persisting] The state stepper. On input, this object must be fully initialized with everything it needs except for the initial condition.
stateIntegrator[inout,persisting] The integrator that will be used to integrate the state equations. This integrator need not be set up with a stepper as that will be done internally. It just needs to have its parameters set to know which algorithm to run. This integrator will be cloned to create the integrator for the state + sens equations if stateAndSensIntegrator==null.
stateAndSensStepper[inout,persisting] The state + sens stepper. On input, this object must be fully initialized with everything it needs except for the initial condition.
stateAndSensIntegrator[inout,persisting] The integrator that will be used to integrate the state + sens equations. This integrator need not be set up with a stepper as that will be done internally. It just needs to have its parameters set to know which algorithm to run. The stateIntegrator will be cloned for this purpose if stateAndSensIntegrator==null.
stateAndSensInitCond[in,persisting] The initial condition for the state + sens system. This object should also include any parameters or other data that is needed to define the initial condition. The underlying parameter subvector will be set and reset again and again in order to satisfy the evaluation p -> g(p) defined by this interface.
responseTimes[in] Array giving the sorted time points for each response function g_k(...). Note that the final time for the integrator will be taken to be finalTime = responseTimes.back().
responseFuncs[in,persisting] Array giving the response functions for each time point responseTimes[k]. It is assumed that the response function will be g(0) in this model and that the same parameter subvector used in the state function is used here also.
responseFuncBasePoints[in,persisting] Array giving the base point of the response functions.
responseType[in] Determines how the response functions are treated and how they are returned. If responseType==RESPONSE_TYPE_SUM, then the response functions are summed and the single summed response function is returned from evalModel(). If responseType==RESPONSE_TYPE_BLOCK, then the response functions are blocked up and returned as a product vector from evalModel().

ToDo: Finish Documentation!

Definition at line 396 of file Rythmos_ForwardSensitivityIntegratorAsModelEvaluator.hpp.

template<class Scalar >
const Array< Scalar > & Rythmos::ForwardSensitivityIntegratorAsModelEvaluator< Scalar >::getResponseTimes ( ) const
template<class Scalar >
void Rythmos::ForwardSensitivityIntegratorAsModelEvaluator< Scalar >::setParameterList ( RCP< Teuchos::ParameterList > const &  paramList)

Note that observationTargetIO() and parameterBaseIO() must be set before calling this function in order to use the parameter sublist to read in the vectors observationTarget() and parameterBase().

Definition at line 528 of file Rythmos_ForwardSensitivityIntegratorAsModelEvaluator.hpp.

template<class Scalar >
RCP< Teuchos::ParameterList > Rythmos::ForwardSensitivityIntegratorAsModelEvaluator< Scalar >::getNonconstParameterList ( )
template<class Scalar >
RCP< Teuchos::ParameterList > Rythmos::ForwardSensitivityIntegratorAsModelEvaluator< Scalar >::unsetParameterList ( )
template<class Scalar >
RCP< const Teuchos::ParameterList > Rythmos::ForwardSensitivityIntegratorAsModelEvaluator< Scalar >::getParameterList ( ) const
template<class Scalar >
RCP< const Teuchos::ParameterList > Rythmos::ForwardSensitivityIntegratorAsModelEvaluator< Scalar >::getValidParameters ( ) const

Note that observationTargetIO() and parameterBaseIO() must be set before calling this function in order to have the sublists added that will allow the vectors observationTarget() and parameterBase() to be read in latter when the parameter list is set..

Definition at line 573 of file Rythmos_ForwardSensitivityIntegratorAsModelEvaluator.hpp.

template<class Scalar >
int Rythmos::ForwardSensitivityIntegratorAsModelEvaluator< Scalar >::Np ( ) const
template<class Scalar >
int Rythmos::ForwardSensitivityIntegratorAsModelEvaluator< Scalar >::Ng ( ) const
template<class Scalar >
RCP< const Thyra::VectorSpaceBase< Scalar > > Rythmos::ForwardSensitivityIntegratorAsModelEvaluator< Scalar >::get_p_space ( int  l) const
template<class Scalar >
RCP< const Thyra::VectorSpaceBase< Scalar > > Rythmos::ForwardSensitivityIntegratorAsModelEvaluator< Scalar >::get_g_space ( int  j) const
template<class Scalar >
Thyra::ModelEvaluatorBase::InArgs< Scalar > Rythmos::ForwardSensitivityIntegratorAsModelEvaluator< Scalar >::createInArgs ( ) const

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