Thyra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Related Functions | List of all members
Thyra::EpetraModelEvaluator Class Reference

Concrete Adapter subclass that takes an EpetraExt::ModelEvaluator object and wraps it as a Thyra::ModelEvaluator object. More...

#include <Thyra_EpetraModelEvaluator.hpp>

Inheritance diagram for Thyra::EpetraModelEvaluator:
Inheritance graph
[legend]

Related Functions

(Note that these are not member functions.)

RCP< EpetraModelEvaluatorepetraModelEvaluator (const RCP< const EpetraExt::ModelEvaluator > &epetraModel, const RCP< LinearOpWithSolveFactoryBase< double > > &W_factory)
 
ModelEvaluatorBase::EDerivativeMultiVectorOrientation convert (const EpetraExt::ModelEvaluator::EDerivativeMultiVectorOrientation &mvOrientation)
 
EpetraExt::ModelEvaluator::EDerivativeMultiVectorOrientation convert (const ModelEvaluatorBase::EDerivativeMultiVectorOrientation &mvOrientation)
 
ModelEvaluatorBase::DerivativeProperties convert (const EpetraExt::ModelEvaluator::DerivativeProperties &derivativeProperties)
 
ModelEvaluatorBase::DerivativeSupport convert (const EpetraExt::ModelEvaluator::DerivativeSupport &derivativeSupport)
 
EpetraExt::ModelEvaluator::Derivative convert (const ModelEvaluatorBase::Derivative< double > &derivative, const RCP< const Epetra_Map > &fnc_map, const RCP< const Epetra_Map > &var_map)
 

Constructors/initializers/accessors/utilities.

 EpetraModelEvaluator ()
 
 EpetraModelEvaluator (const RCP< const EpetraExt::ModelEvaluator > &epetraModel, const RCP< LinearOpWithSolveFactoryBase< double > > &W_factory)
 
void initialize (const RCP< const EpetraExt::ModelEvaluator > &epetraModel, const RCP< LinearOpWithSolveFactoryBase< double > > &W_factory)
 
RCP< const
EpetraExt::ModelEvaluator > 
getEpetraModel () const
 
void setNominalValues (const ModelEvaluatorBase::InArgs< double > &nominalValues)
 Set the nominal values. More...
 
void setStateVariableScalingVec (const RCP< const Epetra_Vector > &stateVariableScalingVec)
 Set the state variable scaling vector s_x (see above). More...
 
RCP< const Epetra_VectorgetStateVariableInvScalingVec () const
 Get the state variable scaling vector s_x (see above). More...
 
RCP< const Epetra_VectorgetStateVariableScalingVec () const
 Get the inverse state variable scaling vector inv_s_x (see above). More...
 
void setStateFunctionScalingVec (const RCP< const Epetra_Vector > &stateFunctionScalingVec)
 Set the state function scaling vector s_f (see above). More...
 
RCP< const Epetra_VectorgetStateFunctionScalingVec () const
 Get the state function scaling vector s_f (see above). More...
 
void uninitialize (RCP< const EpetraExt::ModelEvaluator > *epetraModel=NULL, RCP< LinearOpWithSolveFactoryBase< double > > *W_factory=NULL)
 
const
ModelEvaluatorBase::InArgs
< double > & 
getFinalPoint () const
 
bool finalPointWasSolved () const
 

Public functions overridden from Teuchos::Describable.

std::string description () const
 

Overridden from ParameterListAcceptor

void setParameterList (RCP< Teuchos::ParameterList > const &paramList)
 
RCP< Teuchos::ParameterListgetNonconstParameterList ()
 
RCP< Teuchos::ParameterListunsetParameterList ()
 
RCP< const Teuchos::ParameterListgetParameterList () const
 
RCP< const Teuchos::ParameterListgetValidParameters () const
 

Public functions overridden from ModelEvaulator.

int Np () const
 
int Ng () const
 
RCP< const VectorSpaceBase
< double > > 
get_x_space () const
 
RCP< const VectorSpaceBase
< double > > 
get_f_space () const
 
RCP< const VectorSpaceBase
< double > > 
get_p_space (int l) const
 
RCP< const Teuchos::Array
< std::string > > 
get_p_names (int l) const
 
RCP< const VectorSpaceBase
< double > > 
get_g_space (int j) const
 
Teuchos::ArrayView< const
std::string > 
get_g_names (int j) const
 
ModelEvaluatorBase::InArgs
< double > 
getNominalValues () const
 
ModelEvaluatorBase::InArgs
< double > 
getLowerBounds () const
 
ModelEvaluatorBase::InArgs
< double > 
getUpperBounds () const
 
RCP< LinearOpBase< double > > create_W_op () const
 
RCP< PreconditionerBase< double > > create_W_prec () const
 Returns null currently. More...
 
RCP< const
LinearOpWithSolveFactoryBase
< double > > 
get_W_factory () const
 
ModelEvaluatorBase::InArgs
< double > 
createInArgs () const
 
void reportFinalPoint (const ModelEvaluatorBase::InArgs< double > &finalPoint, const bool wasSolved)
 

Additional Inherited Members

- Public Types inherited from Thyra::ModelEvaluator< double >
typedef Teuchos::ScalarTraits
< double >::magnitudeType 
ScalarMag
 
- Public Types inherited from Thyra::ModelEvaluatorBase
enum  EInArgsMembers {
  IN_ARG_x_dot_dot, IN_ARG_x_dot, IN_ARG_x, IN_ARG_x_dot_poly,
  IN_ARG_x_poly, IN_ARG_x_dot_mp, IN_ARG_x_mp, IN_ARG_t,
  IN_ARG_alpha, IN_ARG_beta, IN_ARG_W_x_dot_dot_coeff, IN_ARG_step_size,
  IN_ARG_stage_number
}
 
enum  EInArgs_p_mp { IN_ARG_p_mp }
 
enum  EEvalType { , EVAL_TYPE_APPROX_DERIV, EVAL_TYPE_VERY_APPROX_DERIV }
 The type of an evaluation. More...
 
enum  EDerivativeMultiVectorOrientation { DERIV_MV_JACOBIAN_FORM, DERIV_MV_GRADIENT_FORM, DERIV_MV_BY_COL = DERIV_MV_JACOBIAN_FORM, DERIV_TRANS_MV_BY_ROW = DERIV_MV_GRADIENT_FORM }
 
enum  EDerivativeLinearOp { DERIV_LINEAR_OP }
 
enum  EDerivativeLinearity { DERIV_LINEARITY_UNKNOWN, DERIV_LINEARITY_CONST, DERIV_LINEARITY_NONCONST }
 
enum  ERankStatus { DERIV_RANK_UNKNOWN, DERIV_RANK_FULL, DERIV_RANK_DEFICIENT }
 
enum  EOutArgsMembers {
  OUT_ARG_f, OUT_ARG_W, OUT_ARG_f_mp, OUT_ARG_W_mp,
  OUT_ARG_W_op, OUT_ARG_W_prec, OUT_ARG_f_poly
}
 
enum  EOutArgsDfDp { OUT_ARG_DfDp }
 
enum  EOutArgsDgDx_dot { OUT_ARG_DgDx_dot }
 
enum  EOutArgsDgDx { OUT_ARG_DgDx }
 
enum  EOutArgsDgDp { OUT_ARG_DgDp }
 
enum  EOutArgs_hess_vec_prod_f_xx { OUT_ARG_hess_vec_prod_f_xx }
 
enum  EOutArgs_hess_vec_prod_f_xp { OUT_ARG_hess_vec_prod_f_xp }
 
enum  EOutArgs_hess_vec_prod_f_px { OUT_ARG_hess_vec_prod_f_px }
 
enum  EOutArgs_hess_vec_prod_f_pp { OUT_ARG_hess_vec_prod_f_pp }
 
enum  EOutArgs_hess_vec_prod_g_xx { OUT_ARG_hess_vec_prod_g_xx }
 
enum  EOutArgs_hess_vec_prod_g_xp { OUT_ARG_hess_vec_prod_g_xp }
 
enum  EOutArgs_hess_vec_prod_g_px { OUT_ARG_hess_vec_prod_g_px }
 
enum  EOutArgs_hess_vec_prod_g_pp { OUT_ARG_hess_vec_prod_g_pp }
 
enum  EOutArgs_hess_f_xx { OUT_ARG_hess_f_xx }
 
enum  EOutArgs_hess_f_xp { OUT_ARG_hess_f_xp }
 
enum  EOutArgs_hess_f_pp { OUT_ARG_hess_f_pp }
 
enum  EOutArgs_hess_g_xx { OUT_ARG_hess_g_xx }
 
enum  EOutArgs_hess_g_xp { OUT_ARG_hess_g_xp }
 
enum  EOutArgs_hess_g_pp { OUT_ARG_hess_g_pp }
 
enum  EOutArgs_H_xx { OUT_ARG_H_xx }
 
enum  EOutArgs_H_xp { OUT_ARG_H_xp }
 
enum  EOutArgs_H_pp { OUT_ARG_H_pp }
 
enum  EOutArgsDfDp_mp { OUT_ARG_DfDp_mp }
 
enum  EOutArgs_g_mp { OUT_ARG_g_mp }
 
enum  EOutArgsDgDx_dot_mp { OUT_ARG_DgDx_dot_mp }
 
enum  EOutArgsDgDx_mp { OUT_ARG_DgDx_mp }
 
enum  EOutArgsDgDp_mp { OUT_ARG_DgDp_mp }
 
- Public Member Functions inherited from Thyra::ModelEvaluatorDefaultBase< double >
int Np () const
 
int Ng () const
 
RCP< LinearOpBase< double > > create_DfDp_op (int l) const
 
RCP< LinearOpBase< double > > create_DgDx_dot_op (int j) const
 
RCP< LinearOpBase< double > > create_DgDx_op (int j) const
 
RCP< LinearOpBase< double > > create_DgDp_op (int j, int l) const
 
RCP< LinearOpWithSolveBase
< double > > 
create_W () const
 
ModelEvaluatorBase::OutArgs
< double > 
createOutArgs () const
 
void evalModel (const ModelEvaluatorBase::InArgs< double > &inArgs, const ModelEvaluatorBase::OutArgs< double > &outArgs) const
 
virtual RCP< const
VectorSpaceBase< double > > 
get_f_multiplier_space () const
 
virtual RCP< const
VectorSpaceBase< double > > 
get_g_multiplier_space (int j) const
 
virtual RCP< LinearOpBase
< double > > 
create_hess_f_xx () const
 
virtual RCP< LinearOpBase
< double > > 
create_hess_f_xp (int l) const
 
virtual RCP< LinearOpBase
< double > > 
create_hess_f_pp (int l1, int l2) const
 
virtual RCP< LinearOpBase
< double > > 
create_hess_g_xx (int j) const
 
virtual RCP< LinearOpBase
< double > > 
create_hess_g_xp (int j, int l) const
 
virtual RCP< LinearOpBase
< double > > 
create_hess_g_pp (int j, int l1, int l2) const
 
- Public Member Functions inherited from Thyra::ModelEvaluator< double >
- Public Member Functions inherited from Thyra::ModelEvaluatorBase
 ModelEvaluatorBase ()
 constructor More...
 
- Static Public Attributes inherited from Thyra::ModelEvaluatorBase
static const int NUM_E_IN_ARGS_MEMBERS =13
 
static const int NUM_E_OUT_ARGS_MEMBERS =7
 
- Protected Member Functions inherited from Thyra::ModelEvaluatorDefaultBase< double >
 ModelEvaluatorDefaultBase ()
 
void initializeDefaultBase ()
 Function called by subclasses to fully initialize this object on any important change. More...
 
void resetDefaultBase ()
 Sets the the DefaultBase to an uninitialized state, forcing lazy initialization when needed. More...
 

Detailed Description

Concrete Adapter subclass that takes an EpetraExt::ModelEvaluator object and wraps it as a Thyra::ModelEvaluator object.

This class takes care of the basic details of wrapping and unwrapping Epetra from Thyra objects. This class is highly configurable and will be maintained and modified in the future as the basic link between the Epetra world and the Thyra world for nonlinear models and nonlinear algorithms.

Scaling

This class can handle scaling of the state function f(...) and of the state variables x and all of the affected derivatives.

The scaling for the state function can be set manually using setStateFunctionScalingVec() or can be computed automatically using the parameter "State Function Scaling" (see documentation output from this->getValidParameters()->print(...)) in the input parameter list set by setParameterList(). Reguardless of how the state function scaling is computed, it will compute a positive vector s_f that defines a diagonal matrix S_f = diag(s_f) that transforms the state function:

   f(...) = S_f * f_orig(...)

where f_orig(...) is the original state function as computed by the underlying EpetraExt::ModelEvaluator object and f(...) is the state function as computed by evalModel().

The scaling for the state variables must be set manually using Thyra::setStateVariableScalingVec(). The vector that is set s_x>/tt> defines a diagonal scaling matrix S_x = diag(s_x) that transforms the variables as:

   x = S_x * x_orig

where x_orig is the original unscaled state variable vector as defined by the underlying EpetraExt::ModelEvaluator object and x is the scaled state varaible vector as returned from getNominalValues() and as accepted by evalModel(). Note that when the scaled variables x are passed into evalModel that they are unscaled as:

   x_orig = inv(S_x) * x

where inv(S_x) is the inverse of the diagonals of S_x which is stored as a vector inv_s_x.

Note how these scalings affect the state function:

   f(x,...) = S_f * f_orig( inv(S_x)*x...)

which as the state/state Jacobian:

   W = d(f)/d(x) = S_f * d(f_orig)/d(x_orig) * inv(S_x)

Currently, this class does not handle scalings of the parameters p(l) or of the auxilary response functions g(j)(...).

The state varaible and state function scaling gives the following scaled quantities:

   f = S_f * f_orig

   W = S_f * W_orig * inv(S_x)

   DfDp(l) = S_f * DfDp_orig(l)

   g(j) = g_orig(j)

   DgDx_dot(j) = DgDx_dot_orig(j) * inv(S_x)

   DgDx(j) = DgDx_orig(j) * inv(S_x)
   
   DgDp(j,l) = DgDp_orig(j,l)

Since the scaling is done explicitly, the client never even sees the orginal scaling and the linear solver (and contained preconditioner) are computed from the scaled W shown above.

ToDo: Describe how scaling affects the Hessian-vector products an how you just need to scale the Lagrange mutipliers as:

 u^T * f(...) = u^T * (S_f * f_orig(...)) = u_f^T * f_orig(...)

where u_f = S_f * u.

ToDo: Finish documentation!

Definition at line 175 of file Thyra_EpetraModelEvaluator.hpp.

Constructor & Destructor Documentation

Thyra::EpetraModelEvaluator::EpetraModelEvaluator ( )

Definition at line 123 of file Thyra_EpetraModelEvaluator.cpp.

Thyra::EpetraModelEvaluator::EpetraModelEvaluator ( const RCP< const EpetraExt::ModelEvaluator > &  epetraModel,
const RCP< LinearOpWithSolveFactoryBase< double > > &  W_factory 
)

Definition at line 129 of file Thyra_EpetraModelEvaluator.cpp.

Member Function Documentation

void Thyra::EpetraModelEvaluator::initialize ( const RCP< const EpetraExt::ModelEvaluator > &  epetraModel,
const RCP< LinearOpWithSolveFactoryBase< double > > &  W_factory 
)

Definition at line 140 of file Thyra_EpetraModelEvaluator.cpp.

RCP< const EpetraExt::ModelEvaluator > Thyra::EpetraModelEvaluator::getEpetraModel ( ) const

Definition at line 197 of file Thyra_EpetraModelEvaluator.cpp.

void Thyra::EpetraModelEvaluator::setNominalValues ( const ModelEvaluatorBase::InArgs< double > &  nominalValues)

Set the nominal values.

Warning, if scaling is being used, these must be according to the scaled values, not the original unscaled values.

Definition at line 203 of file Thyra_EpetraModelEvaluator.cpp.

void Thyra::EpetraModelEvaluator::setStateVariableScalingVec ( const RCP< const Epetra_Vector > &  stateVariableScalingVec)

Set the state variable scaling vector s_x (see above).

This function must be called after intialize() or the constructur in order to set the scaling vector correctly!

ToDo: Move this into an external strategy class object!

Definition at line 212 of file Thyra_EpetraModelEvaluator.cpp.

RCP< const Epetra_Vector > Thyra::EpetraModelEvaluator::getStateVariableInvScalingVec ( ) const

Get the state variable scaling vector s_x (see above).

Definition at line 234 of file Thyra_EpetraModelEvaluator.cpp.

RCP< const Epetra_Vector > Thyra::EpetraModelEvaluator::getStateVariableScalingVec ( ) const

Get the inverse state variable scaling vector inv_s_x (see above).

Definition at line 227 of file Thyra_EpetraModelEvaluator.cpp.

void Thyra::EpetraModelEvaluator::setStateFunctionScalingVec ( const RCP< const Epetra_Vector > &  stateFunctionScalingVec)

Set the state function scaling vector s_f (see above).

Definition at line 241 of file Thyra_EpetraModelEvaluator.cpp.

RCP< const Epetra_Vector > Thyra::EpetraModelEvaluator::getStateFunctionScalingVec ( ) const

Get the state function scaling vector s_f (see above).

Definition at line 250 of file Thyra_EpetraModelEvaluator.cpp.

void Thyra::EpetraModelEvaluator::uninitialize ( RCP< const EpetraExt::ModelEvaluator > *  epetraModel = NULL,
RCP< LinearOpWithSolveFactoryBase< double > > *  W_factory = NULL 
)

Definition at line 256 of file Thyra_EpetraModelEvaluator.cpp.

const ModelEvaluatorBase::InArgs< double > & Thyra::EpetraModelEvaluator::getFinalPoint ( ) const

Definition at line 273 of file Thyra_EpetraModelEvaluator.cpp.

bool Thyra::EpetraModelEvaluator::finalPointWasSolved ( ) const

Definition at line 279 of file Thyra_EpetraModelEvaluator.cpp.

std::string Thyra::EpetraModelEvaluator::description ( ) const
virtual

Reimplemented from Teuchos::Describable.

Definition at line 288 of file Thyra_EpetraModelEvaluator.cpp.

void Thyra::EpetraModelEvaluator::setParameterList ( RCP< Teuchos::ParameterList > const &  paramList)
virtual

Implements Teuchos::ParameterListAcceptor.

Definition at line 310 of file Thyra_EpetraModelEvaluator.cpp.

RCP< Teuchos::ParameterList > Thyra::EpetraModelEvaluator::getNonconstParameterList ( )
virtual

Implements Teuchos::ParameterListAcceptor.

Definition at line 331 of file Thyra_EpetraModelEvaluator.cpp.

RCP< Teuchos::ParameterList > Thyra::EpetraModelEvaluator::unsetParameterList ( )
virtual

Implements Teuchos::ParameterListAcceptor.

Definition at line 338 of file Thyra_EpetraModelEvaluator.cpp.

RCP< const Teuchos::ParameterList > Thyra::EpetraModelEvaluator::getParameterList ( ) const
virtual

Reimplemented from Teuchos::ParameterListAcceptor.

Definition at line 347 of file Thyra_EpetraModelEvaluator.cpp.

RCP< const Teuchos::ParameterList > Thyra::EpetraModelEvaluator::getValidParameters ( ) const
virtual

Reimplemented from Teuchos::ParameterListAcceptor.

Definition at line 354 of file Thyra_EpetraModelEvaluator.cpp.

int Thyra::EpetraModelEvaluator::Np ( ) const
virtual

Implements Thyra::ModelEvaluator< double >.

Definition at line 403 of file Thyra_EpetraModelEvaluator.cpp.

int Thyra::EpetraModelEvaluator::Ng ( ) const
virtual

Implements Thyra::ModelEvaluator< double >.

Definition at line 409 of file Thyra_EpetraModelEvaluator.cpp.

RCP< const VectorSpaceBase< double > > Thyra::EpetraModelEvaluator::get_x_space ( ) const
virtual

Implements Thyra::ModelEvaluator< double >.

Definition at line 416 of file Thyra_EpetraModelEvaluator.cpp.

RCP< const VectorSpaceBase< double > > Thyra::EpetraModelEvaluator::get_f_space ( ) const
virtual

Implements Thyra::ModelEvaluator< double >.

Definition at line 423 of file Thyra_EpetraModelEvaluator.cpp.

RCP< const VectorSpaceBase< double > > Thyra::EpetraModelEvaluator::get_p_space ( int  l) const
virtual

Implements Thyra::ModelEvaluator< double >.

Definition at line 430 of file Thyra_EpetraModelEvaluator.cpp.

RCP< const Teuchos::Array< std::string > > Thyra::EpetraModelEvaluator::get_p_names ( int  l) const
virtual

Implements Thyra::ModelEvaluator< double >.

Definition at line 440 of file Thyra_EpetraModelEvaluator.cpp.

RCP< const VectorSpaceBase< double > > Thyra::EpetraModelEvaluator::get_g_space ( int  j) const
virtual

Implements Thyra::ModelEvaluator< double >.

Definition at line 450 of file Thyra_EpetraModelEvaluator.cpp.

Teuchos::ArrayView< const std::string > Thyra::EpetraModelEvaluator::get_g_names ( int  j) const
virtual

Implements Thyra::ModelEvaluator< double >.

Definition at line 458 of file Thyra_EpetraModelEvaluator.cpp.

ModelEvaluatorBase::InArgs< double > Thyra::EpetraModelEvaluator::getNominalValues ( ) const
virtual

Implements Thyra::ModelEvaluator< double >.

Definition at line 468 of file Thyra_EpetraModelEvaluator.cpp.

ModelEvaluatorBase::InArgs< double > Thyra::EpetraModelEvaluator::getLowerBounds ( ) const
virtual

Implements Thyra::ModelEvaluator< double >.

Definition at line 476 of file Thyra_EpetraModelEvaluator.cpp.

ModelEvaluatorBase::InArgs< double > Thyra::EpetraModelEvaluator::getUpperBounds ( ) const
virtual

Implements Thyra::ModelEvaluator< double >.

Definition at line 484 of file Thyra_EpetraModelEvaluator.cpp.

RCP< LinearOpBase< double > > Thyra::EpetraModelEvaluator::create_W_op ( ) const
virtual

Implements Thyra::ModelEvaluator< double >.

Definition at line 492 of file Thyra_EpetraModelEvaluator.cpp.

RCP< PreconditionerBase< double > > Thyra::EpetraModelEvaluator::create_W_prec ( ) const
virtual

Returns null currently.

Implements Thyra::ModelEvaluator< double >.

Definition at line 499 of file Thyra_EpetraModelEvaluator.cpp.

RCP< const LinearOpWithSolveFactoryBase< double > > Thyra::EpetraModelEvaluator::get_W_factory ( ) const
virtual

.

Implements Thyra::ModelEvaluator< double >.

Definition at line 506 of file Thyra_EpetraModelEvaluator.cpp.

ModelEvaluatorBase::InArgs< double > Thyra::EpetraModelEvaluator::createInArgs ( ) const
virtual

Implements Thyra::ModelEvaluator< double >.

Definition at line 512 of file Thyra_EpetraModelEvaluator.cpp.

void Thyra::EpetraModelEvaluator::reportFinalPoint ( const ModelEvaluatorBase::InArgs< double > &  finalPoint,
const bool  wasSolved 
)
virtual

Implements Thyra::ModelEvaluator< double >.

Definition at line 520 of file Thyra_EpetraModelEvaluator.cpp.

Friends And Related Function Documentation

RCP< EpetraModelEvaluator > epetraModelEvaluator ( const RCP< const EpetraExt::ModelEvaluator > &  epetraModel,
const RCP< LinearOpWithSolveFactoryBase< double > > &  W_factory 
)
related

ModelEvaluatorBase::EDerivativeMultiVectorOrientation convert ( const EpetraExt::ModelEvaluator::EDerivativeMultiVectorOrientation &  mvOrientation)
related

EpetraExt::ModelEvaluator::EDerivativeMultiVectorOrientation convert ( const ModelEvaluatorBase::EDerivativeMultiVectorOrientation mvOrientation)
related

ModelEvaluatorBase::DerivativeProperties convert ( const EpetraExt::ModelEvaluator::DerivativeProperties &  derivativeProperties)
related

ModelEvaluatorBase::DerivativeSupport convert ( const EpetraExt::ModelEvaluator::DerivativeSupport &  derivativeSupport)
related

EpetraExt::ModelEvaluator::Derivative convert ( const ModelEvaluatorBase::Derivative< double > &  derivative,
const RCP< const Epetra_Map > &  fnc_map,
const RCP< const Epetra_Map > &  var_map 
)
related


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