42 #ifndef EPETRA_EXT_MODEL_EVALUATOR_HPP
43 #define EPETRA_EXT_MODEL_EVALUATOR_HPP
45 #if defined(EpetraExt_SHOW_DEPRECATED_WARNINGS)
47 #warning "The EpetraExt package is deprecated"
59 #ifdef HAVE_PYTRILINOS
62 typedef _object PyObject;
72 class EpetraVectorOrthogPoly;
73 class EpetraMultiVectorOrthogPoly;
74 class EpetraOperatorOrthogPoly;
76 template <
typename ordinal_type,
typename scalar_type>
class Quadrature;
80 class ProductEpetraVector;
81 class ProductEpetraMultiVector;
82 class ProductEpetraOperator;
213 void set_t(
double t );
227 double get_t()
const;
311 template<
class ObjType>
379 switch(mvOrientation) {
395 switch(mvOrientation) {
484 ) :
dmv_(mv,orientation) {}
516 bool isAlreadyInverted_ )
572 ) :
dmv_(mv,orientation) {}
643 ) :
dmv_(mv,orientation) {}
1282 #ifdef HAVE_PYTRILINOS
1284 friend InArgs convertInArgsFromPython(PyObject * source);
1287 friend OutArgs convertOutArgsFromPython(PyObject * source);
1402 const std::string &modelEvalDescription,
1404 const std::string &derivName
1410 const std::string &modelEvalDescription,
1412 const std::string &derivName,
1476 {
return p_.size(); }
1496 { assert_supports(
IN_ARG_x); x_ = x; }
1500 { assert_supports(
IN_ARG_x);
return x_; }
1585 { assert_l(l); p_[l] = p_l; }
1589 { assert_l(l);
return p_[l]; }
1594 { assert_supports(
IN_ARG_p_sg, l); p_sg_[l] = p_sg_l; }
1599 { assert_supports(
IN_ARG_p_sg, l);
return p_sg_[l]; }
1604 { assert_supports(
IN_ARG_p_mp, l); p_mp_[l] = p_mp_l; }
1609 { assert_supports(
IN_ARG_p_mp, l);
return p_mp_[l]; }
1613 { assert_supports(
IN_ARG_t); t_ = t; }
1617 { assert_supports(
IN_ARG_t);
return t_; }
1689 modelEvalDescription_ = new_modelEvalDescription;
1696 p_sg_.resize(new_Np);
1697 p_mp_.resize(new_Np);
1698 supports_p_sg_.resize(new_Np);
1699 supports_p_mp_.resize(new_Np);
1708 {
return modelEvalDescription_; }
1713 return DfDp_.size();
1786 {
return W_properties_; }
1789 {
return WPrec_properties_; }
1811 return DfDp_properties_[l];
1818 DfDp_sg_[l] = DfDp_sg_l;
1834 return DfDp_sg_properties_[l];
1841 DfDp_mp_[l] = DfDp_mp_l;
1857 return DfDp_mp_properties_[l];
1864 DgDx_dot_[j] = DgDx_dot_j;
1872 return DgDx_dot_[j];
1880 return DgDx_dot_properties_[j];
1887 DgDx_dot_sg_[j] = DgDx_dot_sg_j;
1895 return DgDx_dot_sg_[j];
1903 return DgDx_dot_sg_properties_[j];
1910 DgDx_dot_mp_[j] = DgDx_dot_mp_j;
1918 return DgDx_dot_mp_[j];
1926 return DgDx_dot_mp_properties_[j];
1933 DgDx_dotdot_[j] = DgDx_dotdot_j;
1941 return DgDx_dotdot_[j];
1949 return DgDx_dotdot_properties_[j];
1956 DgDx_dotdot_sg_[j] = DgDx_dotdot_sg_j;
1964 return DgDx_dotdot_sg_[j];
1972 return DgDx_dotdot_sg_properties_[j];
1979 DgDx_dotdot_mp_[j] = DgDx_dotdot_mp_j;
1987 return DgDx_dotdot_mp_[j];
1995 return DgDx_dotdot_mp_properties_[j];
2018 return DgDx_properties_[j];
2025 DgDx_sg_[j] = DgDx_sg_j;
2041 return DgDx_sg_properties_[j];
2048 DgDx_mp_[j] = DgDx_mp_j;
2064 return DgDx_mp_properties_[j];
2071 DgDp_[ j*Np() + l ] = DgDp_j_l;
2079 return DgDp_[ j*Np() + l ];
2087 return DgDp_properties_[ j*Np() + l ];
2094 DgDp_sg_[ j*Np() + l ] = DgDp_sg_j_l;
2102 return DgDp_sg_[ j*Np() + l ];
2110 return DgDp_sg_properties_[ j*Np() + l ];
2117 DgDp_mp_[ j*Np() + l ] = DgDp_mp_j_l;
2125 return DgDp_mp_[ j*Np() + l ];
2133 return DgDp_mp_properties_[ j*Np() + l ];
2138 { f_poly_ = f_poly; }
2182 this->_setModelEvalDescription(new_modelEvalDescription);
2187 { this->_set_Np(new_Np); }
2191 { this->_setSupports(arg,new_supports); }
2195 { this->_setSupports(arg,l,new_supports); }
2199 { this->_setSupports(arg,l,new_supports); }
2208 this->_setModelEvalDescription(new_modelEvalDescription);
2213 { this->_set_Np_Ng(new_Np,new_Ng); }
2217 { this->_setSupports(arg,new_supports); }
2221 { this->_setSupports(arg,l,new_supports); }
2225 { this->_setSupports(arg,j,new_supports); }
2229 { this->_setSupports(arg,j,new_supports); }
2233 { this->_setSupports(arg,j,new_supports); }
2237 { this->_setSupports(arg,j,l,new_supports); }
2241 { this->_setSupports(arg,j,new_supports); }
2245 { this->_setSupports(arg,l,new_supports); }
2249 { this->_setSupports(arg,j,new_supports); }
2253 { this->_setSupports(arg,j,new_supports); }
2257 { this->_setSupports(arg,j,new_supports); }
2261 { this->_setSupports(arg,j,l,new_supports); }
2265 { this->_setSupports(arg,j,new_supports); }
2269 { this->_setSupports(arg,l,new_supports); }
2273 { this->_setSupports(arg,j,new_supports); }
2277 { this->_setSupports(arg,j,new_supports); }
2281 { this->_setSupports(arg,j,new_supports); }
2285 { this->_setSupports(arg,j,l,new_supports); }
2289 { this->_set_W_properties(properties); }
2292 { this->_set_WPrec_properties(properties); }
2297 this->_set_DfDp_properties(l,properties);
2303 this->_set_DgDx_dot_properties(j,properties);
2309 this->_set_DgDx_dotdot_properties(j,properties);
2315 this->_set_DgDx_properties(j,properties);
2321 this->_set_DgDp_properties(j,l,properties);
2327 this->_set_DfDp_sg_properties(l,properties);
2333 this->_set_DgDx_dot_sg_properties(j,properties);
2339 this->_set_DgDx_dotdot_sg_properties(j,properties);
2345 this->_set_DgDx_sg_properties(j,properties);
2351 this->_set_DgDp_sg_properties(j,l,properties);
2357 this->_set_DfDp_mp_properties(l,properties);
2363 this->_set_DgDx_dot_mp_properties(j,properties);
2369 this->_set_DgDx_dotdot_mp_properties(j,properties);
2375 this->_set_DgDx_mp_properties(j,properties);
2381 this->_set_DgDp_mp_properties(j,l,properties);
2386 #endif // EPETRA_EXT_MODEL_EVALUATOR_HPP
Simple aggregate class that stores a derivative object as a general linear operator or as a multi-vec...
DerivativeMultiVector getDerivativeMultiVector() const
void set_stage_number(int stage_number)
deriv_properties_t DgDp_mp_properties_
mp_const_vector_t get_x_mp() const
Get multi-point solution vector.
Teuchos::RCP< const Teuchos::Polynomial< Epetra_Vector > > get_x_dot_poly() const
Get time derivative vector Taylor polynomial.
EDerivativeMultiVectorOrientation getMultiVectorOrientation() const
SGDerivative get_DgDx_dotdot_sg(int j) const
supports_t supports_DgDx_dotdot_sg_
bool supports(EOutArgsMembers arg) const
Residual vector Taylor polynomial.
supports_t supports_DgDx_dot_
MPDerivative(const mp_operator_t &lo)
void set_WPrec_properties(const DerivativeProperties &properties)
Teuchos::RCP< const Epetra_Vector > x_dotdot_
MPDerivativeMultiVector(const mp_multivector_t &mv, const EDerivativeMultiVectorOrientation orientation=DERIV_MV_BY_COL, const Teuchos::Array< int > ¶mIndexes=Teuchos::Array< int >())
virtual Teuchos::RCP< const Epetra_Map > get_g_map(int j) const
.
void set_W_properties(const DerivativeProperties &properties)
DerivativeProperties get_DgDx_sg_properties(int j) const
virtual InArgs createInArgs() const =0
MPDerivativeMultiVector dmv_
void set_step_size(double step_size)
mp_const_vector_t get_p_mp(int l) const
Get multi-point parameter vector.
Teuchos::RCP< const Teuchos::Polynomial< Epetra_Vector > > x_poly_
virtual Teuchos::RCP< Epetra_Operator > create_DgDx_op(int j) const
Evaluation< Epetra_Vector > get_g(int j) const
Get g(j) where 0 <= j && j < this->Ng().
mp_const_vector_t get_x_dotdot_mp() const
void set_DgDx(int j, const Derivative &DgDx_j)
Derivative get_DfDp(int l) const
void setSupports(EOutArgsMembers arg, bool supports=true)
std::string modelEvalDescription_
EEvalType getType() const
supports_g_sg_t supports_g_mp_
A very approximate derivative (i.e. for a preconditioner)
RCP(ObjType *p, ERCPWeakNoDealloc)
Teuchos::RCP< const Epetra_Vector > get_x_dotdot() const
void set_x_dotdot_poly(const Teuchos::RCP< const Teuchos::Polynomial< Epetra_Vector > > &x_dotdot_poly)
DerivativeMultiVector dmv_
mp_vector_t get_g_mp(int j) const
Get multi-point response.
Preconditioner()
Default constructor of null Operatir.
bool supports(EDerivativeLinearOp) const
Time second derivative vector Taylor polynomial.
std::string modelEvalDescription_
void _setSupports(EOutArgsMembers arg, bool supports)
Teuchos::Array< DerivativeSupport > supports_t
void _set_DfDp_properties(int l, const DerivativeProperties &properties)
virtual Teuchos::RCP< const Epetra_Vector > get_p_init(int l) const
Exact function evaluation.
void set_DgDx_sg_properties(int j, const DerivativeProperties &properties)
virtual Teuchos::RCP< const Epetra_Map > get_f_map() const =0
.
supports_t supports_DgDp_
virtual Teuchos::RCP< const Epetra_Map > get_x_map() const =0
.
Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > sg_operator_t
Short-hand for stochastic Galerkin operator type.
MPDerivativeMultiVector()
Teuchos::Array< mp_vector_t > g_mp_t
DerivativeProperties WPrec_properties_
void set_beta(double beta)
virtual ~ModelEvaluator()
static const int NUM_E_IN_ARGS_MEMBERS
deriv_properties_t DfDp_mp_properties_
bool isAlreadyInverted
Bool flag.
Teuchos::RCP< Teuchos::Polynomial< Epetra_Vector > > get_f_poly() const
Get residual vector Taylor polynomial.
Preconditioner(const Teuchos::RCP< Epetra_Operator > &PrecOp_, bool isAlreadyInverted_)
Usable constructor to set the (Epetra_Operator,bool) pair.
void set_x(const Teuchos::RCP< const Epetra_Vector > &x)
void set_x_dot(const Teuchos::RCP< const Epetra_Vector > &x_dot)
virtual Teuchos::RCP< const Epetra_Vector > get_p_upper_bounds(int l) const
void set_x_dotdot_mp(const mp_const_vector_t &x_dotdot_mp)
Teuchos::RCP< Epetra_Operator > getLinearOp(const std::string &modelEvalDescription, const ModelEvaluator::Derivative &deriv, const std::string &derivName)
Multi-point solution vector.
Teuchos::RCP< Epetra_MultiVector > get_DgDp_mv(const int j, const int l, const ModelEvaluator::OutArgs &outArgs, const ModelEvaluator::EDerivativeMultiVectorOrientation mvOrientation)
virtual double get_t_upper_bound() const
mp_multivector_t getMultiVector() const
Teuchos::RCP< const Stokhos::OrthogPolyBasis< int, double > > sg_basis_
Coeff of second derivative term d(x_dotdot)/dx.
virtual Teuchos::RCP< const Teuchos::Array< std::string > > get_p_names(int l) const
Get the names of the parameters associated with parameter subvector l if available.
SGDerivativeMultiVector(const Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > &mv, const EDerivativeMultiVectorOrientation orientation=DERIV_MV_BY_COL, const Teuchos::Array< int > ¶mIndexes=Teuchos::Array< int >())
bool funcOrDerivesAreSet(EOutArgsMembers arg) const
Return true if the function or its derivatives are set.
DerivativeProperties get_DfDp_sg_properties(int l) const
supports_t supports_DfDp_sg_
MPDerivative get_DgDx_dotdot_mp(int j) const
void assert_l(int l) const
EDerivativeMultiVectorOrientation orientation_
Teuchos::RCP< Epetra_Operator > get_W() const
Teuchos::Array< MPDerivative > mp_deriv_t
EDerivativeMultiVectorOrientation orientation_
DerivativeProperties get_DgDx_mp_properties(int j) const
Teuchos::RCP< Epetra_MultiVector > get_DgDx_dot_mv(const int j, const ModelEvaluator::OutArgs &outArgs, const ModelEvaluator::EDerivativeMultiVectorOrientation mvOrientation)
void set_DfDp(int l, const Derivative &DfDp_l)
void set_WPrec(const Teuchos::RCP< Epetra_Operator > &WPrec)
Teuchos::RCP< Stokhos::ProductEpetraMultiVector > mp_multivector_t
Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > getMultiVector() const
Stochastic Galerkin quadrature.
Stochastic Galerkin expansion.
void _set_DgDp_mp_properties(int j, int l, const DerivativeProperties &properties)
Teuchos::RCP< Epetra_Operator > W_
Teuchos::RCP< Stokhos::ProductEpetraOperator > mp_operator_t
deriv_properties_t DfDp_sg_properties_
void set_DgDx_dotdot(int j, const Derivative &DgDx_dotdot_j)
Derivative get_DgDx_dot(int j) const
void set_g_sg(int j, const sg_vector_t &g_sg_j)
Set stochastic Galerkin vector polynomial response.
Teuchos::RCP< const Stokhos::ProductEpetraMultiVector > mp_const_multivector_t
void set_f(const Evaluation< Epetra_Vector > &f)
SGDerivativeMultiVector getDerivativeMultiVector() const
void setSupports(EInArgsMembers arg, bool supports=true)
Teuchos::RCP< Epetra_Operator > PrecOp
Accessor for the Epetra_Operator.
Teuchos::RCP< const Teuchos::Polynomial< Epetra_Vector > > get_x_poly() const
Get solution vector Taylor polynomial.
DerivativeProperties get_WPrec_properties() const
DerivativeSupport(EDerivativeMultiVectorOrientation mvOrientation)
DerivativeProperties get_DgDp_properties(int j, int l) const
const Teuchos::Array< int > & getParamIndexes() const
DerivativeProperties get_DgDx_properties(int j) const
Teuchos::RCP< Epetra_Operator > getLinearOp() const
bool supports(EDerivativeMultiVectorOrientation mvOrientation) const
DerivativeProperties get_DgDp_mp_properties(int j, int l) const
Teuchos::Array< sg_vector_t > g_sg_t
Teuchos::Array< int > paramIndexes_
DerivativeProperties get_W_properties() const
supports_t supports_DgDx_sg_
void setFailed() const
Set that the evaluation as a whole failed.
Teuchos::RCP< const Epetra_Vector > get_p(int l) const
void set_DfDp_sg(int l, const SGDerivative &DfDp_sg_l)
DerivativeProperties get_DgDx_dotdot_sg_properties(int j) const
Simple aggregate class that stores a derivative object as a general linear operator or as a multi-vec...
EDerivativeMultiVectorOrientation getOrientation() const
MPDerivative get_DfDp_mp(int l) const
Evaluation< Epetra_Vector > f_
Teuchos::RCP< const Epetra_Vector > x_dot_
EDerivativeMultiVectorOrientation orientation_
supports_t supports_DgDx_mp_
Teuchos::RCP< const Stokhos::Quadrature< int, double > > get_sg_quadrature() const
SGDerivative get_DgDp_sg(int j, int l) const
Derivative get_DgDx(int j) const
Teuchos::RCP< Epetra_MultiVector > getMultiVector(const std::string &modelEvalDescription, const ModelEvaluator::Derivative &deriv, const std::string &derivName, const ModelEvaluator::EDerivativeMultiVectorOrientation mvOrientation)
virtual Teuchos::RCP< Epetra_Operator > create_DgDx_dot_op(int j) const
void changeOrientation(const EDerivativeMultiVectorOrientation orientation)
deriv_properties_t DgDx_dotdot_mp_properties_
virtual Teuchos::RCP< const Epetra_Vector > get_x_upper_bounds() const
Simple aggregate class that stores a derivative object as a general linear operator or as a multi-vec...
deriv_properties_t DgDx_dot_mp_properties_
supports_t supports_DgDx_dot_sg_
void set_sg_basis(const Teuchos::RCP< const Stokhos::OrthogPolyBasis< int, double > > &basis)
EDerivativeLinearity linearity
Teuchos::RCP< Stokhos::OrthogPolyExpansion< int, double, Stokhos::StandardStorage< int, double > > > sg_exp_
Multi-point "W" operator.
deriv_properties_t DgDx_dot_properties_
DerivativeSupport(EDerivativeLinearOp)
Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > getLinearOp() const
Simple aggregate class for a derivative object represented as a column-wise multi-vector or its trans...
MPDerivativeMultiVector getDerivativeMultiVector() const
Teuchos::RCP< Epetra_Operator > get_WPrec() const
void set_x_dotdot(const Teuchos::RCP< const Epetra_Vector > &x_dotdot)
Stochastic Galerkin solution vector polynomial.
Teuchos::RCP< const Teuchos::Polynomial< Epetra_Vector > > x_dot_poly_
std::string modelEvalDescription() const
SGDerivativeMultiVector dmv_
void set_W_sg(const sg_operator_t &W_sg)
Set stochastic Galerkin W operator polynomial.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > sg_vector_t
Short-hand for stochastic Galerkin vector type.
void set_DgDx_dot_sg(int j, const SGDerivative &DgDx_dot_j)
void set_DgDx_dot(int j, const Derivative &DgDx_dot_j)
DerivativeProperties get_DgDx_dot_sg_properties(int j) const
mp_const_vector_t get_x_dot_mp() const
Get multi-point time derivative vector.
sg_vector_t get_f_sg() const
Get stochastic Galerkin residual vector polynomial.
An approximate derivative (i.e. for a Jacobian)
supports_p_sg_t supports_p_mp_
void set_DgDx_dot_mp_properties(int j, const DerivativeProperties &properties)
void set_DgDp_properties(int j, int l, const DerivativeProperties &properties)
void set_DgDx_dotdot_mp(int j, const MPDerivative &DgDx_dotdot_j)
Teuchos::RCP< Epetra_MultiVector > getMultiVector() const
Evaluation(const Teuchos::RCP< ObjType > &obj)
void _set_DgDx_dotdot_mp_properties(int j, const DerivativeProperties &properties)
void set_DgDx_sg(int j, const SGDerivative &DgDx_j)
virtual Teuchos::ArrayView< const std::string > get_g_names(int j) const
Get the names of the response functions associated with response subvector j if available.
int get_stage_number() const
Teuchos::RCP< const Stokhos::EpetraVectorOrthogPoly > sg_const_vector_t
Short-hand for stochastic Galerkin vector type.
void set_p_mp(int l, const mp_const_vector_t &p_mp_l)
Set multi-point parameter vector.
void set_DgDx_mp(int j, const MPDerivative &DgDx_j)
void set_DgDp_sg(int j, int l, const SGDerivative &DgDp_sg_j_l)
void set_DfDp_properties(int l, const DerivativeProperties &properties)
void assert_l(int l) const
void setModelEvalDescription(const std::string &modelEvalDescription)
Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > lo_
void set_x_dotdot_sg(const sg_const_vector_t &x_dotdot_sg)
virtual Teuchos::RCP< const Epetra_Vector > get_x_lower_bounds() const
Teuchos::RCP< const Teuchos::Polynomial< Epetra_Vector > > get_x_dotdot_poly() const
virtual Teuchos::RCP< const Epetra_Vector > get_x_dotdot_init() const
Teuchos::RCP< Epetra_MultiVector > get_DfDp_mv(const int l, const ModelEvaluator::OutArgs &outArgs)
supports_t supports_DgDx_dotdot_mp_
Derivative(const DerivativeMultiVector &dmv)
void set_p(int l, const Teuchos::RCP< const Epetra_Vector > &p_l)
void _set_DgDx_mp_properties(int j, const DerivativeProperties &properties)
void set_DgDp_sg_properties(int j, int l, const DerivativeProperties &properties)
supports_t supports_DgDp_sg_
MPDerivative get_DgDx_mp(int j) const
Simple aggregate class for a derivative object represented as a column-wise multi-vector or its trans...
Teuchos::RCP< Epetra_MultiVector > getMultiVector() const
Derivative get_DgDp(int j, int l) const
void set_DgDx_dot_properties(int j, const DerivativeProperties &properties)
const Teuchos::Array< int > & getParamIndexes() const
SGDerivativeMultiVector()
Teuchos::RCP< Epetra_Operator > get_DfDp_op(const int l, const ModelEvaluator::OutArgs &outArgs)
sg_operator_t get_W_sg() const
Get stochastic Galerkin W operator polynomial.
void set_x_dot_sg(const sg_const_vector_t &x_dot_sg)
Set stochastic Galerkin time derivative vector polynomial.
Teuchos::RCP< Epetra_MultiVector > mv_
MPDerivative(const MPDerivativeMultiVector &dmv)
DerivativeSupport & plus(EDerivativeMultiVectorOrientation mvOrientation)
Teuchos::RCP< const Epetra_Vector > get_x_dot() const
Multi-point time second derivative vector.
SGDerivative get_DfDp_sg(int l) const
Preconditioner operator (approx Jacobian)
MPDerivative(const mp_multivector_t &mv, const EDerivativeMultiVectorOrientation orientation=DERIV_MV_BY_COL)
Teuchos::RCP< const Stokhos::ProductEpetraVector > mp_const_vector_t
bool supports(EInArgsMembers arg) const
EDerivativeMultiVectorOrientation getOrientation() const
void assert_supports(EInArgsMembers arg) const
deriv_properties_t DgDx_sg_properties_
void assert_j(int j) const
Teuchos::RCP< Epetra_MultiVector > get_DgDx_mv(const int j, const ModelEvaluator::OutArgs &outArgs, const ModelEvaluator::EDerivativeMultiVectorOrientation mvOrientation)
Teuchos::Array< int > paramIndexes_
void _set_WPrec_properties(const DerivativeProperties &WPrec_properties)
Teuchos::RCP< const Stokhos::OrthogPolyBasis< int, double > > get_sg_basis() const
void set_W(const Teuchos::RCP< Epetra_Operator > &W)
void _set_DgDx_dot_mp_properties(int j, const DerivativeProperties &properties)
virtual Teuchos::RCP< const Epetra_Vector > get_x_dot_init() const
virtual Teuchos::RCP< const Epetra_Map > get_p_map(int l) const
.
Derivative(const Teuchos::RCP< Epetra_Operator > &lo)
void _set_DgDx_dot_properties(int j, const DerivativeProperties &properties)
mp_operator_t getLinearOp() const
EDerivativeMultiVectorOrientation
mp_operator_t get_W_mp() const
Get multi-point W.
Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > mv_
Multi-point residual vector.
RCP< ObjType > & operator=(const RCP< ObjType > &r_ptr)
deriv_properties_t DgDx_dot_sg_properties_
void set_sg_expansion(const Teuchos::RCP< Stokhos::OrthogPolyExpansion< int, double, Stokhos::StandardStorage< int, double > > > &exp)
virtual Teuchos::RCP< Epetra_Operator > create_DgDp_op(int j, int l) const
Teuchos::Array< int > paramIndexes_
void setModelEvalDescription(const std::string &modelEvalDescription)
Derivative get_DgDx_dotdot(int j) const
void set_f_sg(const sg_vector_t &f_sg)
Set stochastic Galerkin residual vector polynomial.
DerivativeProperties(EDerivativeLinearity in_linearity, ERankStatus in_rank, bool in_supportsAdjoint)
SGDerivative(const Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > &lo)
mp_const_vector_t x_dot_mp_
Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > getMultiVector() const
Teuchos::RCP< const Stokhos::Quadrature< int, double > > sg_quad_
void set_DgDp(int j, int l, const Derivative &DgDp_j_l)
void set_W_mp(const mp_operator_t &W_sg)
Set multi-point W.
Teuchos::Array< Evaluation< Epetra_Vector > > g_t
void set_f_poly(const Teuchos::RCP< Teuchos::Polynomial< Epetra_Vector > > &f_poly)
Set residual vector Taylor polynomial.
void set_DgDx_dotdot_sg_properties(int j, const DerivativeProperties &properties)
virtual Teuchos::RCP< Epetra_Operator > create_DgDx_dotdot_op(int j) const
void _set_DgDx_dot_sg_properties(int j, const DerivativeProperties &properties)
Multi-point time derivative vector.
supports_t supports_DfDp_
EDerivativeMultiVectorOrientation getMultiVectorOrientation() const
Teuchos::Array< DerivativeProperties > deriv_properties_t
sg_vector_t get_g_sg(int j) const
Get stochastic Galerkin vector polynomial response.
DerivativeProperties get_DgDx_dotdot_mp_properties(int j) const
deriv_properties_t DgDx_properties_
DerivativeProperties get_DfDp_properties(int l) const
void set_x_dot_mp(const mp_const_vector_t &x_dot_mp)
Set multi-point time derivative vector.
bool isFailed() const
Return if the evaluation failed or not.
void set_DfDp_mp_properties(int l, const DerivativeProperties &properties)
void set_DgDx_dot_mp(int j, const MPDerivative &DgDx_dot_j)
Derivative(const Teuchos::RCP< Epetra_MultiVector > &mv, const EDerivativeMultiVectorOrientation orientation=DERIV_MV_BY_COL)
void set_x_poly(const Teuchos::RCP< const Teuchos::Polynomial< Epetra_Vector > > &x_poly)
DerivativeMultiVector(const Teuchos::RCP< Epetra_MultiVector > &mv, const EDerivativeMultiVectorOrientation orientation=DERIV_MV_BY_COL, const Teuchos::Array< int > ¶mIndexes=Teuchos::Array< int >())
virtual Teuchos::RCP< const Epetra_Vector > get_p_lower_bounds(int l) const
void changeOrientation(const EDerivativeMultiVectorOrientation orientation)
void reset(const Teuchos::RCP< ObjType > &obj, EEvalType evalType)
Simple aggregate struct that stores a preconditioner as an Epetra_Operator and a bool, about whether it is inverted or not.
EDerivativeMultiVectorOrientation getOrientation() const
void _set_DfDp_mp_properties(int l, const DerivativeProperties &properties)
Solution vector Taylor polynomial.
Stochastic Galerkin time derivative vector polynomial.
Teuchos::Array< Teuchos::RCP< const Epetra_Vector > > p_t
virtual double getInfBound() const
Return the value of an infinite bound.
void set_omega(double omega)
sg_const_vector_t get_p_sg(int l) const
Get stochastic Galerkin vector polynomial parameter.
supports_t supports_DgDx_dot_mp_
void _set_DgDp_properties(int j, int l, const DerivativeProperties &properties)
Stochastic Galerkin residual vector polynomial.
virtual double get_t_lower_bound() const
Stochastic Galerkin "W" operator polyomial.
deriv_properties_t DgDp_properties_
SGDerivative get_DgDx_dot_sg(int j) const
Teuchos::Array< sg_const_vector_t > p_sg_t
void set_alpha(double alpha)
Teuchos::RCP< Epetra_MultiVector > get_DgDx_dotdot_mv(const int j, const ModelEvaluator::OutArgs &outArgs, const ModelEvaluator::EDerivativeMultiVectorOrientation mvOrientation)
deriv_properties_t DgDx_dotdot_properties_
SGDerivative(const Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > &mv, const EDerivativeMultiVectorOrientation orientation=DERIV_MV_BY_COL)
Teuchos::RCP< const Epetra_Vector > x_
DerivativeSupport(EDerivativeLinearOp, EDerivativeMultiVectorOrientation mvOrientation)
virtual Teuchos::RCP< const Epetra_Vector > get_x_init() const
mp_vector_t get_f_mp() const
Get multi-point residual vector.
void set_DfDp_sg_properties(int l, const DerivativeProperties &properties)
Teuchos::Array< mp_const_vector_t > p_mp_t
virtual double get_t_init() const
DerivativeProperties get_DgDx_dotdot_properties(int j) const
virtual void evalModel(const InArgs &inArgs, const OutArgs &outArgs) const =0
EDerivativeMultiVectorOrientation getMultiVectorOrientation() const
supports_t supports_DgDx_
mp_multivector_t getMultiVector() const
supports_t supports_DgDx_dotdot_
supports_t supports_DgDp_mp_
Teuchos::RCP< Teuchos::Polynomial< Epetra_Vector > > f_poly_
Teuchos::RCP< Stokhos::ProductEpetraVector > mp_vector_t
void set_DfDp_mp(int l, const MPDerivative &DfDp_mp_l)
void set_DgDx_dotdot_mp_properties(int j, const DerivativeProperties &properties)
Teuchos::RCP< Stokhos::OrthogPolyExpansion< int, double, Stokhos::StandardStorage< int, double > > > get_sg_expansion() const
sg_deriv_t DgDx_dotdot_sg_
Time derivative vector Taylor polynomial.
Simple aggregate class for a derivative object represented as a column-wise multi-vector or its trans...
bool supportsTransMVByRow_
void _set_DgDx_properties(int j, const DerivativeProperties &properties)
virtual Teuchos::RCP< EpetraExt::ModelEvaluator::Preconditioner > create_WPrec() const
void set_DgDx_properties(int j, const DerivativeProperties &properties)
deriv_properties_t DgDp_sg_properties_
Teuchos::Array< Derivative > deriv_t
MPDerivative get_DgDp_mp(int j, int l) const
Teuchos::Array< bool > supports_g_sg_t
sg_const_vector_t x_dot_sg_
void set_DgDp_mp(int j, int l, const MPDerivative &DgDp_mp_j_l)
void set_DgDx_dot_sg_properties(int j, const DerivativeProperties &properties)
double get_step_size() const
Teuchos::RCP< Epetra_Operator > lo_
void _set_Np_Ng(int Np, int Ng)
SGDerivative(const SGDerivativeMultiVector &dmv)
DerivativeProperties get_DgDx_dot_mp_properties(int j) const
DerivativeProperties get_DfDp_mp_properties(int l) const
supports_g_sg_t supports_g_sg_
static const int NUM_E_OUT_ARGS_MEMBERS
DerivativeProperties get_DgDx_dot_properties(int j) const
void set_x_dot_poly(const Teuchos::RCP< const Teuchos::Polynomial< Epetra_Vector > > &x_dot_poly)
Set time derivative vector Taylor polynomial.
void _set_DgDp_sg_properties(int j, int l, const DerivativeProperties &properties)
void set_g_mp(int j, const mp_vector_t &g_mp_j)
Set multi-point response.
virtual Teuchos::RCP< Epetra_Operator > create_DfDp_op(int l) const
void set_x_mp(const mp_const_vector_t &x_mp)
Set multi-point solution vector.
Stochastic Galerkin basis.
supports_t supports_DfDp_mp_
SGDerivative get_DgDx_sg(int j) const
Stochastic Galerkin time second derivative vector polynomial.
const Teuchos::Array< int > & getParamIndexes() const
void set_sg_quadrature(const Teuchos::RCP< const Stokhos::Quadrature< int, double > > &quad)
void set_x_sg(const sg_const_vector_t &x_sg)
Set stochastic Galerkin solution vector polynomial.
std::string modelEvalDescription() const
void set_DgDx_dotdot_sg(int j, const SGDerivative &DgDx_dotdot_j)
DerivativeProperties W_properties_
void _set_DgDx_dotdot_properties(int j, const DerivativeProperties &properties)
MPDerivative get_DgDx_dot_mp(int j) const
void _setSupports(EInArgsMembers arg, bool supports)
bool supports_[NUM_E_OUT_ARGS_MEMBERS]
Teuchos::RCP< Epetra_Operator > WPrec_
void set_f_mp(const mp_vector_t &f_mp)
Set multi-point residual vector.
Evaluation< Epetra_Vector > get_f() const
void assert_supports(EOutArgsMembers arg) const
virtual Teuchos::RCP< Epetra_Operator > create_W() const
If supported, create a Epetra_Operator object for W to be evaluated.
sg_const_vector_t get_x_sg() const
Get stochastic Galerkin solution vector polynomial.
void _set_DgDx_sg_properties(int j, const DerivativeProperties &properties)
Teuchos::RCP< const Epetra_Vector > get_x() const
Set solution vector Taylor polynomial.
void _set_W_properties(const DerivativeProperties &W_properties)
void changeOrientation(const EDerivativeMultiVectorOrientation orientation)
mp_deriv_t DgDx_dotdot_mp_
sg_const_vector_t x_dotdot_sg_
deriv_properties_t DfDp_properties_
DerivativeSupport(EDerivativeMultiVectorOrientation mvOrientation1, EDerivativeMultiVectorOrientation mvOrientation2)
DerivativeSupport & plus(EDerivativeLinearOp)
DerivativeProperties get_DgDp_sg_properties(int j, int l) const
virtual OutArgs createOutArgs() const =0
void _set_DgDx_dotdot_sg_properties(int j, const DerivativeProperties &properties)
sg_const_vector_t get_x_dot_sg() const
Get stochastic Galerkin time derivative vector polynomial.
Teuchos::Array< SGDerivative > sg_deriv_t
void set_p_sg(int l, const sg_const_vector_t &p_sg_l)
Set stochastic Galerkin vector polynomial parameter.
Base interface for evaluating a stateless "model".
deriv_properties_t DgDx_dotdot_sg_properties_
void set_DgDx_dotdot_properties(int j, const DerivativeProperties &properties)
supports_p_sg_t supports_p_sg_
void _setModelEvalDescription(const std::string &modelEvalDescription)
Teuchos::RCP< const Teuchos::Polynomial< Epetra_Vector > > x_dotdot_poly_
void set_Np_Ng(int Np, int Ng)
#define TEUCHOS_TEST_FOR_EXCEPT(throw_exception_test)
void _set_DfDp_sg_properties(int l, const DerivativeProperties &properties)
std::string toString(const int &x)
Teuchos::Array< bool > supports_p_sg_t
mp_const_vector_t x_dotdot_mp_
#define TEUCHOS_UNREACHABLE_RETURN(dummyReturnVal)
void _setModelEvalDescription(const std::string &modelEvalDescription)
Teuchos::RCP< const Stokhos::ProductEpetraOperator > mp_const_operator_t
void set_g(int j, const Evaluation< Epetra_Vector > &g_j)
Set g(j) where 0 <= j && j < this->Ng().
bool supports_[NUM_E_IN_ARGS_MEMBERS]
void set_DgDx_mp_properties(int j, const DerivativeProperties &properties)
void set_DgDp_mp_properties(int j, int l, const DerivativeProperties &properties)
deriv_properties_t DgDx_mp_properties_
Evaluation(const Teuchos::RCP< ObjType > &obj, EEvalType evalType)
sg_const_vector_t get_x_dotdot_sg() const