39 InArgs me_inargs =
me->createInArgs();
40 num_p = me_inargs.Np();
41 if (me_inargs.supports(IN_ARG_x_dot))
43 if (me_inargs.supports(IN_ARG_x))
46 for (
int i=0; i<
num_p; i++)
50 OutArgs me_outargs =
me->createOutArgs();
51 num_g = me_outargs.Ng();
54 if (me_outargs.supports(OUT_ARG_f))
58 if (me_outargs.supports(OUT_ARG_W))
63 for (
int i=0; i<
num_p; i++)
64 if (me_outargs.supports(OUT_ARG_DfDp,i).supports(DERIV_MV_BY_COL))
65 dfdp_qp[i] = EpetraExt::ModelEvaluator::Derivative(
68 me->get_p_map(i)->NumGlobalElements())));
69 else if (me_outargs.supports(OUT_ARG_DfDp,i).supports(DERIV_TRANS_MV_BY_ROW))
70 dfdp_qp[i] = EpetraExt::ModelEvaluator::Derivative(
73 me->get_f_map()->NumGlobalElements())));
74 else if (me_outargs.supports(OUT_ARG_DfDp,i).supports(DERIV_LINEAR_OP))
75 dfdp_qp[i] = EpetraExt::ModelEvaluator::Derivative(
76 me->create_DfDp_op(i));
82 for (
int i=0; i<
num_g; i++) {
89 if (me_outargs.supports(OUT_ARG_DgDx, i).supports(DERIV_TRANS_MV_BY_ROW))
90 dgdx_qp[i] = EpetraExt::ModelEvaluator::Derivative(
93 me->get_g_map(i)->NumGlobalElements())));
94 else if (me_outargs.supports(OUT_ARG_DgDx, i).supports(DERIV_MV_BY_COL))
95 dgdx_qp[i] = EpetraExt::ModelEvaluator::Derivative(
98 me->get_x_map()->NumGlobalElements())));
99 else if (me_outargs.supports(OUT_ARG_DgDx, i).supports(DERIV_LINEAR_OP))
100 dgdx_qp[i] = EpetraExt::ModelEvaluator::Derivative(
101 me->create_DgDx_op(i));
104 if (me_outargs.supports(OUT_ARG_DgDx_dot, i).supports(DERIV_TRANS_MV_BY_ROW))
105 dgdx_dot_qp[i] = EpetraExt::ModelEvaluator::Derivative(
108 me->get_g_map(i)->NumGlobalElements())));
109 else if (me_outargs.supports(OUT_ARG_DgDx_dot, i).supports(DERIV_MV_BY_COL))
110 dgdx_dot_qp[i] = EpetraExt::ModelEvaluator::Derivative(
113 me->get_x_map()->NumGlobalElements())));
114 else if (me_outargs.supports(OUT_ARG_DgDx_dot, i).supports(DERIV_LINEAR_OP))
115 dgdx_dot_qp[i] = EpetraExt::ModelEvaluator::Derivative(
116 me->create_DgDx_dot_op(i));
121 if (me_outargs.supports(OUT_ARG_DgDp, i,
j).supports(DERIV_TRANS_MV_BY_ROW))
122 dgdp_qp[i][
j] = EpetraExt::ModelEvaluator::Derivative(
125 me->get_g_map(i)->NumGlobalElements())));
126 else if (me_outargs.supports(OUT_ARG_DgDp, i,
j).supports(DERIV_MV_BY_COL))
127 dgdp_qp[i][
j] = EpetraExt::ModelEvaluator::Derivative(
130 me->get_p_map(
j)->NumGlobalElements())));
131 else if (me_outargs.supports(OUT_ARG_DgDp, i,
j).supports(DERIV_LINEAR_OP))
132 dgdp_qp[i][
j] = EpetraExt::ModelEvaluator::Derivative(
133 me->create_DgDp_op(i,
j));
143 return me->get_x_map();
150 return me->get_f_map();
157 return me->get_p_map(l);
164 return me->get_g_map(l);
171 return me->get_p_names(l);
178 return me->get_x_init();
185 return me->get_p_init(l);
192 return me->create_W();
195 EpetraExt::ModelEvaluator::InArgs
200 InArgs me_inargs = me->createInArgs();
202 inArgs.setModelEvalDescription(this->description());
203 inArgs.set_Np(num_p);
204 inArgs.setSupports(IN_ARG_x_dot, me_inargs.supports(IN_ARG_x_dot));
205 inArgs.setSupports(IN_ARG_x, me_inargs.supports(IN_ARG_x));
206 inArgs.setSupports(IN_ARG_t, me_inargs.supports(IN_ARG_t));
207 inArgs.setSupports(IN_ARG_alpha, me_inargs.supports(IN_ARG_alpha));
208 inArgs.setSupports(IN_ARG_beta, me_inargs.supports(IN_ARG_beta));
210 for (
int i=0; i<num_p; i++)
211 inArgs.setSupports(IN_ARG_p_sg, i,
true);
212 inArgs.setSupports(IN_ARG_x_sg, me_inargs.supports(IN_ARG_x));
213 inArgs.setSupports(IN_ARG_x_dot_sg, me_inargs.supports(IN_ARG_x_dot));
214 inArgs.setSupports(IN_ARG_sg_basis,
true);
215 inArgs.setSupports(IN_ARG_sg_quadrature,
true);
220 EpetraExt::ModelEvaluator::OutArgs
224 OutArgsSetup outArgs;
225 OutArgs me_outargs = me->createOutArgs();
227 outArgs.setModelEvalDescription(this->description());
228 outArgs.set_Np_Ng(num_p, num_g);
229 outArgs.setSupports(OUT_ARG_f, me_outargs.supports(OUT_ARG_f));
230 outArgs.setSupports(OUT_ARG_W, me_outargs.supports(OUT_ARG_W));
231 for (
int j=0;
j<num_p;
j++)
232 outArgs.setSupports(OUT_ARG_DfDp,
j,
233 me_outargs.supports(OUT_ARG_DfDp,
j));
234 for (
int i=0; i<num_g; i++) {
235 outArgs.setSupports(OUT_ARG_DgDx, i,
236 me_outargs.supports(OUT_ARG_DgDx, i));
237 outArgs.setSupports(OUT_ARG_DgDx_dot, i,
238 me_outargs.supports(OUT_ARG_DgDx_dot, i));
239 for (
int j=0;
j<num_p;
j++)
240 outArgs.setSupports(OUT_ARG_DgDp, i,
j,
241 me_outargs.supports(OUT_ARG_DgDp, i,
j));
244 outArgs.setSupports(OUT_ARG_f_sg, me_outargs.supports(OUT_ARG_f));
245 if (me_outargs.supports(OUT_ARG_W)) {
246 outArgs.set_W_properties(me_outargs.get_W_properties());
247 outArgs.setSupports(OUT_ARG_W_sg,
true);
249 for (
int j=0;
j<num_p;
j++) {
250 outArgs.setSupports(OUT_ARG_DfDp_sg,
j,
251 me_outargs.supports(OUT_ARG_DfDp,
j));
253 for (
int i=0; i<num_g; i++) {
254 outArgs.setSupports(OUT_ARG_g_sg, i,
true);
255 outArgs.setSupports(OUT_ARG_DgDx_sg, i,
256 me_outargs.supports(OUT_ARG_DgDx, i));
257 outArgs.setSupports(OUT_ARG_DgDx_dot_sg, i,
258 me_outargs.supports(OUT_ARG_DgDx_dot, i));
259 for (
int j=0;
j<num_p;
j++)
260 outArgs.setSupports(OUT_ARG_DgDp_sg, i,
j,
261 me_outargs.supports(OUT_ARG_DgDp, i,
j));
269 evalModel(
const InArgs& inArgs,
const OutArgs& outArgs)
const
272 InArgs me_inargs = me->createInArgs();
273 if (me_inargs.supports(IN_ARG_x))
274 me_inargs.set_x(inArgs.get_x());
275 if (me_inargs.supports(IN_ARG_x_dot))
276 me_inargs.set_x_dot(inArgs.get_x_dot());
277 if (me_inargs.supports(IN_ARG_alpha))
278 me_inargs.set_alpha(inArgs.get_alpha());
279 if (me_inargs.supports(IN_ARG_beta))
280 me_inargs.set_beta(inArgs.get_beta());
281 if (me_inargs.supports(IN_ARG_t))
282 me_inargs.set_t(inArgs.get_t());
283 for (
int i=0; i<num_p; i++)
284 me_inargs.set_p(i, inArgs.get_p(i));
287 OutArgs me_outargs = me->createOutArgs();
288 if (me_outargs.supports(OUT_ARG_f))
289 me_outargs.set_f(outArgs.get_f());
290 if (me_outargs.supports(OUT_ARG_W))
291 me_outargs.set_W(outArgs.get_W());
292 for (
int j=0;
j<num_p;
j++)
293 if (!outArgs.supports(OUT_ARG_DfDp,
j).none())
294 me_outargs.set_DfDp(
j, outArgs.get_DfDp(
j));
295 for (
int i=0; i<num_g; i++) {
296 me_outargs.set_g(i, outArgs.get_g(i));
297 if (!outArgs.supports(OUT_ARG_DgDx, i).none())
298 me_outargs.set_DgDx(i, outArgs.get_DgDx(i));
299 if (!outArgs.supports(OUT_ARG_DgDx_dot, i).none())
300 me_outargs.set_DgDx(i, outArgs.get_DgDx_dot(i));
301 for (
int j=0;
j<num_p;
j++)
302 if (!outArgs.supports(OUT_ARG_DgDp, i,
j).none())
303 me_outargs.set_DgDp(i,
j, outArgs.get_DgDp(i,
j));
306 bool do_quad =
false;
307 InArgs::sg_const_vector_t x_sg;
308 InArgs::sg_const_vector_t x_dot_sg;
310 OutArgs::sg_vector_t f_sg;
311 OutArgs::sg_operator_t W_sg;
319 "Error! Stokhos::SGQuadModelEvaluator::evalModel(): " <<
320 "SG basis inArg cannot be null!");
323 "Error! Stokhos::SGQuadModelEvaluator::evalModel(): " <<
324 "SG quadrature inArg cannot be null!");
326 inArgs.get_sg_basis();
328 inArgs.get_sg_quadrature();
329 if (inArgs.supports(IN_ARG_x_sg)) {
330 x_sg = inArgs.get_x_sg();
335 if (inArgs.supports(IN_ARG_x_dot_sg)) {
336 x_dot_sg = inArgs.get_x_dot_sg();
341 for (
int i=0; i<num_p; i++) {
342 p_sg[i] = inArgs.get_p_sg(i);
347 if (outArgs.supports(OUT_ARG_f_sg)) {
348 f_sg = outArgs.get_f_sg();
352 if (outArgs.supports(OUT_ARG_W_sg)) {
353 W_sg = outArgs.get_W_sg();
357 for (
int i=0; i<num_p; i++) {
358 if (!outArgs.supports(OUT_ARG_DfDp_sg, i).none()) {
359 dfdp_sg[i] = outArgs.get_DfDp_sg(i);
361 dfdp_sg[i].getMultiVector()->init(0.0);
363 dfdp_sg[i].getLinearOp()->init(0.0);
367 for (
int i=0; i<num_g; i++) {
368 g_sg[i] = outArgs.get_g_sg(i);
372 if (!outArgs.supports(OUT_ARG_DgDx_sg, i).none()) {
373 dgdx_sg[i] = outArgs.get_DgDx_sg(i);
375 dgdx_sg[i].getMultiVector()->init(0.0);
377 dgdx_sg[i].getLinearOp()->init(0.0);
380 if (!outArgs.supports(OUT_ARG_DgDx_dot_sg, i).none()) {
381 dgdx_dot_sg[i] = outArgs.get_DgDx_dot_sg(i);
383 dgdx_dot_sg[i].getMultiVector()->init(0.0);
385 dgdx_dot_sg[i].getLinearOp()->init(0.0);
389 for (
int j=0;
j<num_p;
j++) {
390 if (!outArgs.supports(OUT_ARG_DgDp_sg, i,
j).none()) {
391 dgdp_sg[i][
j] = outArgs.get_DgDp_sg(i,
j);
393 dgdp_sg[i][
j].getMultiVector()->init(0.0);
395 dgdp_sg[i][
j].getLinearOp()->init(0.0);
403 quad->getQuadPoints();
405 quad->getQuadWeights();
407 quad->getBasisAtQuadPoints();
411 for (
int qp=0; qp<quad_points.
size(); qp++) {
416 if (quad_weights[qp] == 0.0)
420 #ifdef STOKHOS_TEUCHOS_TIME_MONITOR
427 #ifdef STOKHOS_TEUCHOS_TIME_MONITOR
430 x_sg->evaluate(quad_values[qp], *x_qp);
431 me_inargs.set_x(x_qp);
434 #ifdef STOKHOS_TEUCHOS_TIME_MONITOR
437 x_dot_sg->evaluate(quad_values[qp], *x_dot_qp);
438 me_inargs.set_x_dot(x_qp);
440 for (
int i=0; i<num_p; i++) {
442 #ifdef STOKHOS_TEUCHOS_TIME_MONITOR
445 p_sg[i]->evaluate(quad_values[qp], *(p_qp[i]));
446 me_inargs.set_p(i, p_qp[i]);
450 me_outargs.set_f(f_qp);
452 me_outargs.set_W(W_qp);
453 for (
int i=0; i<num_p; i++) {
454 if (!dfdp_sg[i].isEmpty())
455 me_outargs.set_DfDp(i, dfdp_qp[i]);
457 for (
int i=0; i<num_g; i++) {
459 me_outargs.set_g(i, g_qp[i]);
460 if (!dgdx_dot_sg[i].isEmpty())
461 me_outargs.set_DgDx_dot(i, dgdx_dot_qp[i]);
462 if (!dgdx_sg[i].isEmpty())
463 me_outargs.set_DgDx(i, dgdx_qp[i]);
464 for (
int j=0;
j<num_p;
j++)
465 if (!dgdp_sg[i][
j].isEmpty())
466 me_outargs.set_DgDp(i,
j, dgdp_qp[i][
j]);
472 #ifdef STOKHOS_TEUCHOS_TIME_MONITOR
477 me->evalModel(me_inargs, me_outargs);
482 #ifdef STOKHOS_TEUCHOS_TIME_MONITOR
484 "Stokhos: SGQuadModelEvaluator -- Polynomial Integration", Integration);
489 #ifdef STOKHOS_TEUCHOS_TIME_MONITOR
492 f_sg->sumIntoAllTerms(quad_weights[qp], quad_values[qp], basis_norms,
496 #ifdef STOKHOS_TEUCHOS_TIME_MONITOR
499 W_sg->sumIntoAllTerms(quad_weights[qp], quad_values[qp], basis_norms,
502 for (
int j=0;
j<num_p;
j++) {
503 if (!dfdp_sg[
j].isEmpty()) {
504 #ifdef STOKHOS_TEUCHOS_TIME_MONITOR
506 "Stokhos: SGQuadModelEvaluator -- df/dp Integration");
509 dfdp_sg[
j].getMultiVector()->sumIntoAllTerms(
510 quad_weights[qp], quad_values[qp], basis_norms,
511 *(dfdp_qp[
j].getMultiVector()));
514 dfdp_sg[
j].getLinearOp()->sumIntoAllTerms(
515 quad_weights[qp], quad_values[qp], basis_norms,
516 *(dfdp_qp[
j].getLinearOp()));
520 for (
int i=0; i<num_g; i++) {
522 #ifdef STOKHOS_TEUCHOS_TIME_MONITOR
525 g_sg[i]->sumIntoAllTerms(quad_weights[qp], quad_values[qp],
526 basis_norms, *g_qp[i]);
528 if (!dgdx_dot_sg[i].isEmpty()) {
529 #ifdef STOKHOS_TEUCHOS_TIME_MONITOR
531 "Stokhos: SGQuadModelEvaluator -- dg/dx_dot Integration");
534 dgdx_dot_sg[i].getMultiVector()->sumIntoAllTerms(
535 quad_weights[qp], quad_values[qp], basis_norms,
536 *(dgdx_dot_qp[i].getMultiVector()));
539 dgdx_dot_sg[i].getLinearOp()->sumIntoAllTerms(
540 quad_weights[qp], quad_values[qp], basis_norms,
541 *(dgdx_dot_qp[i].getLinearOp()));
544 if (!dgdx_sg[i].isEmpty()) {
545 #ifdef STOKHOS_TEUCHOS_TIME_MONITOR
547 "Stokhos: SGQuadModelEvaluator -- dg/dx Integration");
550 dgdx_sg[i].getMultiVector()->sumIntoAllTerms(
551 quad_weights[qp], quad_values[qp], basis_norms,
552 *(dgdx_qp[i].getMultiVector()));
555 dgdx_sg[i].getLinearOp()->sumIntoAllTerms(
556 quad_weights[qp], quad_values[qp], basis_norms,
557 *(dgdx_qp[i].getLinearOp()));
560 for (
int j=0;
j<num_p;
j++) {
561 if (!dgdp_sg[i][
j].isEmpty()) {
562 #ifdef STOKHOS_TEUCHOS_TIME_MONITOR
564 "Stokhos: SGQuadModelEvaluator -- dg/dp Integration");
567 dgdp_sg[i][
j].getMultiVector()->sumIntoAllTerms(
568 quad_weights[qp], quad_values[qp], basis_norms,
569 *(dgdp_qp[i][
j].getMultiVector()));
572 dgdp_sg[i][
j].getLinearOp()->sumIntoAllTerms(
573 quad_weights[qp], quad_values[qp], basis_norms,
574 *(dgdp_qp[i][
j].getLinearOp()));
585 me->evalModel(me_inargs, me_outargs);
int num_p
Number of parameter vectors.
#define TEUCHOS_FUNC_TIME_MONITOR(FUNCNAME)
Teuchos::RCP< const Epetra_Map > get_x_map() const
Return solution vector map.
int num_g
Number of response vectors.
Teuchos::RCP< const Epetra_Vector > get_x_init() const
Return initial solution.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Teuchos::Array< EpetraExt::ModelEvaluator::Derivative > dgdx_dot_qp
Response derivative.
Teuchos::Array< EpetraExt::ModelEvaluator::Derivative > dgdx_qp
Response derivative.
Teuchos::RCP< const Teuchos::Array< std::string > > get_p_names(int l) const
Return array of parameter names.
InArgs createInArgs() const
Create InArgs.
Teuchos::Array< Teuchos::RCP< Epetra_Vector > > g_qp
Response vectors.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Teuchos::Array< Teuchos::RCP< Epetra_Vector > > p_qp
Parameter vectors.
Teuchos::RCP< const Epetra_Map > get_p_map(int l) const
Return parameter vector map.
Teuchos::Array< EpetraExt::ModelEvaluator::Derivative > dfdp_qp
Residual derivatives.
Teuchos::RCP< Epetra_Vector > x_dot_qp
Time derivative vector.
void resize(size_type new_size, const value_type &x=value_type())
Teuchos::RCP< Epetra_Operator > create_W() const
Create W = alpha*M + beta*J matrix.
Teuchos::RCP< Epetra_Vector > f_qp
Residual vector.
Teuchos::RCP< Epetra_Vector > x_qp
Solution vector.
Teuchos::RCP< Epetra_Operator > W_qp
W operator.
Teuchos::Array< Teuchos::Array< EpetraExt::ModelEvaluator::Derivative > > dgdp_qp
Response sensitivities.
Teuchos::RCP< EpetraExt::ModelEvaluator > me
Underlying model evaluator.
Teuchos::RCP< const Epetra_Vector > get_p_init(int l) const
Return initial parameters.
Teuchos::RCP< const Epetra_Map > get_f_map() const
Return residual vector map.
OutArgs createOutArgs() const
Create OutArgs.
SGQuadModelEvaluator(const Teuchos::RCP< EpetraExt::ModelEvaluator > &me)
void evalModel(const InArgs &inArgs, const OutArgs &outArgs) const
Evaluate model on InArgs.
#define TEUCHOS_FUNC_TIME_MONITOR_DIFF(FUNCNAME, DIFF)
Teuchos::RCP< const Epetra_Map > get_g_map(int l) const
Return observation vector map.