42 #ifndef TWOD_DIFFUSION_ME_HPP 43 #define TWOD_DIFFUSION_ME_HPP 45 #include "Teuchos_RCP.hpp" 46 #include "Teuchos_Array.hpp" 47 #include "Teuchos_ParameterList.hpp" 49 #include "EpetraExt_ModelEvaluator.h" 51 #include "Epetra_Map.h" 52 #include "Epetra_LocalMap.h" 53 #include "Epetra_Import.h" 54 #include "Epetra_CrsGraph.h" 55 #include "Epetra_CrsMatrix.h" 67 const Teuchos::RCP<const Epetra_Comm>& comm,
int n,
int d,
68 double s = 0.1,
double mu = 0.2,
73 const Teuchos::RCP<Teuchos::ParameterList>&
precParams = Teuchos::null);
82 Teuchos::RCP<const Epetra_Map>
get_x_map()
const;
85 Teuchos::RCP<const Epetra_Map>
get_f_map()
const;
88 Teuchos::RCP<const Epetra_Map>
get_p_map(
int l)
const;
91 Teuchos::RCP<const Epetra_Map>
get_g_map(
int j)
const;
94 Teuchos::RCP<const Teuchos::Array<std::string> >
98 Teuchos::RCP<const Epetra_Vector>
get_x_init()
const;
101 Teuchos::RCP<const Epetra_Vector>
get_p_init(
int l)
const;
104 Teuchos::RCP<Epetra_Operator>
create_W()
const;
107 Teuchos::RCP<EpetraExt::ModelEvaluator::Preconditioner>
create_WPrec()
const;
116 void evalModel(
const InArgs& inArgs,
const OutArgs& outArgs)
const;
121 Teuchos::RCP<Epetra_CrsMatrix>
get_mean()
const {
return A_k[0]; }
126 template <
typename FuncT>
138 Teuchos::Array<MeshPoint>
mesh;
141 Teuchos::RCP<const Stokhos::OrthogPolyBasis<int,double> >
basis;
167 Teuchos::RCP< Teuchos::Array<std::string> >
p_names;
170 Teuchos::RCP<Epetra_CrsGraph>
graph;
173 Teuchos::Array<Teuchos::RCP<Epetra_CrsMatrix> >
A_k;
176 Teuchos::RCP<Epetra_Vector>
b;
182 Teuchos::RCP<Epetra_CrsMatrix>
A;
185 mutable Teuchos::Array<double>
point;
192 template <
typename FuncT>
197 int NumMyElements =
x_map->NumMyElements();
198 int *MyGlobalElements =
x_map->MyGlobalElements();
203 for (
int k=0; k<sz; k++) {
204 A_k[k] = Teuchos::rcp(
new Epetra_CrsMatrix(Copy, *
graph));
205 for(
int i=0 ; i<NumMyElements; ++i ) {
208 int global_idx = MyGlobalElements[i];
209 if (
mesh[global_idx].boundary) {
214 A_k[k]->ReplaceGlobalValues(global_idx, 1, &
val, &global_idx);
218 -func(
mesh[global_idx].
x,
mesh[global_idx].
y-
h/2.0, k)/h2;
220 -func(
mesh[global_idx].
x-
h/2.0,
mesh[global_idx].
y, k)/h2;
222 -func(
mesh[global_idx].
x+
h/2.0,
mesh[global_idx].
y, k)/h2;
224 -func(
mesh[global_idx].
x,
mesh[global_idx].
y+
h/2.0, k)/h2;
227 val = -(a_down + a_left + a_right + a_up);
228 A_k[k]->ReplaceGlobalValues(global_idx, 1, &
val, &global_idx);
232 A_k[k]->ReplaceGlobalValues(global_idx, 1, &a_down,
233 &
mesh[global_idx].down);
237 A_k[k]->ReplaceGlobalValues(global_idx, 1, &a_left,
238 &
mesh[global_idx].left);
242 A_k[k]->ReplaceGlobalValues(global_idx, 1, &a_right,
243 &
mesh[global_idx].right);
247 A_k[k]->ReplaceGlobalValues(global_idx, 1, &a_up,
248 &
mesh[global_idx].up);
251 A_k[k]->FillComplete();
252 A_k[k]->OptimizeStorage();
256 #endif // TWOD_DIFFUSION_ME_HPP Teuchos::Array< int > bcIndices
Teuchos::RCP< Epetra_CrsMatrix > get_mean() const
Get mean matrix.
Teuchos::Array< Teuchos::RCP< Epetra_CrsMatrix > > A_k
KL coefficients of operator.
Teuchos::RCP< Epetra_Vector > x_init
Initial guess.
Teuchos::RCP< Epetra_Map > p_map
Parameter vector map.
Teuchos::RCP< const Epetra_Map > get_f_map() const
Return residual vector map.
Teuchos::RCP< const Epetra_Vector > get_p_init(int l) const
Return initial parameters.
Teuchos::RCP< Stokhos::AbstractPreconditionerFactory > precFactory
Teuchos::RCP< Epetra_CrsGraph > graph
Jacobian graph.
Teuchos::RCP< const Epetra_Map > get_g_map(int j) const
Return response function map.
const IndexType const IndexType const IndexType const IndexType const ValueType const ValueType * x
ModelEvaluator for a linear 2-D diffusion problem.
Teuchos::Array< MeshPoint > mesh
Teuchos::RCP< Teuchos::ParameterList > precParams
Teuchos::Array< double > basis_vals
Array to store values of basis at a point.
Teuchos::RCP< const Epetra_Map > get_x_map() const
Return solution vector map.
~twoD_diffusion_ME()
Destructor.
Teuchos::RCP< Epetra_Map > x_map
Solution vector map.
Teuchos::RCP< Epetra_CrsMatrix > A
Matrix to store deterministic operator.
twoD_diffusion_ME(const Teuchos::RCP< const Epetra_Comm > &comm, int n, int d, double s=0.1, double mu=0.2, const Teuchos::RCP< const Stokhos::OrthogPolyBasis< int, double > > &basis=Teuchos::null, bool log_normal=false, bool eliminate_bcs=false, const Teuchos::RCP< Teuchos::ParameterList > &precParams=Teuchos::null)
Constructor.
void evalModel(const InArgs &inArgs, const OutArgs &outArgs) const
Evaluate model on InArgs.
Teuchos::RCP< Epetra_Vector > p_init
Initial parameters.
Teuchos::RCP< Epetra_Operator > create_W() const
Create W = alpha*M + beta*J matrix.
InArgs createInArgs() const
Create InArgs.
Teuchos::RCP< const Epetra_Vector > get_x_init() const
Return initial solution.
Teuchos::RCP< const Teuchos::Array< std::string > > get_p_names(int l) const
Return array of parameter names.
Teuchos::RCP< const Epetra_Map > get_p_map(int l) const
Return parameter vector map.
Teuchos::RCP< Epetra_Vector > b
Deterministic RHS.
Teuchos::RCP< Teuchos::Array< std::string > > p_names
Parameter names.
OutArgs createOutArgs() const
Create OutArgs.
Teuchos::Array< double > point
Array to store a point for basis evaluation.
void fillMatrices(const FuncT &func, int sz)
Fill coefficient matrix given function to evaluate diffusion coefficient.
Teuchos::RCP< Epetra_Map > g_map
Response vector map.
const IndexType const IndexType const IndexType const IndexType const ValueType const ValueType ValueType * y
Teuchos::RCP< Epetra_Import > importer
Importer to overlapped distribution.
Teuchos::RCP< EpetraExt::ModelEvaluator::Preconditioner > create_WPrec() const
Create preconditioner for W.
Teuchos::RCP< const Stokhos::OrthogPolyBasis< int, double > > basis
Teuchos::Array< Teuchos::RCP< Epetra_Vector > > sg_kx_vec_all
Vectors to store matrix-vector products in SG residual calculation.