45 #include "Teuchos_TimeMonitor.hpp" 49 const Teuchos::RCP<const EpetraExt::MultiComm>& sg_comm_,
51 const Teuchos::RCP<const Stokhos::EpetraSparse3Tensor>& epetraCijk_,
52 const Teuchos::RCP<const Epetra_CrsGraph>& base_graph,
53 const Teuchos::RCP<const Epetra_Map>& domain_sg_map_,
54 const Teuchos::RCP<const Epetra_Map>& range_sg_map_,
55 const Teuchos::RCP<Teuchos::ParameterList>& params) :
57 *(epetraCijk_->getStochasticGraph()),
61 epetraCijk(epetraCijk_),
62 domain_sg_map(domain_sg_map_),
63 range_sg_map(range_sg_map_),
64 Cijk(epetraCijk->getParallelCijk()),
68 only_use_linear(false)
70 scale_op = params->get(
"Scale Operator by Inverse Basis Norms",
true);
83 const Teuchos::RCP<Stokhos::EpetraOperatorOrthogPoly >& ops)
85 #ifdef STOKHOS_TEUCHOS_TIME_MONITOR 86 TEUCHOS_FUNC_TIME_MONITOR(
"Stokhos: SG Fully Assembled Operator Assembly");
97 if (!include_mean && index(k_begin) == 0)
99 if (only_use_linear) {
100 int dim = sg_basis->dimension();
101 k_end = Cijk->find_k(dim+1);
105 const Teuchos::Array<double>& norms = sg_basis->norm_squared();
108 Teuchos::RCP<Epetra_RowMatrix> block =
112 j_it != Cijk->j_end(k_it); ++j_it) {
113 int j = epetraCijk->GCID(index(j_it));
115 i_it != Cijk->i_end(j_it); ++i_it) {
116 int i = epetraCijk->GRID(index(i_it));
117 double c = value(i_it);
120 this->SumIntoGlobalBlock(c, *block, i,
j);
125 this->FillComplete(*domain_sg_map, *range_sg_map);
128 Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly >
135 Teuchos::RCP<const Stokhos::EpetraOperatorOrthogPoly >
146 #ifdef STOKHOS_TEUCHOS_TIME_MONITOR 147 TEUCHOS_FUNC_TIME_MONITOR(
"Stokhos: SG Operator Apply");
virtual void setupOperator(const Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > &poly)
Setup operator.
virtual Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > getSGPolynomial()
Get SG polynomial.
Bi-directional iterator for traversing a sparse array.
virtual int Apply(const Epetra_MultiVector &Input, Epetra_MultiVector &Result) const
Wrap Apply() to add a timer.
bool only_use_linear
Flag indicating whether to only use linear terms.
bool scale_op
Flag indicating whether operator be scaled with <^2>
bool include_mean
Flag indicating whether to include mean term.
CRS matrix of dense blocks.
int Apply(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
FullyAssembledOperator(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_CrsGraph > &base_graph, const Teuchos::RCP< const Epetra_Map > &domain_sg_map, const Teuchos::RCP< const Epetra_Map > &range_sg_map, const Teuchos::RCP< Teuchos::ParameterList > ¶ms)
Constructor.
virtual ~FullyAssembledOperator()
Destructor.