ModelEvaluator
subclass.
More...
#include <Rythmos_ForwardSensitivityModelEvaluator.hpp>
Inheritance diagram for Rythmos::ForwardSensitivityModelEvaluator< Scalar >:
Constructors/Intializers/Accessors  
ForwardSensitivityModelEvaluator ()  
 
void  initializeStructure (const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &stateModel, const int p_index) 
Intialize the with the model structure.  
Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > >  getStateModel () const 
 
int  get_p_index () const 
 
void  initializeState (const Thyra::ModelEvaluatorBase::InArgs< Scalar > &stateBasePoint, const Teuchos::RCP< const Thyra::LinearOpWithSolveBase< Scalar > > &W_tilde, const Scalar &coeff_x_dot, const Scalar &coeff_x, const Teuchos::RCP< const Thyra::LinearOpBase< Scalar > > &DfDx_dot=Teuchos::null, const Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > &DfDp=Teuchos::null) 
Initialize full state.  
Public functions overridden from ModelEvaulator.  
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 

ModelEvaluator
subclass.
This class provides a very general implemenation of a linear forward sensitivity model evaluator for a DAE.
f(x_dot(t),x(t),{p_l},t) = 0, over t = [t0,tf] x(t0) = x_init(p)
As shown above, the parameters are assumed to be steady state and can enter through the intial condition and/or through the DAE equation itself.
The forward sensitivity equations written in multivector form are:
F_sens(S_dot,S,t) = 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 multivector 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. The sensitivity parameter subvector p
here is really just one of the parameter subvectors in the state equation. This index of the parameter subvector for which the sensitivity equations are written for is given by p_index
. Note that above d(f)/d(x_dot)
, d(f)/d(x)
and d(f)/d(p
are all evaluated at the solution to the state equation (x_dot(t),x(t),t)
and are not functions of S_dot
or S
.
Since the model evaluator interface must be expressed in vector form, the multivector form of the forward sensitivity equations is flattened out into:
f_sens(s_bar_dot(t),s_bar(t),{p_l},t) = 0, over t = [t0,tf] s_bar(t0) = s_bar_init
where
s_bar = [ S(:,0); S(:,0); ...; S(:,np1) ] [ d(f)/d(x_dot)*S_dot(:,0) + d(f)/d(x)*S(:,0) + d(f)/d(p(0)) ] [ d(f)/d(x_dot)*S_dot(:,1) + d(f)/d(x)*S(:,1) + d(f)/d(p(1)) ] f_sens = [ ... ] [ d(f)/d(x_dot)*S_dot(:,np1) + d(f)/d(x)*S(:,np1) + d(f)/d(p(np1)) ] s_bar_init = [ d(x_init)/d(p(0)); d(x_init)/d(p(1)); ...; d(x_init)/d(p(np1)) ]
The product vector s_bar
is represented as a specialized Thyra::ProductVectorBase
subclass object with np
"blocks" in terms of a single Thyra::MultiVectorBase
object (which has np
columns).
W = alpha*d(f)/d(x_dot) + beta*d(f)/d(x)
which is the Jacobian for the nonlinear timestep equation in many methods.
First, let us consider how to represent the forward sensitivity equations using a precomputed LOWSB object
W_tilde = coeff_x_dot * d(f)/d(x_dot) + coeff_x * d(f)/d(x_dot)
computed at some point (x_dot_tilde,x_tilde,t_tilde,...)
.
Here is how to evaluate the forward sensitivity equations F_sens
using W_tilde
:
F_sens = d(f)/d(x_dot)*S_dot + ( (1/coeff_x) * ( W_tilde  coeff_x_dot * d(f)/d(x_dot) ) ) * S + d(f)/d(p) ==> F_sens = (1/coeff_x) * W_tilde * S + d(f)/d(x_dot) * ( S_dot  (coeff_x_dot/coeff_x)*S ) + d(f)/d(p)
Above, we see that in order to compute the multivector form of the residual for the forward sensitivity equations F_sens
, we must be able to compute something like:
d(f)/d(x_dot) * V + d(f)/d(p)
This can be done by computing d(f)/d(x_dot)
as W
with alpha=1.0
and beta=0.0
by calling the underlying state model. Or, a special sensitivity computation could be added to the model evaluator that would compute a generalization of:
F_sens_var = x_dot_scalar * d(f)/d(x_dot) * V_x_dot + x_scalar * d(f)/d(x) * V_x + p_scalar * d(f)/d(p) * V_p
We could then compute what we need using x_dot_scalar=1.0
, x_scalar=0.0
, p_scalar=1.0
, V_x_dot=V
, and V_p=I
. For explicit timestepping methods, this operation could compute the needed sensitivity in one shot. Such an addition to the Thyra::ModelEvaluator
interface would be handled through additions to the InArgs and OutArgs classes and the details of which are yet to be worked out.
Up to this point, the only assumption that we have made about the time stepping algorithm is that the sensitivities are only needed and only represented at a single point (x_dot_tilde,x_tilde,t_tilde,...)
. It is up to the client to ensure that this is indeed the case and it will be the case for many different types of timestepping methods.
In order to get maximum reuse out of a precomputed W_tilde
matrix we also have to put in special logic for the evaluation of W_hat
with respect to the forward sensitivity residual. When using a precomputed W_tilde
from the state timestep computation, the matrix W_hat
for the sensitivity residual becomes:
W_hat = alpha * d(f)/d(x_dot) + beta * (1/coeff_x) * ( W_tilde  coeff_x_dot * d(f)/d(x_dot) )
The special logic is that when alpha=coeff_x_dot
and beta=coeff_x
, then W_hat
above simplifies to W_hat=W_tilde
as:
W_hat = coeff_x_dot * d(f)/d(x_dot) + coeff_x * (1/coeff_x) * ( W_tilde  coeff_x_dot * d(f)/d(x_dot) ) ==> W_hat = coeff_x_dot * d(f)/d(x_dot)  coeff_x_dot * d(f)/d(x_dot) + W_tilde ==> W_hat = W_tilde
For all of the timestepping methods currently implemented in Rythmos at the time of this writing, alpha=coeff_x_dot
and beta=coeff_x
will always be true and this is checked for in the code implementation. Because of this, only a single W
LOWSB object will be given to a client and only it will be returned. Any other values of alpha
and beta
requested by the client will thown exceptions. In the future, other values of alpha
and beta
might be allowed but for now this would only be an error.
Note that there are also a few simplifications in cases where single residual timestepping methods are used. In the case of BDF methods (of which backward Euler is one type), we have:
S_dot = coeff_x_dot * S_tilde + B_x_dot S_dot = coeff_x * S_tilde
In this type of method we see that :
S_dot  (coeff_x_dot/coeff_x) * S ==> coeff_x_dot * S_tilde + B_x_dot  (coeff_x_dot/coeff_x) * ( coeff_x * S_tilde ) ==> coeff_x_dot * S_tilde  coeff_x_dot * S_tilde + B_x_dot ==> B_x_dot
Therefore, in this types of method, the term involving d(f)/d(x_dot)
becomes:
d(f)/d(x_dot) * ( S_dot  (coeff_x_dot/coeff_x)*S ) ==> d(f)/d(x_dot) * B_x_dot
and is independent of the unknown quantity S_tilde
. What this means is that if the residual for the sensitivity equaitions is to be computed multiple times for different values of S_tilde
, the term d(f)/d(x_dot) * B_x_dot
need only be computed once and can then be reused each time.
ToDo: Finish documention!
Definition at line 297 of file Rythmos_ForwardSensitivityModelEvaluator.hpp.
Rythmos::ForwardSensitivityModelEvaluator< Scalar >::ForwardSensitivityModelEvaluator  (  ) 
void Rythmos::ForwardSensitivityModelEvaluator< Scalar >::initializeStructure  (  const Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > &  stateModel,  
const int  p_index  
) 
Intialize the with the model structure.
stateModel  [in,persisting] The ModelEvaluator that defines the parameterized state model f(x_dot,x,p) .  
p_index  [in] The index of the parameter subvector in stateModel for which sensitivities will be computed for. 
*this
model object is not fully initialized at this point in that evalModel()
will not work yet and will thrown exceptions if called. The function initalizeState()
must be called later in order to fully initalize the model.
Definition at line 471 of file Rythmos_ForwardSensitivityModelEvaluator.hpp.
Teuchos::RCP< const Thyra::ModelEvaluator< Scalar > > Rythmos::ForwardSensitivityModelEvaluator< Scalar >::getStateModel  (  )  const 
int Rythmos::ForwardSensitivityModelEvaluator< Scalar >::get_p_index  (  )  const 
void Rythmos::ForwardSensitivityModelEvaluator< Scalar >::initializeState  (  const Thyra::ModelEvaluatorBase::InArgs< Scalar > &  stateBasePoint,  
const Teuchos::RCP< const Thyra::LinearOpWithSolveBase< Scalar > > &  W_tilde,  
const Scalar &  coeff_x_dot,  
const Scalar &  coeff_x,  
const Teuchos::RCP< const Thyra::LinearOpBase< Scalar > > &  DfDx_dot = Teuchos::null , 

const Teuchos::RCP< const Thyra::MultiVectorBase< Scalar > > &  DfDp = Teuchos::null  
) 
Initialize full state.
stateBasePoint  [in] The base point (x_dot_tilde,x_tilde,t_tilde,...) for which the sensitivities will be computed and represented at. Any sensitivities that are needed are computed at this point. The value of alpha and beta here are ignored.  
W_tilde  [in,persisting] The composite state derivative matrix computed at the base point stateBasePoint with alpha=coeff_x_dot and beta=coeff_x . This argument can be null, in which case it can be computed internally if needed or not at all.  
coeff_x_dot  [in] The value of alpha for which W_tilde was computed.  
coeff_x  [in] The value of beta for which W_tilde was computed.  
DfDx_dot  [in] The matrix d(f)/d(x_dot) computed at stateBasePoint if available. If this argument is null, then it will be computed internally if needed. The value is Teuchos::null .  
DfDp  [in] The matrix d(f)/d(p(p_index)) computed at stateBasePoint if available. If this argument is null, then it will be computed internally if needed. The value is Teuchos::null . 
!is_null(W_tilde)
This function must be called after intializeStructure()
and before evalModel()
is called. After this function is called, *this
model is fully initialized and ready to compute the requested outputs.
Definition at line 530 of file Rythmos_ForwardSensitivityModelEvaluator.hpp.
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > Rythmos::ForwardSensitivityModelEvaluator< Scalar >::get_x_space  (  )  const [virtual] 
Implements Thyra::ModelEvaluator< Scalar >.
Definition at line 588 of file Rythmos_ForwardSensitivityModelEvaluator.hpp.
Teuchos::RCP< const Thyra::VectorSpaceBase< Scalar > > Rythmos::ForwardSensitivityModelEvaluator< Scalar >::get_f_space  (  )  const [virtual] 
Implements Thyra::ModelEvaluator< Scalar >.
Definition at line 596 of file Rythmos_ForwardSensitivityModelEvaluator.hpp.
Thyra::ModelEvaluatorBase::InArgs< Scalar > Rythmos::ForwardSensitivityModelEvaluator< Scalar >::getNominalValues  (  )  const [virtual] 
Reimplemented from Thyra::StateFuncModelEvaluatorBase< Scalar >.
Definition at line 604 of file Rythmos_ForwardSensitivityModelEvaluator.hpp.
Teuchos::RCP< Thyra::LinearOpWithSolveBase< Scalar > > Rythmos::ForwardSensitivityModelEvaluator< Scalar >::create_W  (  )  const [virtual] 
Reimplemented from Thyra::StateFuncModelEvaluatorBase< Scalar >.
Definition at line 612 of file Rythmos_ForwardSensitivityModelEvaluator.hpp.
Thyra::ModelEvaluatorBase::InArgs< Scalar > Rythmos::ForwardSensitivityModelEvaluator< Scalar >::createInArgs  (  )  const [virtual] 
Implements Thyra::ModelEvaluator< Scalar >.
Definition at line 620 of file Rythmos_ForwardSensitivityModelEvaluator.hpp.