411 const OutArgs& outArgs)
const
414 Teuchos::RCP<const Epetra_Vector> x;
415 if (inArgs.supports(IN_ARG_x)) {
417 if (x != Teuchos::null)
420 Teuchos::RCP<const Epetra_Vector> x_dot;
421 if (inArgs.supports(IN_ARG_x_dot))
422 x_dot = inArgs.get_x_dot();
425 EpetraExt::ModelEvaluator::Evaluation<Epetra_Vector> f_out;
426 if (outArgs.supports(OUT_ARG_f))
427 f_out = outArgs.get_f();
428 Teuchos::RCP<Epetra_Operator> W_out;
429 if (outArgs.supports(OUT_ARG_W))
430 W_out = outArgs.get_W();
433 InArgs me_inargs = me->createInArgs();
434 if (x != Teuchos::null) {
435 Teuchos::RCP<Epetra_Vector> overlapped_x
436 = Teuchos::rcp(
new Epetra_Vector(*interlace_overlapped_x_map));
437 overlapped_x->Import(*x,*interlace_overlapped_x_importer,
Insert);
442 copyToPolyOrthogVector(*overlapped_x,*x_sg_blocks);
443 me_inargs.set_x_sg(x_sg_blocks);
445 if (x_dot != Teuchos::null) {
446 Teuchos::RCP<Epetra_Vector> overlapped_x_dot
447 = Teuchos::rcp(
new Epetra_Vector(*interlace_overlapped_x_map));
448 overlapped_x_dot->Import(*x_dot,*interlace_overlapped_x_importer,
Insert);
453 copyToPolyOrthogVector(*overlapped_x_dot,*x_dot_sg_blocks);
454 me_inargs.set_x_dot_sg(x_dot_sg_blocks);
456 if (me_inargs.supports(IN_ARG_alpha))
457 me_inargs.set_alpha(inArgs.get_alpha());
458 if (me_inargs.supports(IN_ARG_beta))
459 me_inargs.set_beta(inArgs.get_beta());
460 if (me_inargs.supports(IN_ARG_t))
461 me_inargs.set_t(inArgs.get_t());
462 if (me_inargs.supports(IN_ARG_sg_basis)) {
463 if (inArgs.get_sg_basis() != Teuchos::null)
464 me_inargs.set_sg_basis(inArgs.get_sg_basis());
466 me_inargs.set_sg_basis(sg_basis);
468 if (me_inargs.supports(IN_ARG_sg_quadrature)) {
469 if (inArgs.get_sg_quadrature() != Teuchos::null)
470 me_inargs.set_sg_quadrature(inArgs.get_sg_quadrature());
472 me_inargs.set_sg_quadrature(sg_quad);
474 if (me_inargs.supports(IN_ARG_sg_expansion)) {
475 if (inArgs.get_sg_expansion() != Teuchos::null)
476 me_inargs.set_sg_expansion(inArgs.get_sg_expansion());
478 me_inargs.set_sg_expansion(sg_exp);
482 for (
int i=0; i<num_p; i++)
483 me_inargs.set_p(i, inArgs.get_p(i));
484 for (
int i=0; i<num_p_sg; i++) {
485 Teuchos::RCP<const Epetra_Vector> p = inArgs.get_p(i+num_p);
489 if (p == Teuchos::null)
490 p = sg_p_init[i]->getBlockVector();
493 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> p_sg =
494 create_p_sg(sg_p_index_map[i],
View, p.get());
495 me_inargs.set_p_sg(sg_p_index_map[i], p_sg);
499 OutArgs me_outargs = me->createOutArgs();
502 if (f_out != Teuchos::null) {
503 me_outargs.set_f_sg(f_sg_blocks);
505 me_outargs.set_W_sg(W_sg_blocks);
509 if (W_out != Teuchos::null && !eval_W_with_f )
510 me_outargs.set_W_sg(W_sg_blocks);
513 for (
int i=0; i<outArgs.Np(); i++) {
514 if (!outArgs.supports(OUT_ARG_DfDp, i).none()) {
515 Derivative dfdp = outArgs.get_DfDp(i);
516 if (dfdp.getMultiVector() != Teuchos::null) {
517 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dfdp_sg;
518 if (dfdp.getMultiVectorOrientation() == DERIV_MV_BY_COL)
521 sg_basis, overlapped_stoch_row_map,
522 me->get_f_map(), interlace_overlapped_f_map, sg_comm,
523 me->get_p_map(i)->NumMyElements()));
524 else if (dfdp.getMultiVectorOrientation() == DERIV_TRANS_MV_BY_ROW)
527 sg_basis, overlapped_stoch_row_map,
528 me->get_p_map(i), sg_comm,
529 me->get_f_map()->NumMyElements()));
530 me_outargs.set_DfDp_sg(i,
531 SGDerivative(dfdp_sg,
532 dfdp.getMultiVectorOrientation()));
534 TEUCHOS_TEST_FOR_EXCEPTION(dfdp.getLinearOp() != Teuchos::null, std::logic_error,
535 "Error! Stokhos::SGModelEvaluator_Interlaced::evalModel " <<
536 "cannot handle operator form of df/dp!");
541 for (
int i=0; i<num_g_sg; i++) {
542 int ii = sg_g_index_map[i];
545 Teuchos::RCP<Epetra_Vector>
g = outArgs.get_g(i);
546 if (
g != Teuchos::null) {
547 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> g_sg =
548 create_g_sg(sg_g_index_map[i],
View,
g.get());
549 me_outargs.set_g_sg(i, g_sg);
553 if (outArgs.supports(OUT_ARG_DgDx_dot, i).supports(DERIV_LINEAR_OP)) {
554 Derivative dgdx_dot = outArgs.get_DgDx_dot(i);
555 if (dgdx_dot.getLinearOp() != Teuchos::null) {
556 Teuchos::RCP<Stokhos::SGOperator> op =
557 Teuchos::rcp_dynamic_cast<Stokhos::SGOperator>(
558 dgdx_dot.getLinearOp(),
true);
559 Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > sg_blocks =
560 op->getSGPolynomial();
561 if (me_outargs.supports(OUT_ARG_DgDx, ii).supports(DERIV_LINEAR_OP))
562 me_outargs.set_DgDx_dot_sg(ii, sg_blocks);
564 for (
unsigned int k=0; k<num_sg_blocks; k++) {
565 Teuchos::RCP<Epetra_MultiVector> mv =
566 Teuchos::rcp_dynamic_cast<Stokhos::EpetraMultiVectorOperator>(
567 sg_blocks->getCoeffPtr(k),
true)->getMultiVector();
568 dgdx_dot_sg_blocks[i]->setCoeffPtr(k, mv);
570 if (me_outargs.supports(OUT_ARG_DgDx_dot_sg, ii).supports(DERIV_MV_BY_COL))
571 me_outargs.set_DgDx_dot_sg(ii, SGDerivative(dgdx_dot_sg_blocks[i],
574 me_outargs.set_DgDx_dot_sg(ii, SGDerivative(dgdx_dot_sg_blocks[i],
575 DERIV_TRANS_MV_BY_ROW));
578 TEUCHOS_TEST_FOR_EXCEPTION(dgdx_dot.getLinearOp() == Teuchos::null &&
579 dgdx_dot.isEmpty() ==
false,
581 "Error! Stokhos::SGModelEvaluator_Interlaced::evalModel: " <<
582 "Operator form of dg/dxdot is required!");
586 if (outArgs.supports(OUT_ARG_DgDx, i).supports(DERIV_LINEAR_OP)) {
587 Derivative dgdx = outArgs.get_DgDx(i);
588 if (dgdx.getLinearOp() != Teuchos::null) {
589 Teuchos::RCP<Stokhos::SGOperator> op =
590 Teuchos::rcp_dynamic_cast<Stokhos::SGOperator>(
591 dgdx.getLinearOp(),
true);
592 Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > sg_blocks =
593 op->getSGPolynomial();
594 if (me_outargs.supports(OUT_ARG_DgDx, ii).supports(DERIV_LINEAR_OP))
595 me_outargs.set_DgDx_sg(i, sg_blocks);
597 for (
unsigned int k=0; k<num_sg_blocks; k++) {
598 Teuchos::RCP<Epetra_MultiVector> mv =
599 Teuchos::rcp_dynamic_cast<Stokhos::EpetraMultiVectorOperator>(
600 sg_blocks->getCoeffPtr(k),
true)->getMultiVector();
601 dgdx_sg_blocks[i]->setCoeffPtr(k, mv);
603 if (me_outargs.supports(OUT_ARG_DgDx_sg, ii).supports(DERIV_MV_BY_COL))
604 me_outargs.set_DgDx_sg(ii, SGDerivative(dgdx_sg_blocks[i],
607 me_outargs.set_DgDx_sg(ii, SGDerivative(dgdx_sg_blocks[i],
608 DERIV_TRANS_MV_BY_ROW));
611 TEUCHOS_TEST_FOR_EXCEPTION(dgdx.getLinearOp() == Teuchos::null &&
612 dgdx.isEmpty() ==
false,
614 "Error! Stokhos::SGModelEvaluator_Interlaced::evalModel: " <<
615 "Operator form of dg/dxdot is required!");
620 for (
int j=0;
j<num_p;
j++) {
621 if (!outArgs.supports(OUT_ARG_DgDp, i,
j).none()) {
622 Derivative dgdp = outArgs.get_DgDp(i,
j);
623 if (dgdp.getMultiVector() != Teuchos::null) {
624 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dgdp_sg;
625 if (dgdp.getMultiVectorOrientation() == DERIV_MV_BY_COL)
628 sg_basis, overlapped_stoch_row_map,
629 me->get_g_map(ii), sg_g_map[i], sg_comm,
630 View, *(dgdp.getMultiVector())));
631 else if (dgdp.getMultiVectorOrientation() == DERIV_TRANS_MV_BY_ROW) {
632 Teuchos::RCP<const Epetra_BlockMap> product_map =
633 Teuchos::rcp(&(dgdp.getMultiVector()->Map()),
false);
636 sg_basis, overlapped_stoch_row_map,
637 me->get_p_map(
j), product_map, sg_comm,
638 View, *(dgdp.getMultiVector())));
640 me_outargs.set_DgDp_sg(ii,
j,
641 SGDerivative(dgdp_sg,
642 dgdp.getMultiVectorOrientation()));
644 TEUCHOS_TEST_FOR_EXCEPTION(dgdp.getLinearOp() != Teuchos::null,
646 "Error! Stokhos::SGModelEvaluator_Interlaced::evalModel " <<
647 "cannot handle operator form of dg/dp!");
654 me->evalModel(me_inargs, me_outargs);
657 if ((W_out != Teuchos::null || (eval_W_with_f && f_out != Teuchos::null)) ) {
659 Teuchos::RCP<Epetra_Operator> W;
660 if (W_out != Teuchos::null)
664 Teuchos::RCP<Stokhos::SGOperator> W_sg =
665 Teuchos::rcp_dynamic_cast<Stokhos::SGOperator>(W,
true);
666 W_sg->setupOperator(W_sg_blocks);
670 if (f_out!=Teuchos::null){
672 for (
int i=0; i<sg_basis->size(); i++)
673 (*f_sg_blocks)[i].Scale(sg_basis->norm_squared(i));
675 Teuchos::RCP<Epetra_Vector> overlapped_f
676 = Teuchos::rcp(
new Epetra_Vector(*interlace_overlapped_f_map));
677 copyToInterlacedVector(*f_sg_blocks,*overlapped_f);
678 f_out->Export(*overlapped_f,*interlace_overlapped_f_exporter,
Insert);
682 for (
int i=0; i<num_p; i++) {
683 if (!outArgs.supports(OUT_ARG_DfDp, i).none()) {
684 Derivative dfdp = outArgs.get_DfDp(i);
685 SGDerivative dfdp_sg = me_outargs.get_DfDp_sg(i);
686 if (dfdp.getMultiVector() != Teuchos::null) {
687 dfdp.getMultiVector()->Export(
688 *(dfdp_sg.getMultiVector()->getBlockMultiVector()),
689 *interlace_overlapped_f_exporter,
Insert);