45#include "Teuchos_Assert.hpp"
46#include "EpetraExt_BlockUtility.h"
47#include "EpetraExt_BlockMultiVector.h"
57 const Teuchos::RCP<EpetraExt::ModelEvaluator>& me_,
61 const Teuchos::RCP<const Stokhos::ParallelData>& sg_parallel_data_,
62 const Teuchos::RCP<Teuchos::ParameterList>& params_,
69 num_sg_blocks(sg_basis->size()),
70 num_W_blocks(sg_basis->size()),
71 num_p_blocks(sg_basis->size()),
73 x_map(me->get_x_map()),
74 f_map(me->get_f_map()),
75 sg_parallel_data(sg_parallel_data_),
76 sg_comm(sg_parallel_data->getMultiComm()),
77 epetraCijk(sg_parallel_data->getEpetraCijk()),
80 sg_overlapped_x_map(),
82 sg_overlapped_f_map(),
83 sg_overlapped_x_importer(),
84 sg_overlapped_f_exporter(),
100 eval_W_with_f(false),
103 if (
x_map != Teuchos::null)
117 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
120 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
124 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
127 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
144 (*sg_x_init)[0] = *(
me->get_x_init());
150 std::string W_expansion_type =
151 params->get(
"Jacobian Expansion Type",
"Full");
152 if (W_expansion_type ==
"Linear")
156 Teuchos::RCP<Epetra_BlockMap> W_overlap_map =
168 Teuchos::RCP<Teuchos::ParameterList> sgPrecParams =
169 Teuchos::rcp(&(
params->sublist(
"SG Preconditioner")),
false);
177 InArgs me_inargs =
me->createInArgs();
178 OutArgs me_outargs =
me->createOutArgs();
179 num_p = me_inargs.Np();
182 for (
int i=0; i<
num_p; i++) {
183 if (me_inargs.supports(IN_ARG_p_sg, i))
193 std::string p_expansion_type =
194 params->get(
"Parameter Expansion Type",
"Full");
195 if (p_expansion_type ==
"Linear")
208 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
211 Teuchos::RCP<const Teuchos::Array<std::string> > p_names =
213 if (p_names != Teuchos::null) {
215 Teuchos::rcp(
new Teuchos::Array<std::string>(
num_sg_blocks*(p_names->size())));
216 for (
int j=0;
j<p_names->size();
j++) {
217 std::stringstream ss;
218 ss << (*p_names)[
j] <<
" -- SG Coefficient " << i;
231 num_g = me_outargs.Ng();
234 for (
int i=0; i<
num_g; i++) {
235 if (me_outargs.supports(OUT_ARG_g_sg, i))
246 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
263Teuchos::RCP<const Epetra_Map>
269Teuchos::RCP<const Epetra_Map>
275Teuchos::RCP<const Epetra_Map>
278 TEUCHOS_TEST_FOR_EXCEPTION(l < 0 || l >= num_p + num_p_sg, std::logic_error,
279 "Error! Invalid p map index " << l);
281 return me->get_p_map(l);
283 return sg_p_map[l-num_p];
285 return Teuchos::null;
288Teuchos::RCP<const Epetra_Map>
291 TEUCHOS_TEST_FOR_EXCEPTION(j < 0 || j >= num_g_sg, std::logic_error,
292 "Error! Invalid g map index " <<
j);
296Teuchos::RCP<const Teuchos::Array<std::string> >
299 TEUCHOS_TEST_FOR_EXCEPTION(l < 0 || l >= num_p + num_p_sg, std::logic_error,
300 "Error! Invalid p map index " << l);
302 return me->get_p_names(l);
304 return sg_p_names[l-num_p];
306 return Teuchos::null;
309Teuchos::RCP<const Epetra_Vector>
312 return sg_x_init->getBlockVector();
315Teuchos::RCP<const Epetra_Vector>
318 TEUCHOS_TEST_FOR_EXCEPTION(l < 0 || l >= num_p + num_p_sg, std::logic_error,
319 "Error! Invalid p map index " << l);
321 return me->get_p_init(l);
323 return sg_p_init[l-num_p]->getBlockVector();
325 return Teuchos::null;
328Teuchos::RCP<Epetra_Operator>
332 Teuchos::RCP<Teuchos::ParameterList> sgOpParams =
333 Teuchos::rcp(&(params->sublist(
"SG Operator")),
false);
334 Teuchos::RCP<Epetra_CrsMatrix> W_crs;
335 if (sgOpParams->get(
"Operator Method",
"Matrix Free") ==
337 W_crs = Teuchos::rcp_dynamic_cast<Epetra_CrsMatrix>(me->create_W(),
true);
338 Teuchos::RCP<const Epetra_CrsGraph> W_graph =
339 Teuchos::rcp(&(W_crs->Graph()),
false);
340 sgOpParams->set< Teuchos::RCP<const Epetra_CrsGraph> >(
"Base Graph",
344 my_W = sg_op_factory.
build(sg_comm, sg_basis, epetraCijk, x_map, f_map,
346 my_W->setupOperator(W_sg_blocks);
351 return Teuchos::null;
354Teuchos::RCP<EpetraExt::ModelEvaluator::Preconditioner>
359 Teuchos::RCP<Epetra_Operator> precOp =
360 sg_prec_factory->build(sg_comm, sg_basis, epetraCijk, x_map, sg_x_map);
361 return Teuchos::rcp(
new EpetraExt::ModelEvaluator::Preconditioner(precOp,
365 return Teuchos::null;
368Teuchos::RCP<Epetra_Operator>
371 TEUCHOS_TEST_FOR_EXCEPTION(j < 0 || j >= num_g_sg || !supports_x,
373 "Error: dg/dx index " <<
j <<
" is not supported!");
375 int jj = sg_g_index_map[
j];
376 Teuchos::RCP<const Epetra_Map> g_map = me->get_g_map(jj);
377 Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > sg_blocks;
378 OutArgs me_outargs = me->createOutArgs();
379 DerivativeSupport ds = me_outargs.supports(OUT_ARG_DgDx_sg, jj);
380 if (ds.supports(DERIV_LINEAR_OP)) {
383 sg_basis, overlapped_stoch_row_map, x_map,
384 g_map, sg_g_map[
j], sg_comm));
385 for (
unsigned int i=0; i<num_sg_blocks; i++)
386 sg_blocks->setCoeffPtr(i, me->create_DgDx_op(jj));
388 else if (ds.supports(DERIV_MV_BY_COL)) {
389 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_mv_blocks =
391 sg_basis, overlapped_stoch_row_map, g_map, sg_g_map[
j],
392 sg_comm, x_map->NumMyElements()));
395 sg_mv_blocks,
false));
397 else if (ds.supports(DERIV_TRANS_MV_BY_ROW)) {
398 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_mv_blocks =
400 sg_basis, overlapped_stoch_row_map, x_map, sg_x_map,
401 sg_comm, g_map->NumMyElements()));
404 sg_mv_blocks,
true));
407 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
408 "Error! me_outargs.supports(OUT_ARG_DgDx_sg, " << jj
409 <<
").none() is true!");
411 Teuchos::RCP<Teuchos::ParameterList> pl =
412 Teuchos::rcp(
new Teuchos::ParameterList);
413 Teuchos::RCP<Stokhos::SGOperator> dgdx_sg =
415 sg_comm, sg_basis, serialCijk, x_map, g_map,
416 sg_x_map, sg_g_map[
j], pl));
417 dgdx_sg->setupOperator(sg_blocks);
421Teuchos::RCP<Epetra_Operator>
424 TEUCHOS_TEST_FOR_EXCEPTION(j < 0 || j >= num_g_sg || !supports_x,
426 "Error: dg/dx_dot index " <<
j <<
" is not supported!");
428 int jj = sg_g_index_map[
j];
429 Teuchos::RCP<const Epetra_Map> g_map = me->get_g_map(jj);
430 Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > sg_blocks;
431 OutArgs me_outargs = me->createOutArgs();
432 DerivativeSupport ds = me_outargs.supports(OUT_ARG_DgDx_dot_sg, jj);
433 if (ds.supports(DERIV_LINEAR_OP)) {
436 sg_basis, overlapped_stoch_row_map, x_map,
437 g_map, sg_g_map[
j], sg_comm));
438 for (
unsigned int i=0; i<num_sg_blocks; i++)
439 sg_blocks->setCoeffPtr(i, me->create_DgDx_dot_op(jj));
441 else if (ds.supports(DERIV_MV_BY_COL)) {
442 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_mv_blocks =
444 sg_basis, overlapped_stoch_row_map, g_map, sg_g_map[
j],
445 sg_comm, x_map->NumMyElements()));
448 sg_mv_blocks,
false));
450 else if (ds.supports(DERIV_TRANS_MV_BY_ROW)) {
451 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_mv_blocks =
453 sg_basis, overlapped_stoch_row_map, x_map, sg_x_map,
454 sg_comm, g_map->NumMyElements()));
457 sg_mv_blocks,
true));
460 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
461 "Error! me_outargs.supports(OUT_ARG_DgDx_dot_sg, "
462 << jj <<
").none() is true!");
464 Teuchos::RCP<Teuchos::ParameterList> pl =
465 Teuchos::rcp(
new Teuchos::ParameterList);
466 Teuchos::RCP<Stokhos::SGOperator> dgdx_dot_sg =
468 sg_comm, sg_basis, serialCijk, x_map, g_map,
469 sg_x_map, sg_g_map[
j], pl));
470 dgdx_dot_sg->setupOperator(sg_blocks);
474Teuchos::RCP<Epetra_Operator>
477 TEUCHOS_TEST_FOR_EXCEPTION(
478 j < 0 || j >= num_g_sg || i < 0 || i >= num_p+num_p_sg,
480 "Error: dg/dp index " <<
j <<
"," << i <<
" is not supported!");
482 int jj = sg_g_index_map[
j];
483 Teuchos::RCP<const Epetra_Map> g_map = me->get_g_map(jj);
484 OutArgs me_outargs = me->createOutArgs();
486 if (me_outargs.supports(OUT_ARG_DgDp_sg,jj,i).supports(DERIV_LINEAR_OP)) {
487 Teuchos::RCP<Stokhos::EpetraOperatorOrthogPoly> sg_blocks =
489 sg_basis, overlapped_stoch_row_map,
490 me->get_p_map(i), me->get_g_map(jj),
491 sg_g_map[
j], sg_comm));
492 for (
unsigned int l=0; l<num_sg_blocks; l++)
493 sg_blocks->setCoeffPtr(l, me->create_DgDp_op(jj,i));
497 TEUCHOS_TEST_FOR_EXCEPTION(
498 true, std::logic_error,
499 "Error: Underlying model evaluator must support DERIV_LINER_OP " <<
500 "to create operator for dg/dp index " <<
j <<
"," << i <<
"!");
503 int ii = sg_p_index_map[i-num_p];
504 Teuchos::RCP<const Epetra_Map> p_map = me->get_p_map(ii);
505 Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > sg_blocks;
506 DerivativeSupport ds = me_outargs.supports(OUT_ARG_DgDp_sg,jj,ii);
507 if (ds.supports(DERIV_LINEAR_OP)) {
510 sg_basis, overlapped_stoch_row_map,
511 p_map, g_map, sg_g_map[
j], sg_comm));
512 for (
unsigned int l=0; l<num_sg_blocks; l++)
513 sg_blocks->setCoeffPtr(l, me->create_DgDp_op(jj,ii));
515 else if (ds.supports(DERIV_MV_BY_COL)) {
516 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_mv_blocks =
518 sg_basis, overlapped_stoch_row_map, g_map, sg_g_map[
j],
519 sg_comm, p_map->NumMyElements()));
522 sg_mv_blocks,
false));
524 else if (ds.supports(DERIV_TRANS_MV_BY_ROW)) {
525 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_mv_blocks =
527 sg_basis, overlapped_stoch_row_map, p_map,
529 sg_comm, g_map->NumMyElements()));
532 sg_mv_blocks,
true));
535 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
536 "Error! me_outargs.supports(OUT_ARG_DgDp_sg, " << jj
537 <<
"," << ii <<
").none() is true!");
539 Teuchos::RCP<Teuchos::ParameterList> pl =
540 Teuchos::rcp(
new Teuchos::ParameterList);
541 Teuchos::RCP<Stokhos::SGOperator> dgdp_sg =
543 sg_comm, sg_basis, serialCijk,
544 p_map, g_map, sg_p_map[i-num_p], sg_g_map[
j], pl));
545 dgdp_sg->setupOperator(sg_blocks);
549 return Teuchos::null;
552Teuchos::RCP<Epetra_Operator>
555 TEUCHOS_TEST_FOR_EXCEPTION(i < 0 || i >= num_p+num_p_sg,
557 "Error: df/dp index " << i <<
" is not supported!");
559 OutArgs me_outargs = me->createOutArgs();
561 if (me_outargs.supports(OUT_ARG_DfDp_sg,i).supports(DERIV_LINEAR_OP)) {
562 Teuchos::RCP<Stokhos::EpetraOperatorOrthogPoly> sg_blocks =
564 sg_basis, overlapped_stoch_row_map,
565 me->get_p_map(i), f_map,
566 sg_overlapped_f_map, sg_comm));
567 for (
unsigned int l=0; l<num_sg_blocks; l++)
568 sg_blocks->setCoeffPtr(l, me->create_DfDp_op(i));
572 TEUCHOS_TEST_FOR_EXCEPTION(
573 true, std::logic_error,
574 "Error: Underlying model evaluator must support DERIV_LINER_OP " <<
575 "to create operator for df/dp index " << i <<
"!");
578 int ii = sg_p_index_map[i-num_p];
579 Teuchos::RCP<const Epetra_Map> p_map = me->get_p_map(ii);
580 Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > sg_blocks;
581 DerivativeSupport ds = me_outargs.supports(OUT_ARG_DfDp_sg,ii);
582 if (ds.supports(DERIV_LINEAR_OP)) {
585 sg_basis, overlapped_stoch_row_map,
586 p_map, f_map, sg_overlapped_f_map, sg_comm));
587 for (
unsigned int l=0; l<num_sg_blocks; l++)
588 sg_blocks->setCoeffPtr(l, me->create_DfDp_op(ii));
590 else if (ds.supports(DERIV_MV_BY_COL)) {
591 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_mv_blocks =
593 sg_basis, overlapped_stoch_row_map, f_map, sg_f_map,
594 sg_comm, p_map->NumMyElements()));
597 sg_mv_blocks,
false));
599 else if (ds.supports(DERIV_TRANS_MV_BY_ROW)) {
600 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_mv_blocks =
602 sg_basis, overlapped_stoch_row_map, p_map,
604 sg_comm, f_map->NumMyElements()));
607 sg_mv_blocks,
true));
610 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::logic_error,
611 "Error! me_outargs.supports(OUT_ARG_DfDp_sg, " << ii
612 <<
").none() is true!");
614 Teuchos::RCP<Teuchos::ParameterList> pl =
615 Teuchos::rcp(
new Teuchos::ParameterList);
616 Teuchos::RCP<Stokhos::SGOperator> dfdp_sg =
618 sg_comm, sg_basis, epetraCijk,
619 p_map, f_map, sg_p_map[i-num_p], sg_f_map, pl));
620 dfdp_sg->setupOperator(sg_blocks);
624 return Teuchos::null;
627EpetraExt::ModelEvaluator::InArgs
631 InArgs me_inargs = me->createInArgs();
633 inArgs.setModelEvalDescription(this->description());
634 inArgs.set_Np(num_p+num_p_sg);
635 inArgs.setSupports(IN_ARG_x_dot, me_inargs.supports(IN_ARG_x_dot_sg));
636 inArgs.setSupports(IN_ARG_x, me_inargs.supports(IN_ARG_x_sg));
637 inArgs.setSupports(IN_ARG_t, me_inargs.supports(IN_ARG_t));
638 inArgs.setSupports(IN_ARG_alpha, me_inargs.supports(IN_ARG_alpha));
639 inArgs.setSupports(IN_ARG_beta, me_inargs.supports(IN_ARG_beta));
640 inArgs.setSupports(IN_ARG_sg_basis, me_inargs.supports(IN_ARG_sg_basis));
641 inArgs.setSupports(IN_ARG_sg_quadrature,
642 me_inargs.supports(IN_ARG_sg_quadrature));
643 inArgs.setSupports(IN_ARG_sg_expansion,
644 me_inargs.supports(IN_ARG_sg_expansion));
649EpetraExt::ModelEvaluator::OutArgs
652 OutArgsSetup outArgs;
653 OutArgs me_outargs = me->createOutArgs();
655 outArgs.setModelEvalDescription(this->description());
656 outArgs.set_Np_Ng(num_p+num_p_sg, num_g_sg);
657 outArgs.setSupports(OUT_ARG_f, me_outargs.supports(OUT_ARG_f_sg));
658 outArgs.setSupports(OUT_ARG_W, me_outargs.supports(OUT_ARG_W_sg));
661 me_outargs.supports(OUT_ARG_W_sg) && sg_prec_factory->isPrecSupported());
662 for (
int j=0;
j<num_p;
j++)
663 outArgs.setSupports(OUT_ARG_DfDp,
j,
664 me_outargs.supports(OUT_ARG_DfDp_sg,
j));
665 for (
int j=0;
j<num_p_sg;
j++)
666 if (!me_outargs.supports(OUT_ARG_DfDp_sg, sg_p_index_map[
j]).none())
667 outArgs.setSupports(OUT_ARG_DfDp,
j+num_p, DERIV_LINEAR_OP);
668 for (
int i=0; i<num_g_sg; i++) {
669 int ii = sg_g_index_map[i];
670 if (!me_outargs.supports(OUT_ARG_DgDx_dot_sg, ii).none())
671 outArgs.setSupports(OUT_ARG_DgDx_dot, i, DERIV_LINEAR_OP);
672 if (!me_outargs.supports(OUT_ARG_DgDx_sg, ii).none())
673 outArgs.setSupports(OUT_ARG_DgDx, i, DERIV_LINEAR_OP);
674 for (
int j=0;
j<num_p;
j++)
675 outArgs.setSupports(OUT_ARG_DgDp, i,
j,
676 me_outargs.supports(OUT_ARG_DgDp_sg, ii,
j));
677 for (
int j=0;
j<num_p_sg;
j++)
678 if (!me_outargs.supports(OUT_ARG_DgDp_sg, ii, sg_p_index_map[
j]).none())
679 outArgs.setSupports(OUT_ARG_DgDp, i,
j+num_p, DERIV_LINEAR_OP);
687 const OutArgs& outArgs)
const
690 Teuchos::RCP<const Epetra_Vector> x;
691 if (inArgs.supports(IN_ARG_x)) {
693 if (x != Teuchos::null)
696 Teuchos::RCP<const Epetra_Vector> x_dot;
697 if (inArgs.supports(IN_ARG_x_dot))
698 x_dot = inArgs.get_x_dot();
701 EpetraExt::ModelEvaluator::Evaluation<Epetra_Vector> f_out;
702 if (outArgs.supports(OUT_ARG_f))
703 f_out = outArgs.get_f();
704 Teuchos::RCP<Epetra_Operator> W_out;
705 if (outArgs.supports(OUT_ARG_W))
706 W_out = outArgs.get_W();
707 Teuchos::RCP<Epetra_Operator> WPrec_out;
708 if (outArgs.supports(OUT_ARG_WPrec))
709 WPrec_out = outArgs.get_WPrec();
713 bool eval_prec = (W_out == Teuchos::null && WPrec_out != Teuchos::null);
719 Teuchos::RCP<Stokhos::SGPreconditioner> W_prec =
720 Teuchos::rcp_dynamic_cast<Stokhos::SGPreconditioner>(WPrec_out,
true);
721 W_prec->setupPreconditioner(my_W, *my_x);
724 bool done = (f_out == Teuchos::null);
725 for (
int i=0; i<outArgs.Ng(); i++) {
726 done = done && (outArgs.get_g(i) == Teuchos::null);
727 for (
int j=0;
j<outArgs.Np();
j++)
728 if (!outArgs.supports(OUT_ARG_DgDp, i,
j).none())
729 done = done && (outArgs.get_DgDp(i,
j).isEmpty());
736 InArgs me_inargs = me->createInArgs();
737 if (x != Teuchos::null) {
738 x_sg_blocks->getBlockVector()->Import(*x, *sg_overlapped_x_importer,
740 me_inargs.set_x_sg(x_sg_blocks);
742 if (x_dot != Teuchos::null) {
743 x_dot_sg_blocks->getBlockVector()->Import(*x_dot, *sg_overlapped_x_importer,
745 me_inargs.set_x_dot_sg(x_dot_sg_blocks);
747 if (me_inargs.supports(IN_ARG_alpha))
748 me_inargs.set_alpha(inArgs.get_alpha());
749 if (me_inargs.supports(IN_ARG_beta))
750 me_inargs.set_beta(inArgs.get_beta());
751 if (me_inargs.supports(IN_ARG_t))
752 me_inargs.set_t(inArgs.get_t());
753 if (me_inargs.supports(IN_ARG_sg_basis)) {
754 if (inArgs.get_sg_basis() != Teuchos::null)
755 me_inargs.set_sg_basis(inArgs.get_sg_basis());
757 me_inargs.set_sg_basis(sg_basis);
759 if (me_inargs.supports(IN_ARG_sg_quadrature)) {
760 if (inArgs.get_sg_quadrature() != Teuchos::null)
761 me_inargs.set_sg_quadrature(inArgs.get_sg_quadrature());
763 me_inargs.set_sg_quadrature(sg_quad);
765 if (me_inargs.supports(IN_ARG_sg_expansion)) {
766 if (inArgs.get_sg_expansion() != Teuchos::null)
767 me_inargs.set_sg_expansion(inArgs.get_sg_expansion());
769 me_inargs.set_sg_expansion(sg_exp);
773 for (
int i=0; i<num_p; i++)
774 me_inargs.set_p(i, inArgs.get_p(i));
775 for (
int i=0; i<num_p_sg; i++) {
776 Teuchos::RCP<const Epetra_Vector> p = inArgs.get_p(i+num_p);
780 if (p == Teuchos::null)
781 p = sg_p_init[i]->getBlockVector();
784 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> p_sg =
785 create_p_sg(sg_p_index_map[i],
View, p.get());
786 me_inargs.set_p_sg(sg_p_index_map[i], p_sg);
790 OutArgs me_outargs = me->createOutArgs();
793 if (f_out != Teuchos::null) {
794 me_outargs.set_f_sg(f_sg_blocks);
796 me_outargs.set_W_sg(W_sg_blocks);
800 if (W_out != Teuchos::null && !eval_W_with_f && !eval_prec)
801 me_outargs.set_W_sg(W_sg_blocks);
804 for (
int i=0; i<num_p; i++) {
805 if (!outArgs.supports(OUT_ARG_DfDp, i).none()) {
806 Derivative dfdp = outArgs.get_DfDp(i);
807 if (dfdp.getMultiVector() != Teuchos::null) {
808 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dfdp_sg;
809 if (dfdp.getMultiVectorOrientation() == DERIV_MV_BY_COL)
812 sg_basis, overlapped_stoch_row_map,
813 me->get_f_map(), sg_overlapped_f_map, sg_comm,
814 me->get_p_map(i)->NumMyElements()));
815 else if (dfdp.getMultiVectorOrientation() == DERIV_TRANS_MV_BY_ROW)
818 sg_basis, overlapped_stoch_row_map,
819 me->get_p_map(i), sg_comm,
820 me->get_f_map()->NumMyElements()));
821 me_outargs.set_DfDp_sg(
822 i, SGDerivative(dfdp_sg, dfdp.getMultiVectorOrientation()));
824 else if (dfdp.getLinearOp() != Teuchos::null) {
825 Teuchos::RCP<Stokhos::EpetraOperatorOrthogPoly> dfdp_sg =
826 Teuchos::rcp_dynamic_cast<Stokhos::EpetraOperatorOrthogPoly>(dfdp.getLinearOp(),
true);
827 me_outargs.set_DfDp_sg(i, SGDerivative(dfdp_sg));
833 for (
int i=0; i<num_p_sg; i++) {
834 if (!outArgs.supports(OUT_ARG_DfDp, i+num_p).none()) {
835 Derivative dfdp = outArgs.get_DfDp(i+num_p);
836 if (dfdp.getLinearOp() != Teuchos::null) {
837 Teuchos::RCP<Stokhos::MatrixFreeOperator> dfdp_op =
838 Teuchos::rcp_dynamic_cast<Stokhos::MatrixFreeOperator>(dfdp.getLinearOp(),
true);
839 Teuchos::RCP<Stokhos::EpetraOperatorOrthogPoly > dfdp_op_sg =
840 dfdp_op->getSGPolynomial();
841 int ii = sg_p_index_map[i];
842 if (me_outargs.supports(OUT_ARG_DfDp_sg,ii).supports(DERIV_LINEAR_OP))
843 me_outargs.set_DfDp_sg(ii, SGDerivative(dfdp_op_sg));
845 Teuchos::RCP<Stokhos::EpetraMultiVectorOperatorOrthogPoly> sg_mv_op =
846 Teuchos::rcp_dynamic_cast<Stokhos::EpetraMultiVectorOperatorOrthogPoly>(dfdp_op_sg,
true);
847 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dfdp_sg =
848 sg_mv_op->multiVectorOrthogPoly();
849 if (me_outargs.supports(OUT_ARG_DfDp_sg,ii).supports(DERIV_MV_BY_COL))
850 me_outargs.set_DfDp_sg(
851 ii, SGDerivative(dfdp_sg, DERIV_MV_BY_COL));
853 me_outargs.set_DfDp_sg(
854 ii, SGDerivative(dfdp_sg, DERIV_TRANS_MV_BY_ROW));
857 TEUCHOS_TEST_FOR_EXCEPTION(
858 dfdp.getLinearOp() == Teuchos::null && dfdp.isEmpty() ==
false,
860 "Error! Stokhos::SGModelEvaluator::evalModel: " <<
861 "Operator form of df/dp(" << i+num_p <<
") is required!");
866 for (
int i=0; i<num_g_sg; i++) {
867 int ii = sg_g_index_map[i];
870 Teuchos::RCP<Epetra_Vector>
g = outArgs.get_g(i);
871 if (
g != Teuchos::null) {
872 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> g_sg =
873 create_g_sg(sg_g_index_map[i],
View,
g.get());
874 me_outargs.set_g_sg(ii, g_sg);
878 if (outArgs.supports(OUT_ARG_DgDx_dot, i).supports(DERIV_LINEAR_OP)) {
879 Derivative dgdx_dot = outArgs.get_DgDx_dot(i);
880 if (dgdx_dot.getLinearOp() != Teuchos::null) {
881 Teuchos::RCP<Stokhos::SGOperator> op =
882 Teuchos::rcp_dynamic_cast<Stokhos::SGOperator>(
883 dgdx_dot.getLinearOp(),
true);
884 Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > sg_blocks =
885 op->getSGPolynomial();
886 if (me_outargs.supports(OUT_ARG_DgDx, ii).supports(DERIV_LINEAR_OP))
887 me_outargs.set_DgDx_dot_sg(ii, sg_blocks);
889 Teuchos::RCP<Stokhos::EpetraMultiVectorOperatorOrthogPoly> sg_mv_op =
890 Teuchos::rcp_dynamic_cast<Stokhos::EpetraMultiVectorOperatorOrthogPoly>(sg_blocks,
true);
891 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dgdx_dot_sg =
892 sg_mv_op->multiVectorOrthogPoly();
893 if (me_outargs.supports(OUT_ARG_DgDx_dot_sg, ii).supports(DERIV_MV_BY_COL))
894 me_outargs.set_DgDx_dot_sg(ii, SGDerivative(dgdx_dot_sg,
897 me_outargs.set_DgDx_dot_sg(ii, SGDerivative(dgdx_dot_sg,
898 DERIV_TRANS_MV_BY_ROW));
901 TEUCHOS_TEST_FOR_EXCEPTION(
902 dgdx_dot.getLinearOp() == Teuchos::null && dgdx_dot.isEmpty() ==
false,
904 "Error! Stokhos::SGModelEvaluator::evalModel: " <<
905 "Operator form of dg/dxdot(" << i <<
") is required!");
909 if (outArgs.supports(OUT_ARG_DgDx, i).supports(DERIV_LINEAR_OP)) {
910 Derivative dgdx = outArgs.get_DgDx(i);
911 if (dgdx.getLinearOp() != Teuchos::null) {
912 Teuchos::RCP<Stokhos::SGOperator> op =
913 Teuchos::rcp_dynamic_cast<Stokhos::SGOperator>(
914 dgdx.getLinearOp(),
true);
915 Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > sg_blocks =
916 op->getSGPolynomial();
917 if (me_outargs.supports(OUT_ARG_DgDx, ii).supports(DERIV_LINEAR_OP))
918 me_outargs.set_DgDx_sg(ii, sg_blocks);
920 Teuchos::RCP<Stokhos::EpetraMultiVectorOperatorOrthogPoly> sg_mv_op =
921 Teuchos::rcp_dynamic_cast<Stokhos::EpetraMultiVectorOperatorOrthogPoly>(sg_blocks,
true);
922 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dgdx_sg =
923 sg_mv_op->multiVectorOrthogPoly();
924 if (me_outargs.supports(OUT_ARG_DgDx_sg, ii).supports(DERIV_MV_BY_COL))
925 me_outargs.set_DgDx_sg(ii, SGDerivative(dgdx_sg,
928 me_outargs.set_DgDx_sg(ii, SGDerivative(dgdx_sg,
929 DERIV_TRANS_MV_BY_ROW));
932 TEUCHOS_TEST_FOR_EXCEPTION(
933 dgdx.getLinearOp() == Teuchos::null && dgdx.isEmpty() ==
false,
935 "Error! Stokhos::SGModelEvaluator::evalModel: " <<
936 "Operator form of dg/dx(" << i <<
") is required!");
940 for (
int j=0;
j<num_p;
j++) {
941 if (!outArgs.supports(OUT_ARG_DgDp, i,
j).none()) {
942 Derivative dgdp = outArgs.get_DgDp(i,
j);
943 if (dgdp.getMultiVector() != Teuchos::null) {
944 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dgdp_sg;
945 if (dgdp.getMultiVectorOrientation() == DERIV_MV_BY_COL)
948 sg_basis, overlapped_stoch_row_map,
949 me->get_g_map(ii), sg_g_map[i], sg_comm,
950 View, *(dgdp.getMultiVector())));
951 else if (dgdp.getMultiVectorOrientation() == DERIV_TRANS_MV_BY_ROW) {
952 Teuchos::RCP<const Epetra_BlockMap> product_map =
953 Teuchos::rcp(&(dgdp.getMultiVector()->Map()),
false);
956 sg_basis, overlapped_stoch_row_map,
957 me->get_p_map(
j), product_map, sg_comm,
958 View, *(dgdp.getMultiVector())));
960 me_outargs.set_DgDp_sg(
961 ii,
j, SGDerivative(dgdp_sg, dgdp.getMultiVectorOrientation()));
963 else if (dgdp.getLinearOp() != Teuchos::null) {
964 Teuchos::RCP<Stokhos::EpetraOperatorOrthogPoly> dgdp_sg =
965 Teuchos::rcp_dynamic_cast<Stokhos::EpetraOperatorOrthogPoly>(dgdp.getLinearOp(),
true);
966 me_outargs.set_DgDp_sg(ii,
j, SGDerivative(dgdp_sg));
972 for (
int j=0;
j<num_p_sg;
j++) {
973 if (!outArgs.supports(OUT_ARG_DgDp, i,
j+num_p).none()) {
974 Derivative dgdp = outArgs.get_DgDp(i,
j+num_p);
975 if (dgdp.getLinearOp() != Teuchos::null) {
976 Teuchos::RCP<Stokhos::MatrixFreeOperator> dgdp_op =
977 Teuchos::rcp_dynamic_cast<Stokhos::MatrixFreeOperator>(dgdp.getLinearOp(),
true);
978 Teuchos::RCP<Stokhos::EpetraOperatorOrthogPoly > dgdp_op_sg =
979 dgdp_op->getSGPolynomial();
980 int jj = sg_p_index_map[
j];
981 if (me_outargs.supports(OUT_ARG_DgDp_sg,ii,jj).supports(DERIV_LINEAR_OP))
982 me_outargs.set_DgDp_sg(ii, jj, SGDerivative(dgdp_op_sg));
984 Teuchos::RCP<Stokhos::EpetraMultiVectorOperatorOrthogPoly> sg_mv_op =
985 Teuchos::rcp_dynamic_cast<Stokhos::EpetraMultiVectorOperatorOrthogPoly>(dgdp_op_sg,
true);
986 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dgdp_sg =
987 sg_mv_op->multiVectorOrthogPoly();
988 if (me_outargs.supports(OUT_ARG_DgDp_sg,ii,jj).supports(DERIV_MV_BY_COL))
989 me_outargs.set_DgDp_sg(
990 ii, jj, SGDerivative(dgdp_sg, DERIV_MV_BY_COL));
992 me_outargs.set_DgDp_sg(
993 ii, jj, SGDerivative(dgdp_sg, DERIV_TRANS_MV_BY_ROW));
996 TEUCHOS_TEST_FOR_EXCEPTION(
997 dgdp.getLinearOp() == Teuchos::null && dgdp.isEmpty() ==
false,
999 "Error! Stokhos::SGModelEvaluator::evalModel: " <<
1000 "Operator form of dg/dp(" << i <<
"," <<
j+num_p <<
") is required!");
1007 me->evalModel(me_inargs, me_outargs);
1010 if ((W_out != Teuchos::null || (eval_W_with_f && f_out != Teuchos::null))
1012 Teuchos::RCP<Epetra_Operator> W;
1013 if (W_out != Teuchos::null)
1017 Teuchos::RCP<Stokhos::SGOperator> W_sg =
1018 Teuchos::rcp_dynamic_cast<Stokhos::SGOperator>(W,
true);
1019 W_sg->setupOperator(W_sg_blocks);
1021 if (WPrec_out != Teuchos::null) {
1022 Teuchos::RCP<Stokhos::SGPreconditioner> W_prec =
1023 Teuchos::rcp_dynamic_cast<Stokhos::SGPreconditioner>(WPrec_out,
true);
1024 W_prec->setupPreconditioner(W_sg, *my_x);
1029 if (f_out!=Teuchos::null){
1031 for (
int i=0; i<sg_basis->size(); i++)
1032 (*f_sg_blocks)[i].Scale(sg_basis->norm_squared(i));
1033 f_out->Export(*(f_sg_blocks->getBlockVector()), *sg_overlapped_f_exporter,
1038 for (
int i=0; i<num_p; i++) {
1039 if (!outArgs.supports(OUT_ARG_DfDp, i).none()) {
1040 Derivative dfdp = outArgs.get_DfDp(i);
1041 SGDerivative dfdp_sg = me_outargs.get_DfDp_sg(i);
1042 if (dfdp.getMultiVector() != Teuchos::null) {
1043 dfdp.getMultiVector()->Export(
1044 *(dfdp_sg.getMultiVector()->getBlockMultiVector()),
1045 *sg_overlapped_f_exporter,
Insert);
1062 *sg_x_init = x_sg_in;
1065Teuchos::RCP<const Stokhos::EpetraVectorOrthogPoly>
1075 Teuchos::Array<int>::iterator it = std::find(sg_p_index_map.begin(),
1076 sg_p_index_map.end(),
1078 TEUCHOS_TEST_FOR_EXCEPTION(it == sg_p_index_map.end(), std::logic_error,
1079 "Error! Invalid p map index " << i);
1080 int ii = it - sg_p_index_map.begin();
1081 *sg_p_init[ii] = p_sg_in;
1084Teuchos::RCP<const Stokhos::EpetraVectorOrthogPoly>
1087 Teuchos::Array<int>::const_iterator it = std::find(sg_p_index_map.begin(),
1088 sg_p_index_map.end(),
1090 TEUCHOS_TEST_FOR_EXCEPTION(it == sg_p_index_map.end(), std::logic_error,
1091 "Error! Invalid p map index " << l);
1092 int ll = it - sg_p_index_map.begin();
1093 return sg_p_init[ll];
1099 return sg_p_index_map;
1105 return sg_g_index_map;
1108Teuchos::Array< Teuchos::RCP<const Epetra_Map> >
1111 Teuchos::Array< Teuchos::RCP<const Epetra_Map> > base_maps(num_g);
1112 for (
int i=0; i<num_g; i++)
1113 base_maps[i] = me->get_g_map(i);
1117Teuchos::RCP<const Epetra_BlockMap>
1120 return overlapped_stoch_row_map;
1123Teuchos::RCP<const Epetra_BlockMap>
1126 return sg_overlapped_x_map;
1129Teuchos::RCP<const Epetra_Import>
1132 return sg_overlapped_x_importer;
1135Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
1139 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> sg_x;
1142 sg_basis, stoch_row_map, x_map, sg_x_map, sg_comm));
1145 sg_basis, stoch_row_map, x_map, sg_x_map, sg_comm,
1150Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
1154 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> sg_x;
1157 sg_basis, overlapped_stoch_row_map, x_map,
1158 sg_overlapped_x_map, sg_comm));
1161 sg_basis, overlapped_stoch_row_map, x_map,
1162 sg_overlapped_x_map, sg_comm, CV, *v));
1166Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly>
1170 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_x;
1173 sg_basis, stoch_row_map, x_map, sg_x_map, sg_comm,
1177 sg_basis, stoch_row_map, x_map, sg_x_map, sg_comm,
1182Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly>
1188 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_x;
1191 sg_basis, overlapped_stoch_row_map, x_map,
1192 sg_overlapped_x_map, sg_comm, num_vecs));
1195 sg_basis, overlapped_stoch_row_map, x_map,
1196 sg_overlapped_x_map, sg_comm, CV, *v));
1200Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
1204 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> sg_p;
1205 Teuchos::Array<int>::const_iterator it = std::find(sg_p_index_map.begin(),
1206 sg_p_index_map.end(),
1208 TEUCHOS_TEST_FOR_EXCEPTION(it == sg_p_index_map.end(), std::logic_error,
1209 "Error! Invalid p map index " << l);
1210 int ll = it - sg_p_index_map.begin();
1213 sg_basis, overlapped_stoch_p_map, me->get_p_map(l),
1214 sg_p_map[ll], sg_comm));
1217 sg_basis, overlapped_stoch_p_map, me->get_p_map(l),
1218 sg_p_map[ll], sg_comm, CV, *v));
1222Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
1226 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> sg_f;
1229 sg_basis, stoch_row_map, f_map, sg_f_map, sg_comm));
1232 sg_basis, stoch_row_map, f_map, sg_f_map, sg_comm,
1237Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
1241 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> sg_f;
1244 sg_basis, overlapped_stoch_row_map, f_map,
1245 sg_overlapped_f_map, sg_comm));
1248 sg_basis, overlapped_stoch_row_map, f_map,
1249 sg_overlapped_f_map, sg_comm, CV, *v));
1253Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly>
1259 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_f;
1262 sg_basis, stoch_row_map, f_map, sg_f_map, sg_comm,
1266 sg_basis, stoch_row_map, f_map, sg_f_map, sg_comm,
1271Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly>
1277 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_f;
1280 sg_basis, overlapped_stoch_row_map, f_map,
1281 sg_overlapped_f_map, sg_comm, num_vecs));
1284 sg_basis, overlapped_stoch_row_map, f_map,
1285 sg_overlapped_f_map, sg_comm, CV, *v));
1289Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
1293 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> sg_g;
1294 Teuchos::Array<int>::const_iterator it = std::find(sg_g_index_map.begin(),
1295 sg_g_index_map.end(),
1297 TEUCHOS_TEST_FOR_EXCEPTION(it == sg_g_index_map.end(), std::logic_error,
1298 "Error! Invalid g map index " << l);
1299 int ll = it - sg_g_index_map.begin();
1302 sg_basis, overlapped_stoch_row_map,
1304 sg_g_map[ll], sg_comm));
1307 sg_basis, overlapped_stoch_row_map,
1309 sg_g_map[ll], sg_comm, CV, *v));
1313Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly>
1318 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_g;
1319 Teuchos::Array<int>::const_iterator it = std::find(sg_g_index_map.begin(),
1320 sg_g_index_map.end(),
1322 TEUCHOS_TEST_FOR_EXCEPTION(it == sg_g_index_map.end(), std::logic_error,
1323 "Error! Invalid g map index " << l);
1324 int ll = it - sg_g_index_map.begin();
1327 sg_basis, overlapped_stoch_row_map,
1329 sg_g_map[ll], sg_comm, num_vecs));
1332 sg_basis, overlapped_stoch_row_map,
1334 sg_g_map[ll], sg_comm, CV, *v));
1338Teuchos::RCP<EpetraExt::BlockVector>
1341 Teuchos::RCP<EpetraExt::BlockVector> x_overlapped =
1342 Teuchos::rcp(
new EpetraExt::BlockVector(*x_map, *sg_overlapped_x_map));
1343 x_overlapped->Import(x, *sg_overlapped_x_importer,
Insert);
1344 return x_overlapped;
1347Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
1350 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> sg_x_overlapped =
1352 sg_basis, overlapped_stoch_row_map, x_map,
1353 sg_overlapped_x_map, sg_comm));
1354 sg_x_overlapped->getBlockVector()->Import(x, *sg_overlapped_x_importer,
1356 return sg_x_overlapped;
1359Teuchos::RCP<EpetraExt::BlockVector>
1362 Teuchos::RCP<EpetraExt::BlockVector> x =
1363 Teuchos::rcp(
new EpetraExt::BlockVector(*x_map, *sg_x_map));
1364 x->Export(x_overlapped, *sg_overlapped_x_importer,
Insert);
1368Teuchos::RCP<EpetraExt::BlockVector>
1371 Teuchos::RCP<EpetraExt::BlockVector> f_overlapped =
1372 Teuchos::rcp(
new EpetraExt::BlockVector(*f_map, *sg_overlapped_f_map));
1373 f_overlapped->Import(
f, *sg_overlapped_f_exporter,
Insert);
1374 return f_overlapped;
1377Teuchos::RCP<EpetraExt::BlockVector>
1380 Teuchos::RCP<EpetraExt::BlockVector>
f =
1381 Teuchos::rcp(
new EpetraExt::BlockVector(*f_map, *sg_f_map));
1382 f->Export(f_overlapped, *sg_overlapped_f_exporter,
Insert);
A container class storing an orthogonal polynomial whose coefficients are vectors,...
A container class storing an orthogonal polynomial whose coefficients are vectors,...
A container class storing an orthogonal polynomial whose coefficients are vectors,...
A container class storing an orthogonal polynomial whose coefficients are vectors,...
An Epetra operator representing the block stochastic Galerkin operator.
Abstract base class for multivariate orthogonal polynomials.
Abstract base class for orthogonal polynomial-based expansions.
Abstract base class for quadrature methods.
Teuchos::RCP< const Stokhos::OrthogPolyBasis< int, double > > sg_basis
Stochastic Galerkin basis.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > create_x_sg_overlap(Epetra_DataAccess CV=Copy, const Epetra_Vector *v=NULL) const
Create vector orthog poly using x map and overlap sg map.
Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > create_f_mv_sg_overlap(int num_vecs, Epetra_DataAccess CV=Copy, const Epetra_MultiVector *v=NULL) const
Create multi-vector orthog poly using f map and overlap sg map.
OutArgs createOutArgs() const
Create OutArgs.
Teuchos::RCP< EpetraExt::BlockVector > export_residual(const Epetra_Vector &f_overlapped) const
Export parallel residual vector.
Teuchos::RCP< const Stokhos::ParallelData > sg_parallel_data
Parallel SG data.
Teuchos::RCP< const Epetra_Map > sg_f_map
Block SG residual map.
int num_p_sg
Number of stochastic parameter vectors.
Teuchos::RCP< Epetra_Operator > create_DgDx_dot_op(int j) const
Create SG operator representing dg/dxdot.
Teuchos::RCP< EpetraExt::ModelEvaluator > me
Underlying model evaluator.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > sg_x_init
SG initial x.
int num_g_sg
Number of stochastic response vectors.
Teuchos::RCP< const Teuchos::Array< std::string > > get_p_names(int l) const
Return array of parameter names.
Teuchos::RCP< Stokhos::SGPreconditionerFactory > sg_prec_factory
Preconditioner factory.
Teuchos::RCP< const Stokhos::EpetraSparse3Tensor > serialCijk
Serial Epetra Cijk for dgdx*.
Teuchos::RCP< const Epetra_Map > get_g_map(int l) const
Return response map.
Teuchos::Array< int > get_g_sg_map_indices() const
Get indices of SG responses.
Teuchos::RCP< const Epetra_Map > f_map
Underlying residual map.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > create_f_sg_overlap(Epetra_DataAccess CV=Copy, const Epetra_Vector *v=NULL) const
Create vector orthog poly using f map and overlap sg map.
Teuchos::RCP< const Epetra_Map > get_p_map(int l) const
Return parameter vector map.
Teuchos::RCP< const Epetra_Vector > get_x_init() const
Return initial solution.
unsigned int num_sg_blocks
Number of stochastic blocks.
Teuchos::RCP< const Epetra_Map > sg_overlapped_x_map
Block SG overlapped unknown map.
Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > create_f_mv_sg(int num_vecs, Epetra_DataAccess CV=Copy, const Epetra_MultiVector *v=NULL) const
Create multi-vector orthog poly using f map and owned sg map.
Teuchos::RCP< const Epetra_BlockMap > overlapped_stoch_row_map
Overlapped map for stochastic blocks (local map)
InArgs createInArgs() const
Create InArgs.
Teuchos::RCP< const Epetra_Map > sg_x_map
Block SG unknown map.
Teuchos::RCP< EpetraExt::BlockVector > import_solution(const Epetra_Vector &x) const
Import parallel solution vector.
Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > create_x_mv_sg_overlap(int num_vecs, Epetra_DataAccess CV=Copy, const Epetra_MultiVector *v=NULL) const
Create vector orthog poly using x map and overlap sg map.
Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > create_x_mv_sg(int num_vecs, Epetra_DataAccess CV=Copy, const Epetra_MultiVector *v=NULL) const
Create vector orthog poly using x map and owned sg map.
Teuchos::RCP< Epetra_Operator > create_W() const
Create W = alpha*M + beta*J matrix.
void set_x_sg_init(const Stokhos::EpetraVectorOrthogPoly &x_sg_in)
Set initial solution polynomial.
Teuchos::Array< Teuchos::RCP< const Epetra_Map > > sg_g_map
Block SG response map.
Teuchos::RCP< Epetra_Import > sg_overlapped_x_importer
Importer from SG to SG-overlapped maps.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > create_x_sg(Epetra_DataAccess CV=Copy, const Epetra_Vector *v=NULL) const
Create vector orthog poly using x map and owned sg map.
Teuchos::RCP< const Epetra_BlockMap > overlapped_stoch_p_map
Overlapped map for p stochastic blocks (local map)
bool eval_W_with_f
Whether to always evaluate W with f.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > create_g_sg(int l, Epetra_DataAccess CV=Copy, const Epetra_Vector *v=NULL) const
Create vector orthog poly using g map.
Teuchos::RCP< Teuchos::ParameterList > params
Algorithmic parameters.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > x_sg_blocks
x stochastic Galerkin components
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > create_p_sg(int l, Epetra_DataAccess CV=Copy, const Epetra_Vector *v=NULL) const
Create vector orthog poly using p map.
Teuchos::Array< int > sg_g_index_map
Index map between block-g and g_sg maps.
Teuchos::Array< Teuchos::RCP< const Epetra_Map > > sg_p_map
Block SG parameter map.
Teuchos::RCP< EpetraExt::BlockVector > export_solution(const Epetra_Vector &x_overlapped) const
Export parallel solution vector.
Teuchos::RCP< const Epetra_Vector > get_p_init(int l) const
Return initial parameters.
Teuchos::Array< Teuchos::RCP< Teuchos::Array< std::string > > > sg_p_names
SG coefficient parameter names.
Teuchos::RCP< Epetra_Vector > my_x
x pointer for evaluating preconditioner
Teuchos::RCP< Epetra_Export > sg_overlapped_f_exporter
Exporter from SG-overlapped to SG maps.
Teuchos::Array< int > sg_p_index_map
Index map between block-p and p_sg maps.
Teuchos::RCP< const Epetra_Import > get_x_sg_importer() const
Return x sg importer.
Teuchos::RCP< Epetra_Operator > create_DgDx_op(int j) const
Create SG operator representing dg/dx.
Teuchos::RCP< const EpetraExt::MultiComm > sg_comm
Parallel SG communicator.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > import_solution_poly(const Epetra_Vector &x) const
Import parallel solution vector.
Teuchos::RCP< const Epetra_Map > get_x_map() const
Return solution vector map.
Teuchos::RCP< Epetra_Operator > create_DgDp_op(int j, int i) const
Create SG operator representing dg/dp.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > x_dot_sg_blocks
x_dot stochastic Galerkin components
Teuchos::RCP< EpetraExt::ModelEvaluator::Preconditioner > create_WPrec() const
Create preconditioner operator.
Teuchos::RCP< const Epetra_BlockMap > get_overlap_stochastic_map() const
Return overlap stochastic map.
Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > W_sg_blocks
W stochastic Galerkin components.
Teuchos::Array< int > get_p_sg_map_indices() const
Get indices of SG parameters.
Teuchos::RCP< const Epetra_BlockMap > stoch_row_map
Map for stochastic blocks.
SGModelEvaluator(const Teuchos::RCP< EpetraExt::ModelEvaluator > &me, const Teuchos::RCP< const Stokhos::OrthogPolyBasis< int, double > > &sg_basis, const Teuchos::RCP< const Stokhos::Quadrature< int, double > > &sg_quad, const Teuchos::RCP< Stokhos::OrthogPolyExpansion< int, double > > &sg_exp, const Teuchos::RCP< const Stokhos::ParallelData > &sg_parallel_data, const Teuchos::RCP< Teuchos::ParameterList > ¶ms, bool scaleOP=true)
Teuchos::Array< Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > > sg_p_init
SG initial p.
Teuchos::RCP< const Stokhos::EpetraVectorOrthogPoly > get_x_sg_init() const
Return initial SG x.
Teuchos::RCP< const Epetra_Map > get_f_map() const
Return residual vector map.
void set_p_sg_init(int i, const Stokhos::EpetraVectorOrthogPoly &p_sg_in)
Set initial parameter polynomial.
Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > create_g_mv_sg(int l, int num_vecs, Epetra_DataAccess CV=Copy, const Epetra_MultiVector *v=NULL) const
Create multi-vector orthog poly using g map.
Teuchos::RCP< const Epetra_BlockMap > get_x_sg_overlap_map() const
Return x sg overlap map.
int num_p
Number of parameter vectors of underlying model evaluator.
void evalModel(const InArgs &inArgs, const OutArgs &outArgs) const
Evaluate model on InArgs.
Teuchos::RCP< const Stokhos::EpetraSparse3Tensor > epetraCijk
Epetra Cijk.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > create_f_sg(Epetra_DataAccess CV=Copy, const Epetra_Vector *v=NULL) const
Create vector orthog poly using f map and owned sg map.
bool supports_x
Whether we support x (and thus f and W)
Teuchos::Array< Teuchos::RCP< const Epetra_Map > > get_g_sg_base_maps() const
Get base maps of SG responses.
unsigned int num_W_blocks
Number of W stochastic blocks (may be smaller than num_sg_blocks)
Teuchos::RCP< EpetraExt::BlockVector > import_residual(const Epetra_Vector &f) const
Import parallel residual vector.
Teuchos::RCP< const Epetra_Map > x_map
Underlying unknown map.
Teuchos::RCP< Epetra_Operator > create_DfDp_op(int i) const
Create SG operator representing df/dp.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > f_sg_blocks
f stochastic Galerkin components
Teuchos::RCP< const Stokhos::EpetraVectorOrthogPoly > get_p_sg_init(int l) const
Return initial SG parameters.
unsigned int num_p_blocks
Number of p stochastic blocks (may be smaller than num_sg_blocks)
int num_g
Number of response vectors of underlying model evaluator.
Teuchos::RCP< const Epetra_Map > sg_overlapped_f_map
Block SG overlapped residual map.
Factory for generating stochastic Galerkin preconditioners.
virtual Teuchos::RCP< Stokhos::SGOperator > build(const Teuchos::RCP< const EpetraExt::MultiComm > &sg_comm, const Teuchos::RCP< const Stokhos::OrthogPolyBasis< int, double > > &sg_basis, const Teuchos::RCP< const Stokhos::EpetraSparse3Tensor > &epetraCijk, const Teuchos::RCP< const Epetra_Map > &domain_base_map, const Teuchos::RCP< const Epetra_Map > &range_base_map, const Teuchos::RCP< const Epetra_Map > &domain_sg_map, const Teuchos::RCP< const Epetra_Map > &range_sg_map)
Build preconditioner operator.
Factory for generating stochastic Galerkin preconditioners.
ScalarType g(const Teuchos::Array< ScalarType > &x, const ScalarType &y)
ScalarType f(const Teuchos::Array< ScalarType > &x, double a, double b)