Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
Stokhos_SGModelEvaluator_Adaptive.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Stokhos Package
5 // Copyright (2009) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Eric T. Phipps (etphipp@sandia.gov).
38 //
39 // ***********************************************************************
40 // @HEADER
41 
42 #ifndef STOKHOS_SGMODELEVALUATOR_ADAPTIVE_HPP
43 #define STOKHOS_SGMODELEVALUATOR_ADAPTIVE_HPP
44 
45 #include <vector>
46 
47 #include "EpetraExt_ModelEvaluator.h"
48 #include "EpetraExt_MultiComm.h"
49 #include "EpetraExt_BlockVector.h"
50 
51 #include "Teuchos_RCP.hpp"
52 #include "Teuchos_Array.hpp"
53 #include "Teuchos_ParameterList.hpp"
54 #include "Stokhos_ParallelData.hpp"
60 #include "Stokhos_SGOperator.hpp"
63 
64 namespace Stokhos {
65 
67 
81  class SGModelEvaluator_Adaptive : public EpetraExt::ModelEvaluator {
82  public:
83 
84  // Constructor
86  const Teuchos::RCP<EpetraExt::ModelEvaluator>& me_,
87  const Teuchos::RCP<Stokhos::AdaptivityManager> & am,
88  const Teuchos::RCP<const Stokhos::Quadrature<int,double> >& sg_quad_,
89  const Teuchos::RCP<Stokhos::OrthogPolyExpansion<int,double> >& sg_exp_,
90  const Teuchos::RCP<const Stokhos::ParallelData>& sg_parallel_data_,
91  bool onlyUseLinear_,int kExpOrder_,
92  const Teuchos::RCP<Teuchos::ParameterList>& params_);
93 
94  // Constructor
96  const Teuchos::RCP<EpetraExt::ModelEvaluator>& me,
97  const Teuchos::RCP<const Stokhos::OrthogPolyBasis<int,double> >& sg_master_basis,
98  const std::vector<Teuchos::RCP<const Stokhos::ProductBasis<int,double> > > & sg_row_dof_basis,
99  const Teuchos::RCP<const Stokhos::Quadrature<int,double> >& sg_quad,
100  const Teuchos::RCP<Stokhos::OrthogPolyExpansion<int,double> >& sg_exp,
101  const Teuchos::RCP<const Stokhos::ParallelData>& sg_parallel_data,
102  bool onlyUseLinear,int kExpOrder,
103  const Teuchos::RCP<Teuchos::ParameterList>& params,
104  bool scaleOP = true);
105 
108 
109  // inputs
111 
113  Teuchos::RCP<const Epetra_Map> get_x_map() const;
114 
116  Teuchos::RCP<const Epetra_Map> get_p_map(int l) const;
117 
119  Teuchos::RCP<const Teuchos::Array<std::string> >
120  get_p_names(int l) const;
121 
123  Teuchos::RCP<const Epetra_Vector> get_x_init() const;
124 
126  Teuchos::RCP<const Epetra_Vector> get_p_init(int l) const;
127 
128  // outputs
130 
132  Teuchos::RCP<const Epetra_Map> get_f_map() const;
133 
135  Teuchos::RCP<const Epetra_Map> get_g_map(int l) const;
136 
138  Teuchos::RCP<Epetra_Operator> create_W() const;
139 
140  // ????
142 
144  InArgs createInArgs() const;
145 
147  OutArgs createOutArgs() const;
148 
150  void evalModel(const InArgs& inArgs, const OutArgs& outArgs) const;
151 
153 
155  void set_x_sg_init(const Stokhos::EpetraVectorOrthogPoly& x_sg_in);
156 
158  Teuchos::RCP<const Stokhos::EpetraVectorOrthogPoly> get_x_sg_init() const;
159 
161  void set_p_sg_init(int i, const Stokhos::EpetraVectorOrthogPoly& p_sg_in);
162 
164  Teuchos::RCP<const Stokhos::EpetraVectorOrthogPoly> get_p_sg_init(int l) const;
165 
167 
170  Teuchos::Array<int> get_p_sg_map_indices() const;
171 
173 
176  Teuchos::Array<int> get_g_sg_map_indices() const;
177 
179  Teuchos::Array< Teuchos::RCP<const Epetra_Map> > get_g_sg_base_maps() const;
180 
182  Teuchos::RCP<const Epetra_BlockMap> get_overlap_stochastic_map() const;
183 
185  Teuchos::RCP<const Epetra_BlockMap> get_x_sg_overlap_map() const;
186 
188  Teuchos::RCP<const Epetra_Import> get_x_sg_importer() const;
189 
191  Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
192  create_x_sg() const;
193 
195  Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
196  create_x_sg_overlap() const;
197 
199  Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly>
200  create_x_mv_sg(int num_vecs) const;
201 
203  Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly>
204  create_x_mv_sg_overlap(int num_vecs) const;
205 
207  Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
208  create_p_sg(int l,Epetra_DataAccess CV=Copy,const Epetra_Vector * v=0) const;
209 
211  Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
212  create_f_sg() const;
213 
215  Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
216  create_f_sg_overlap() const;
217 
219  Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly>
220  create_f_mv_sg(int num_vecs) const;
221 
223  Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly>
224  create_f_mv_sg_overlap(int num_vecs) const;
225 
227  Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
228  create_g_sg(int l, Epetra_DataAccess CV = Copy,
229  const Epetra_Vector* v = NULL) const;
230 
232  Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly>
233  create_g_mv_sg(int l, int num_vecs, Epetra_DataAccess CV = Copy,
234  const Epetra_MultiVector* v = NULL) const;
235 
236  Teuchos::RCP<const Stokhos::AdaptivityManager> getAdaptivityManager() const
237  { return adaptMngr; }
238 
239  protected:
240 
242  Teuchos::RCP<EpetraExt::ModelEvaluator> me;
243 
245  Teuchos::RCP<const Stokhos::OrthogPolyBasis<int, double> > sg_basis;
246 
247  std::vector<Teuchos::RCP<const Stokhos::ProductBasis<int,double> > > sg_row_dof_basis;
248 
250  Teuchos::RCP<const Stokhos::Quadrature<int,double> > sg_quad;
251 
253  Teuchos::RCP<Stokhos::OrthogPolyExpansion<int,double> > sg_exp;
254 
256  Teuchos::RCP<Teuchos::ParameterList> params;
257 
259  unsigned int num_sg_blocks;
260 
262  unsigned int num_W_blocks;
263 
265  unsigned int num_p_blocks;
266 
269 
271  Teuchos::RCP<const Epetra_Map> x_map;
272 
274  Teuchos::RCP<const Epetra_Map> f_map;
275 
277  Teuchos::RCP<const Stokhos::ParallelData> sg_parallel_data;
278 
280  Teuchos::RCP<const EpetraExt::MultiComm> sg_comm;
281 
283  Teuchos::RCP<const Stokhos::EpetraSparse3Tensor> epetraCijk;
284 
286  Teuchos::RCP<const Stokhos::EpetraSparse3Tensor> serialCijk;
287 
288  Teuchos::RCP<const Stokhos::Sparse3Tensor<int,double> > Cijk;
289 
291  Teuchos::RCP<const Epetra_BlockMap> stoch_row_map;
292 
294  Teuchos::RCP<const Epetra_BlockMap> overlapped_stoch_row_map;
295 
297  Teuchos::RCP<const Epetra_BlockMap> overlapped_stoch_p_map;
298 
300  Teuchos::RCP<const Epetra_Map> adapted_x_map;
301 
303  Teuchos::RCP<const Epetra_Map> adapted_overlapped_x_map;
304 
306  Teuchos::RCP<const Epetra_Map> adapted_f_map;
307 
309  Teuchos::RCP<const Epetra_Map> adapted_overlapped_f_map;
310 
312  Teuchos::RCP<Epetra_Import> adapted_overlapped_x_importer;
313 
315  Teuchos::RCP<Epetra_Export> adapted_overlapped_f_exporter;
316 
318  int num_p;
319 
321  int num_p_sg;
322 
324  Teuchos::Array<int> sg_p_index_map;
325 
327  Teuchos::Array< Teuchos::RCP<const Epetra_Map> > sg_p_map;
328 
330  Teuchos::Array< Teuchos::RCP< Teuchos::Array<std::string> > > sg_p_names;
331 
333  int num_g;
334 
336  int num_g_sg;
337 
339  Teuchos::Array<int> sg_g_index_map;
340 
342  Teuchos::Array< Teuchos::RCP<const Epetra_Map> > sg_g_map;
343 
345  Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > x_dot_sg_blocks;
346 
348  Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > x_sg_blocks;
349 
351  mutable Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > f_sg_blocks;
352 
354  mutable Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > W_sg_blocks;
355 
356  mutable Teuchos::Array< Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > > dfdp_sg_blocks;
357 
359  mutable Teuchos::Array< Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > > dgdx_dot_sg_blocks;
360 
362  mutable Teuchos::Array< Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > > dgdx_sg_blocks;
363 
365  Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> sg_x_init;
366 
368  Teuchos::Array< Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> > sg_p_init;
369 
372 
375 
377  mutable Teuchos::RCP<Epetra_CrsMatrix> my_W;
378 
380  mutable Teuchos::RCP<Epetra_Vector> my_x;
381 
382  bool scaleOP;
383 
384  mutable Teuchos::RCP<Stokhos::AdaptivityManager> adaptMngr;
385  };
386 
387 }
388 
389 #endif
void evalModel(const InArgs &inArgs, const OutArgs &outArgs) const
Evaluate model on InArgs.
Teuchos::RCP< const Epetra_Map > get_g_map(int l) const
Return response map.
Teuchos::Array< int > get_p_sg_map_indices() const
Get indices of SG parameters.
bool eval_W_with_f
Whether to always evaluate W with f.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > x_sg_blocks
x stochastic Galerkin components
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > create_x_sg() const
Create vector orthog poly using x map and owned sg map.
int num_g_sg
Number of stochastic response vectors.
Teuchos::RCP< Epetra_Operator > create_W() const
Create W = alpha*M + beta*J matrix.
Teuchos::Array< int > get_g_sg_map_indices() const
Get indices of SG responses.
Teuchos::RCP< const Epetra_Map > get_x_map() const
Return solution vector map.
Teuchos::RCP< const Epetra_BlockMap > stoch_row_map
Map for stochastic blocks.
Teuchos::RCP< const Stokhos::EpetraSparse3Tensor > serialCijk
Serial Epetra Cijk for dgdx*.
Teuchos::RCP< Stokhos::OrthogPolyExpansion< int, double > > sg_exp
Stochastic Galerkin expansion.
Teuchos::RCP< const Stokhos::Sparse3Tensor< int, double > > Cijk
Teuchos::Array< Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > > dgdx_sg_blocks
dg/dx stochastic Galerkin components
Teuchos::RCP< const Epetra_Map > get_f_map() const
Return residual vector map.
Teuchos::RCP< Epetra_Export > adapted_overlapped_f_exporter
Exporter from SG-overlapped to SG maps.
Teuchos::Array< Teuchos::RCP< Teuchos::Array< std::string > > > sg_p_names
SG coefficient parameter names.
Teuchos::Array< Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > > dgdx_dot_sg_blocks
dg/dxdot stochastic Galerkin components
void set_p_sg_init(int i, const Stokhos::EpetraVectorOrthogPoly &p_sg_in)
Set initial parameter polynomial.
Teuchos::RCP< const Epetra_Map > x_map
Underlying unknown map.
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< const Stokhos::ParallelData > sg_parallel_data
Parallel SG data.
Teuchos::RCP< const Epetra_BlockMap > get_overlap_stochastic_map() const
Return overlap stochastic map.
Teuchos::RCP< const Epetra_Map > adapted_overlapped_x_map
Block SG overlapped unknown map.
int num_g
Number of response vectors of underlying model evaluator.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > sg_x_init
SG initial x.
Teuchos::Array< Teuchos::RCP< const Epetra_Map > > sg_g_map
Block SG response map.
Nonlinear, stochastic Galerkin ModelEvaluator that constructs an adapted Jacobian.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > create_f_sg_overlap() 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 Stokhos::EpetraSparse3Tensor > epetraCijk
Epetra Cijk.
Teuchos::RCP< const Teuchos::Array< std::string > > get_p_names(int l) const
Return array of parameter names.
Abstract base class for orthogonal polynomial-based expansions.
unsigned int num_W_blocks
Number of W stochastic blocks (may be smaller than num_sg_blocks)
std::vector< Teuchos::RCP< const Stokhos::ProductBasis< int, double > > > sg_row_dof_basis
Teuchos::RCP< EpetraExt::ModelEvaluator > me
Underlying model evaluator.
Teuchos::Array< Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > > dfdp_sg_blocks
Teuchos::RCP< Epetra_CrsMatrix > my_W
W pointer for evaluating W with f.
Teuchos::RCP< const Stokhos::EpetraVectorOrthogPoly > get_p_sg_init(int l) const
Return initial SG parameters.
bool supports_x
Whether we support x (and thus f and W)
Abstract base class for quadrature methods.
int num_p
Number of parameter vectors of underlying model evaluator.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > create_p_sg(int l, Epetra_DataAccess CV=Copy, const Epetra_Vector *v=0) const
Create vector orthog poly using p map.
A container class storing an orthogonal polynomial whose coefficients are vectors, operators, or in general any type that would have an expensive copy constructor.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > create_x_sg_overlap() const
Create vector orthog poly using x map and overlap sg map.
Top-level namespace for Stokhos classes and functions.
Teuchos::RCP< const Stokhos::Quadrature< int, double > > sg_quad
Stochastic Galerkin quadrature.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > x_dot_sg_blocks
x_dot stochastic Galerkin components
Teuchos::RCP< const Stokhos::AdaptivityManager > getAdaptivityManager() const
Teuchos::RCP< const Epetra_Map > adapted_overlapped_f_map
Block SG overlapped residual map.
Teuchos::RCP< const Epetra_BlockMap > overlapped_stoch_row_map
Overlapped map for stochastic blocks (local map)
Teuchos::Array< Teuchos::RCP< const Epetra_Map > > get_g_sg_base_maps() const
Get base maps of SG responses.
Abstract base class for multivariate orthogonal polynomials generated from tensor products of univari...
Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > create_x_mv_sg(int num_vecs) const
Create vector orthog poly using x map and owned sg map.
Teuchos::RCP< const Epetra_Map > adapted_f_map
Block SG residual map.
Teuchos::RCP< const Epetra_Map > adapted_x_map
Block SG unknown map.
Teuchos::RCP< const Stokhos::EpetraVectorOrthogPoly > get_x_sg_init() const
Return initial SG x.
Teuchos::Array< Teuchos::RCP< const Epetra_Map > > sg_p_map
Block SG parameter map.
Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > create_f_mv_sg_overlap(int num_vecs) const
Create multi-vector orthog poly using f map and overlap sg map.
Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > W_sg_blocks
W stochastic Galerkin components.
Teuchos::RCP< Epetra_Import > adapted_overlapped_x_importer
Importer from SG to SG-overlapped maps.
Teuchos::RCP< const Epetra_BlockMap > get_x_sg_overlap_map() const
Return x sg overlap map.
Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > create_x_mv_sg_overlap(int num_vecs) const
Create vector orthog poly using x map and overlap sg map.
Teuchos::RCP< const Epetra_Vector > get_p_init(int l) const
Return initial parameters.
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< Stokhos::EpetraMultiVectorOrthogPoly > create_f_mv_sg(int num_vecs) const
Create multi-vector orthog poly using f map and owned sg map.
Teuchos::Array< int > sg_p_index_map
Index map between block-p and p_sg maps.
Teuchos::RCP< const Epetra_Vector > get_x_init() const
Return initial solution.
int num_p_sg
Number of stochastic parameter vectors.
void set_x_sg_init(const Stokhos::EpetraVectorOrthogPoly &x_sg_in)
Set initial solution polynomial.
unsigned int num_sg_blocks
Number of stochastic blocks.
Teuchos::Array< int > sg_g_index_map
Index map between block-g and g_sg maps.
Teuchos::RCP< const EpetraExt::MultiComm > sg_comm
Parallel SG communicator.
Teuchos::RCP< const Epetra_Import > get_x_sg_importer() const
Return x sg importer.
Teuchos::RCP< Epetra_Vector > my_x
x pointer for evaluating preconditioner
Teuchos::RCP< const Stokhos::OrthogPolyBasis< int, double > > sg_basis
Stochastic Galerkin basis.
Teuchos::RCP< Teuchos::ParameterList > params
Algorithmic parameters.
unsigned int num_p_blocks
Number of p stochastic blocks (may be smaller than num_sg_blocks)
Teuchos::RCP< Stokhos::AdaptivityManager > adaptMngr
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > create_f_sg() const
Create vector orthog poly using f map and owned sg map.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > f_sg_blocks
f stochastic Galerkin components
Teuchos::RCP< const Epetra_BlockMap > overlapped_stoch_p_map
Overlapped map for p stochastic blocks (local map)
SGModelEvaluator_Adaptive(const Teuchos::RCP< EpetraExt::ModelEvaluator > &me_, const Teuchos::RCP< Stokhos::AdaptivityManager > &am, 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_, bool onlyUseLinear_, int kExpOrder_, const Teuchos::RCP< Teuchos::ParameterList > &params_)
Teuchos::Array< Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > > sg_p_init
SG initial p.
Teuchos::RCP< const Epetra_Map > f_map
Underlying residual map.