46 #include "Epetra_SerialComm.h" 47 #include "Epetra_CrsMatrix.h" 55 ,
const bool showGetInvalidArg
57 :d_(d),showGetInvalidArg_(showGetInvalidArg)
66 x0_ =
rcp(
new Epetra_Vector(*
map_x_)); (*x0_)[0] = x00; (*x0_)[1] = x01;
67 p_ =
rcp(
new Epetra_Vector(*
map_x_)); (*p_)[0] = p0; (*p_)[1] = p1;
72 int indices[nx] = { 0, 1 };
73 for(
int i = 0; i < nx; ++i )
74 W_graph_->InsertGlobalIndices(i,nx,indices);
137 using Teuchos::rcp_dynamic_cast;
141 const Epetra_Vector &x = *inArgs.
get_x();
145 Epetra_Vector *f_out = outArgs.
get_f().
get();
146 Epetra_Operator *W_out = outArgs.
get_W().
get();
148 Epetra_Vector *g_out = outArgs.
get_g(0).
get();
153 const Epetra_Vector &p = *
p_;
155 Epetra_Vector &
f = *f_out;
156 f[0] = x[0] + x[1]*x[1] - p[0];
157 f[1] =
d_ * ( x[0]*x[0] - x[1] - p[1] );
160 Epetra_CrsMatrix &DfDx =
dyn_cast<Epetra_CrsMatrix>(*W_out);
171 values[0] = 1.0; indexes[0] = 0;
172 values[1] = 2.0*x[1]; indexes[1] = 1;
173 DfDx.SumIntoGlobalValues( 0, 2, values, indexes );
175 values[0] = 2.0*
d_*x[0]; indexes[0] = 0;
176 values[1] = -
d_; indexes[1] = 1;
177 DfDx.SumIntoGlobalValues( 1, 2, values, indexes );
Teuchos::RCP< const Epetra_Map > map_x_
void set_W_properties(const DerivativeProperties &properties)
void setSupports(EOutArgsMembers arg, bool supports=true)
Teuchos::RCP< const Epetra_Map > get_x_map() const
Teuchos::RCP< Epetra_CrsGraph > W_graph_
void setSupports(EInArgsMembers arg, bool supports=true)
Teuchos::RCP< const Epetra_Vector > get_x_init() const
virtual std::string description() const
EpetraModelEval2DSim(const double d=10.0, const double p0=2.0, const double p1=0.0, const double x00=1.0, const double x01=1.0, const bool showGetInvalidArg=false)
Teuchos::RCP< const Epetra_Vector > get_x() const
Set solution vector Taylor polynomial.
void setModelEvalDescription(const std::string &modelEvalDescription)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Teuchos::RCP< Epetra_Operator > get_W() const
T_To & dyn_cast(T_From &from)
Teuchos::RCP< const Epetra_Comm > epetra_comm_
Teuchos::RCP< Epetra_Vector > x0_
OutArgs createOutArgs() const
void setModelEvalDescription(const std::string &modelEvalDescription)
InArgs createInArgs() const
Teuchos::RCP< Epetra_Vector > p_
Teuchos::RCP< const Epetra_Map > get_f_map() const
Evaluation< Epetra_Vector > get_g(int j) const
Get g(j) where 0 <= j && j < this->Ng().
void evalModel(const InArgs &inArgs, const OutArgs &outArgs) const
Teuchos::RCP< Epetra_Operator > create_W() const
Evaluation< Epetra_Vector > get_f() const