Stokhos Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
Stokhos_SGModelEvaluator_Adaptive.cpp
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
43
44#include <algorithm>
45#include "Teuchos_Assert.hpp"
46#include "EpetraExt_BlockUtility.h"
47#include "EpetraExt_BlockMultiVector.h"
52#include "Epetra_LocalMap.h"
53#include "Epetra_Export.h"
54#include "Epetra_Import.h"
55
57 const Teuchos::RCP<EpetraExt::ModelEvaluator>& me_,
58 const Teuchos::RCP<Stokhos::AdaptivityManager> & am,
59 const Teuchos::RCP<const Stokhos::Quadrature<int,double> >& sg_quad_,
60 const Teuchos::RCP<Stokhos::OrthogPolyExpansion<int,double> >& sg_exp_,
61 const Teuchos::RCP<const Stokhos::ParallelData>& sg_parallel_data_,
62 bool onlyUseLinear_,int kExpOrder_,
63 const Teuchos::RCP<Teuchos::ParameterList>& params_)
64 : me(me_),
65 sg_basis(am->getMasterStochasticBasis()),
66 sg_row_dof_basis(am->getRowStochasticBasis()),
67 sg_quad(sg_quad_),
68 sg_exp(sg_exp_),
69 params(params_),
70 num_sg_blocks(sg_basis->size()),
71 num_W_blocks(sg_basis->size()),
72 num_p_blocks(sg_basis->size()),
73 supports_x(false),
74 x_map(me->get_x_map()),
75 f_map(me->get_f_map()),
76 sg_parallel_data(sg_parallel_data_),
77 sg_comm(sg_parallel_data->getMultiComm()),
78 epetraCijk(sg_parallel_data->getEpetraCijk()),
79 serialCijk(),
80 Cijk(epetraCijk->getParallelCijk()),
81 num_p(0),
82 num_p_sg(0),
83 sg_p_index_map(),
84 sg_p_map(),
85 sg_p_names(),
86 num_g(0),
87 num_g_sg(0),
88 sg_g_index_map(),
89 sg_g_map(),
90 x_dot_sg_blocks(),
91 x_sg_blocks(),
92 f_sg_blocks(),
93 W_sg_blocks(),
94 dgdx_dot_sg_blocks(),
95 dgdx_sg_blocks(),
96 sg_x_init(),
97 sg_p_init(),
98 eval_W_with_f(false),
99 kExpOrder(kExpOrder_),
100 onlyUseLinear(onlyUseLinear_),
101 scaleOP(am->isScaled()),
102 adaptMngr(am)
103{
104 if (x_map != Teuchos::null)
105 supports_x = true;
106
108 Teuchos::rcp(new Epetra_LocalMap(static_cast<int>(num_sg_blocks), 0, *(sg_parallel_data->getStochasticComm())));
110
111 if (epetraCijk != Teuchos::null)
112 stoch_row_map = epetraCijk->getStochasticRowMap();
113
114 if (supports_x) {
115
116 // Create block SG x and f maps
117 sg_x_map =
118 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
121 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
123
124 sg_f_map =
125 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
128 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
130
131 // Create interlace SG x and f maps
132 adapted_x_map = adaptMngr->getAdaptedMap();
134
135 adapted_f_map = adapted_x_map; // these are the same! No overlap possible in refinement!
137
138 // Create importer/exporter from/to overlapped distribution
140 Teuchos::rcp(new Epetra_Import(*adapted_overlapped_x_map, *get_x_map()));
143
144 // now we create the underlying Epetra block vectors
145 // that will be used by the model evaluator to construct
146 // the solution of the deterministic problem.
148
149 // Create vector blocks
153
154 // Create default sg_x_init
156 if (sg_x_init->myGID(0))
157 (*sg_x_init)[0] = *(me->get_x_init());
158
159 // Preconditioner needs an x: This is adapted
160 my_x = Teuchos::rcp(new Epetra_Vector(*get_x_map()));
161
162 // setup storage for W, these are blocked in Stokhos
163 // format
165
166 // Determine W expansion type
167 std::string W_expansion_type =
168 params->get("Jacobian Expansion Type", "Full");
169 if (W_expansion_type == "Linear")
170 num_W_blocks = sg_basis->dimension() + 1;
171 else
173
174 Teuchos::RCP<Epetra_BlockMap> W_overlap_map =
175 Teuchos::rcp(new Epetra_LocalMap(
176 static_cast<int>(num_W_blocks), 0,
177 *(sg_parallel_data->getStochasticComm())));
179 Teuchos::rcp(new Stokhos::EpetraOperatorOrthogPoly(
180 sg_basis, W_overlap_map, x_map, f_map, adapted_f_map,
181 sg_comm));
182 for (unsigned int i=0; i<num_W_blocks; i++)
183 W_sg_blocks->setCoeffPtr(i, me->create_W()); // allocate a bunch of matrices
184
185 eval_W_with_f = params->get("Evaluate W with F", false);
186 }
187
188 // Parameters -- The idea here is to add new parameter vectors
189 // for the stochastic Galerkin components of the parameters
190
191 InArgs me_inargs = me->createInArgs();
192 OutArgs me_outargs = me->createOutArgs();
193 num_p = me_inargs.Np();
194
195 // Get the p_sg's supported and build index map
196 for (int i=0; i<num_p; i++) {
197 if (me_inargs.supports(IN_ARG_p_sg, i))
198 sg_p_index_map.push_back(i);
199 }
200 num_p_sg = sg_p_index_map.size();
201
202 sg_p_map.resize(num_p_sg);
203 sg_p_names.resize(num_p_sg);
204 sg_p_init.resize(num_p_sg);
205
206 // Determine parameter expansion type
207 std::string p_expansion_type =
208 params->get("Parameter Expansion Type", "Full");
209 if (p_expansion_type == "Linear")
210 num_p_blocks = sg_basis->dimension() + 1;
211 else
213
214 // Create parameter maps, names, and initial values
216 Teuchos::rcp(new Epetra_LocalMap(
217 static_cast<int>(num_p_blocks), 0,
218 *(sg_parallel_data->getStochasticComm())));
219 for (int i=0; i<num_p_sg; i++) {
220 Teuchos::RCP<const Epetra_Map> p_map = me->get_p_map(sg_p_index_map[i]);
221 sg_p_map[i] =
222 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
223 *p_map, *overlapped_stoch_p_map, *sg_comm));
224
225 Teuchos::RCP<const Teuchos::Array<std::string> > p_names =
226 me->get_p_names(sg_p_index_map[i]);
227 if (p_names != Teuchos::null) {
228 sg_p_names[i] =
229 Teuchos::rcp(new Teuchos::Array<std::string>(num_sg_blocks*(p_names->size())));
230 for (int j=0; j<p_names->size(); j++) {
231 std::stringstream ss;
232 ss << (*p_names)[j] << " -- SG Coefficient " << i;
233 (*sg_p_names[i])[j] = ss.str();
234 }
235 }
236
237 // Create default sg_p_init
239 sg_p_init[i]->init(0.0);
240 }
241
242 // Responses -- The idea here is to add new parameter vectors
243 // for the stochastic Galerkin components of the respones
244 num_g = me_outargs.Ng();
245
246 // Get the g_sg's supported and build index map
247 for (int i=0; i<num_g; i++) {
248 if (me_outargs.supports(OUT_ARG_g_sg, i))
249 sg_g_index_map.push_back(i);
250 }
251 num_g_sg = sg_g_index_map.size();
252
253 sg_g_map.resize(num_g_sg);
255 dgdx_sg_blocks.resize(num_g_sg);
256
257 // Create response maps
258 for (int i=0; i<num_g_sg; i++) {
259 Teuchos::RCP<const Epetra_Map> g_map = me->get_g_map(sg_g_index_map[i]);
260 sg_g_map[i] =
261 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
263
264 // Create dg/dxdot, dg/dx SG blocks
265 if (supports_x) {
269 dgdx_sg_blocks[i] =
272 }
273 }
274
275 // We don't support parallel for dgdx yet, so build a new EpetraCijk
276 if (supports_x) {
277 serialCijk =
278 Teuchos::rcp(new Stokhos::EpetraSparse3Tensor(sg_basis,
279 epetraCijk->getCijk(),
280 sg_comm,
282 }
283
284}
285
287 const Teuchos::RCP<EpetraExt::ModelEvaluator>& me_,
288 const Teuchos::RCP<const Stokhos::OrthogPolyBasis<int,double> >& sg_basis_,
289 const std::vector<Teuchos::RCP<const Stokhos::ProductBasis<int,double> > > & sg_row_dof_basis_,
290 const Teuchos::RCP<const Stokhos::Quadrature<int,double> >& sg_quad_,
291 const Teuchos::RCP<Stokhos::OrthogPolyExpansion<int,double> >& sg_exp_,
292 const Teuchos::RCP<const Stokhos::ParallelData>& sg_parallel_data_,
293 bool onlyUseLinear_,int kExpOrder_,
294 const Teuchos::RCP<Teuchos::ParameterList>& params_,
295 bool scaleOP_)
296 : me(me_),
297 sg_basis(sg_basis_),
298 sg_row_dof_basis(sg_row_dof_basis_),
299 sg_quad(sg_quad_),
300 sg_exp(sg_exp_),
301 params(params_),
302 num_sg_blocks(sg_basis->size()),
303 num_W_blocks(sg_basis->size()),
304 num_p_blocks(sg_basis->size()),
305 supports_x(false),
306 x_map(me->get_x_map()),
307 f_map(me->get_f_map()),
308 sg_parallel_data(sg_parallel_data_),
309 sg_comm(sg_parallel_data->getMultiComm()),
310 epetraCijk(sg_parallel_data->getEpetraCijk()),
311 serialCijk(),
312 Cijk(epetraCijk->getParallelCijk()),
313 num_p(0),
314 num_p_sg(0),
315 sg_p_index_map(),
316 sg_p_map(),
317 sg_p_names(),
318 num_g(0),
319 num_g_sg(0),
320 sg_g_index_map(),
321 sg_g_map(),
322 x_dot_sg_blocks(),
323 x_sg_blocks(),
324 f_sg_blocks(),
325 W_sg_blocks(),
326 dgdx_dot_sg_blocks(),
327 dgdx_sg_blocks(),
328 sg_x_init(),
329 sg_p_init(),
330 eval_W_with_f(false),
331 kExpOrder(kExpOrder_),
332 onlyUseLinear(onlyUseLinear_),
333 scaleOP(scaleOP_)
334{
335 if (x_map != Teuchos::null)
336 supports_x = true;
337
338 adaptMngr = Teuchos::rcp(new Stokhos::AdaptivityManager(Teuchos::rcp_dynamic_cast<const Stokhos::ProductBasis<int,double> >(sg_basis),
339 sg_row_dof_basis,*sg_parallel_data->getStochasticComm(),scaleOP));
340
342 Teuchos::rcp(new Epetra_LocalMap(static_cast<int>(num_sg_blocks), 0, *(sg_parallel_data->getStochasticComm())));
344
345 if (epetraCijk != Teuchos::null)
346 stoch_row_map = epetraCijk->getStochasticRowMap();
347
348 if (supports_x) {
349
350 // Create block SG x and f maps
351 sg_x_map =
352 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
355 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
357
358 sg_f_map =
359 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
362 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
364
365 // Create interlace SG x and f maps
366 adapted_x_map = adaptMngr->getAdaptedMap();
368
369 adapted_f_map = adapted_x_map; // these are the same! No overlap possible in refinement!
371
372 // Create importer/exporter from/to overlapped distribution
374 Teuchos::rcp(new Epetra_Import(*adapted_overlapped_x_map, *get_x_map()));
377
378 // now we create the underlying Epetra block vectors
379 // that will be used by the model evaluator to construct
380 // the solution of the deterministic problem.
382
383 // Create vector blocks
387
388 // Create default sg_x_init
390 if (sg_x_init->myGID(0))
391 (*sg_x_init)[0] = *(me->get_x_init());
392
393 // Preconditioner needs an x: This is adapted
394 my_x = Teuchos::rcp(new Epetra_Vector(*get_x_map()));
395
396 // setup storage for W, these are blocked in Stokhos
397 // format
399
400 // Determine W expansion type
401 std::string W_expansion_type =
402 params->get("Jacobian Expansion Type", "Full");
403 if (W_expansion_type == "Linear")
404 num_W_blocks = sg_basis->dimension() + 1;
405 else
407
408 Teuchos::RCP<Epetra_BlockMap> W_overlap_map =
409 Teuchos::rcp(new Epetra_LocalMap(
410 static_cast<int>(num_W_blocks), 0,
411 *(sg_parallel_data->getStochasticComm())));
413 Teuchos::rcp(new Stokhos::EpetraOperatorOrthogPoly(
414 sg_basis, W_overlap_map, x_map, f_map, adapted_f_map,
415 sg_comm));
416 for (unsigned int i=0; i<num_W_blocks; i++)
417 W_sg_blocks->setCoeffPtr(i, me->create_W()); // allocate a bunch of matrices
418
419 eval_W_with_f = params->get("Evaluate W with F", false);
420 }
421
422 // Parameters -- The idea here is to add new parameter vectors
423 // for the stochastic Galerkin components of the parameters
424
425 InArgs me_inargs = me->createInArgs();
426 OutArgs me_outargs = me->createOutArgs();
427 num_p = me_inargs.Np();
428
429 // Get the p_sg's supported and build index map
430 for (int i=0; i<num_p; i++) {
431 if (me_inargs.supports(IN_ARG_p_sg, i))
432 sg_p_index_map.push_back(i);
433 }
434 num_p_sg = sg_p_index_map.size();
435
436 sg_p_map.resize(num_p_sg);
437 sg_p_names.resize(num_p_sg);
438 sg_p_init.resize(num_p_sg);
439
440 // Determine parameter expansion type
441 std::string p_expansion_type =
442 params->get("Parameter Expansion Type", "Full");
443 if (p_expansion_type == "Linear")
444 num_p_blocks = sg_basis->dimension() + 1;
445 else
447
448 // Create parameter maps, names, and initial values
450 Teuchos::rcp(new Epetra_LocalMap(
451 static_cast<int>(num_p_blocks), 0,
452 *(sg_parallel_data->getStochasticComm())));
453 for (int i=0; i<num_p_sg; i++) {
454 Teuchos::RCP<const Epetra_Map> p_map = me->get_p_map(sg_p_index_map[i]);
455 sg_p_map[i] =
456 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
457 *p_map, *overlapped_stoch_p_map, *sg_comm));
458
459 Teuchos::RCP<const Teuchos::Array<std::string> > p_names =
460 me->get_p_names(sg_p_index_map[i]);
461 if (p_names != Teuchos::null) {
462 sg_p_names[i] =
463 Teuchos::rcp(new Teuchos::Array<std::string>(num_sg_blocks*(p_names->size())));
464 for (int j=0; j<p_names->size(); j++) {
465 std::stringstream ss;
466 ss << (*p_names)[j] << " -- SG Coefficient " << i;
467 (*sg_p_names[i])[j] = ss.str();
468 }
469 }
470
471 // Create default sg_p_init
473 sg_p_init[i]->init(0.0);
474 }
475
476 // Responses -- The idea here is to add new parameter vectors
477 // for the stochastic Galerkin components of the respones
478 num_g = me_outargs.Ng();
479
480 // Get the g_sg's supported and build index map
481 for (int i=0; i<num_g; i++) {
482 if (me_outargs.supports(OUT_ARG_g_sg, i))
483 sg_g_index_map.push_back(i);
484 }
485 num_g_sg = sg_g_index_map.size();
486
487 sg_g_map.resize(num_g_sg);
489 dgdx_sg_blocks.resize(num_g_sg);
490
491 // Create response maps
492 for (int i=0; i<num_g_sg; i++) {
493 Teuchos::RCP<const Epetra_Map> g_map = me->get_g_map(sg_g_index_map[i]);
494 sg_g_map[i] =
495 Teuchos::rcp(EpetraExt::BlockUtility::GenerateBlockMap(
497
498 // Create dg/dxdot, dg/dx SG blocks
499 if (supports_x) {
503 dgdx_sg_blocks[i] =
506 }
507 }
508
509 // We don't support parallel for dgdx yet, so build a new EpetraCijk
510 if (supports_x) {
511 serialCijk =
512 Teuchos::rcp(new Stokhos::EpetraSparse3Tensor(sg_basis,
513 epetraCijk->getCijk(),
514 sg_comm,
516 }
517
518}
519
520// Overridden from EpetraExt::ModelEvaluator
521
522Teuchos::RCP<const Epetra_Map>
524{
525 return adapted_x_map;
526}
527
528Teuchos::RCP<const Epetra_Map>
530{
531 return adapted_f_map;
532}
533
534Teuchos::RCP<const Epetra_Map>
536{
537 TEUCHOS_TEST_FOR_EXCEPTION(l < 0 || l >= num_p + num_p_sg, std::logic_error,
538 "Error! Invalid p map index " << l);
539 if (l < num_p)
540 return me->get_p_map(l);
541 else
542 return sg_p_map[l-num_p];
543
544 return Teuchos::null;
545}
546
547Teuchos::RCP<const Epetra_Map>
549{
550 TEUCHOS_TEST_FOR_EXCEPTION(l < 0 || l >= num_g_sg, std::logic_error,
551 "Error! Invalid g map index " << l);
552 return sg_g_map[l];
553}
554
555Teuchos::RCP<const Teuchos::Array<std::string> >
557{
558 TEUCHOS_TEST_FOR_EXCEPTION(l < 0 || l >= num_p + num_p_sg, std::logic_error,
559 "Error! Invalid p map index " << l);
560 if (l < num_p)
561 return me->get_p_names(l);
562 else
563 return sg_p_names[l-num_p];
564
565 return Teuchos::null;
566}
567
568Teuchos::RCP<const Epetra_Vector>
570{
571 // get stochastic galerking initial condition and write it out to x initial condition
572 Teuchos::RCP<Epetra_Vector> x_init = Teuchos::rcp(new Epetra_Vector(*get_x_map()));
573 adaptMngr->copyToAdaptiveVector(*get_x_sg_init(),*x_init);
574
575 return x_init;
576}
577
578Teuchos::RCP<const Epetra_Vector>
580{
581 TEUCHOS_TEST_FOR_EXCEPTION(l < 0 || l >= num_p + num_p_sg, std::logic_error,
582 "Error! Invalid p map index " << l);
583 if (l < num_p)
584 return me->get_p_init(l);
585 else
586 return sg_p_init[l-num_p]->getBlockVector();
587
588 return Teuchos::null;
589}
590
591Teuchos::RCP<Epetra_Operator>
593{
594 if (supports_x) {
595 Teuchos::RCP<Teuchos::ParameterList> sgOpParams =
596 Teuchos::rcp(&(params->sublist("SG Operator")), false);
597
598 Teuchos::RCP<Epetra_CrsMatrix> W_crs
599 = Teuchos::rcp_dynamic_cast<Epetra_CrsMatrix>(me->create_W(), true);
600 const Epetra_CrsGraph & W_graph = W_crs->Graph();
601
602 adaptMngr->setupWithGraph(W_graph,onlyUseLinear,kExpOrder);
603 my_W = adaptMngr->buildMatrixFromGraph();
604 adaptMngr->setupOperator(*my_W,*Cijk,*W_sg_blocks);
605
606 return my_W;
607 }
608
609 return Teuchos::null;
610}
611
612EpetraExt::ModelEvaluator::InArgs
614{
615 InArgsSetup inArgs;
616 InArgs me_inargs = me->createInArgs();
617
618 inArgs.setModelEvalDescription(this->description());
619 inArgs.set_Np(num_p + num_p_sg);
620 inArgs.setSupports(IN_ARG_x_dot, me_inargs.supports(IN_ARG_x_dot_sg));
621 inArgs.setSupports(IN_ARG_x, me_inargs.supports(IN_ARG_x_sg));
622 inArgs.setSupports(IN_ARG_t, me_inargs.supports(IN_ARG_t));
623 inArgs.setSupports(IN_ARG_alpha, me_inargs.supports(IN_ARG_alpha));
624 inArgs.setSupports(IN_ARG_beta, me_inargs.supports(IN_ARG_beta));
625 inArgs.setSupports(IN_ARG_sg_basis, me_inargs.supports(IN_ARG_sg_basis));
626 inArgs.setSupports(IN_ARG_sg_quadrature,
627 me_inargs.supports(IN_ARG_sg_quadrature));
628 inArgs.setSupports(IN_ARG_sg_expansion,
629 me_inargs.supports(IN_ARG_sg_expansion));
630
631 return inArgs;
632}
633
634EpetraExt::ModelEvaluator::OutArgs
636{
637 OutArgsSetup outArgs;
638 OutArgs me_outargs = me->createOutArgs();
639
640 outArgs.setModelEvalDescription(this->description());
641 outArgs.set_Np_Ng(num_p+num_p_sg, num_g_sg);
642
643 outArgs.setSupports(OUT_ARG_f, me_outargs.supports(OUT_ARG_f_sg));
644 outArgs.setSupports(OUT_ARG_W, me_outargs.supports(OUT_ARG_W_sg));
645 outArgs.setSupports(OUT_ARG_WPrec, false);
646 for (int j=0; j<num_p; j++)
647 outArgs.setSupports(OUT_ARG_DfDp, j,
648 me_outargs.supports(OUT_ARG_DfDp_sg, j));
649
650 for (int i=0; i<num_g_sg; i++) {
651 int ii = sg_g_index_map[i];
652// if (!me_outargs.supports(OUT_ARG_DgDx_dot_sg, ii).none())
653// outArgs.setSupports(OUT_ARG_DgDx_dot, i, DERIV_LINEAR_OP);
654// if (!me_outargs.supports(OUT_ARG_DgDx_sg, ii).none())
655// outArgs.setSupports(OUT_ARG_DgDx, i, DERIV_LINEAR_OP);
656 for (int j=0; j<num_p; j++)
657 outArgs.setSupports(OUT_ARG_DgDp, i, j,
658 me_outargs.supports(OUT_ARG_DgDp_sg, ii, j));
659 }
660
661 // We do not support derivatives w.r.t. the new SG parameters, so their
662 // support defaults to none.
663
664 return outArgs;
665}
666
667void
669 const OutArgs& outArgs) const
670{
671 // Get the input arguments
672 Teuchos::RCP<const Epetra_Vector> x;
673 if (inArgs.supports(IN_ARG_x)) {
674 x = inArgs.get_x();
675 if (x != Teuchos::null)
676 *my_x = *x;
677 }
678 Teuchos::RCP<const Epetra_Vector> x_dot;
679 if (inArgs.supports(IN_ARG_x_dot))
680 x_dot = inArgs.get_x_dot();
681
682 // Get the output arguments
683 EpetraExt::ModelEvaluator::Evaluation<Epetra_Vector> f_out;
684 if (outArgs.supports(OUT_ARG_f))
685 f_out = outArgs.get_f();
686 Teuchos::RCP<Epetra_Operator> W_out;
687 if (outArgs.supports(OUT_ARG_W))
688 W_out = outArgs.get_W();
689
690 // Create underlying inargs
691 InArgs me_inargs = me->createInArgs();
692 if (x != Teuchos::null) {
693 // copy to a poly orthog vector
694 adaptMngr->copyFromAdaptiveVector(*x,*x_sg_blocks);
695 me_inargs.set_x_sg(x_sg_blocks);
696 }
697 if (x_dot != Teuchos::null) {
698 // copy to a poly orthog vector
699 adaptMngr->copyFromAdaptiveVector(*x_dot,*x_dot_sg_blocks);
700 me_inargs.set_x_dot_sg(x_dot_sg_blocks);
701 }
702 if (me_inargs.supports(IN_ARG_alpha))
703 me_inargs.set_alpha(inArgs.get_alpha());
704 if (me_inargs.supports(IN_ARG_beta))
705 me_inargs.set_beta(inArgs.get_beta());
706 if (me_inargs.supports(IN_ARG_t))
707 me_inargs.set_t(inArgs.get_t());
708 if (me_inargs.supports(IN_ARG_sg_basis)) {
709 if (inArgs.get_sg_basis() != Teuchos::null)
710 me_inargs.set_sg_basis(inArgs.get_sg_basis());
711 else
712 me_inargs.set_sg_basis(sg_basis);
713 }
714 if (me_inargs.supports(IN_ARG_sg_quadrature)) {
715 if (inArgs.get_sg_quadrature() != Teuchos::null)
716 me_inargs.set_sg_quadrature(inArgs.get_sg_quadrature());
717 else
718 me_inargs.set_sg_quadrature(sg_quad);
719 }
720 if (me_inargs.supports(IN_ARG_sg_expansion)) {
721 if (inArgs.get_sg_expansion() != Teuchos::null)
722 me_inargs.set_sg_expansion(inArgs.get_sg_expansion());
723 else
724 me_inargs.set_sg_expansion(sg_exp);
725 }
726
727 // Pass parameters
728 for (int i=0; i<num_p; i++)
729 me_inargs.set_p(i, inArgs.get_p(i));
730 for (int i=0; i<num_p_sg; i++) {
731 Teuchos::RCP<const Epetra_Vector> p = inArgs.get_p(i+num_p);
732
733 // We always need to pass in the SG parameters, so just use
734 // the initial parameters if it is NULL
735 if (p == Teuchos::null)
736 p = sg_p_init[i]->getBlockVector();
737
738 // Convert block p to SG polynomial
739 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> p_sg =
740 create_p_sg(sg_p_index_map[i], View, p.get());
741 me_inargs.set_p_sg(sg_p_index_map[i], p_sg);
742 }
743
744 // Create underlying outargs
745 OutArgs me_outargs = me->createOutArgs();
746
747 // f
748 if (f_out != Teuchos::null) {
749 me_outargs.set_f_sg(f_sg_blocks);
750 if (eval_W_with_f)
751 me_outargs.set_W_sg(W_sg_blocks);
752 }
753
754 // W
755 if (W_out != Teuchos::null && !eval_W_with_f )
756 me_outargs.set_W_sg(W_sg_blocks);
757
758 // df/dp
759 for (int i=0; i<num_p; i++) {
760 if (!outArgs.supports(OUT_ARG_DfDp, i).none()) {
761 Derivative dfdp = outArgs.get_DfDp(i);
762 if (dfdp.getMultiVector() != Teuchos::null) {
763 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dfdp_sg;
764 if (dfdp.getMultiVectorOrientation() == DERIV_MV_BY_COL)
765 dfdp_sg =
767 sg_basis, overlapped_stoch_row_map,
768 me->get_f_map(), sg_comm,
769 me->get_p_map(i)->NumMyElements()));
770 else if (dfdp.getMultiVectorOrientation() == DERIV_TRANS_MV_BY_ROW)
771 dfdp_sg =
773 sg_basis, overlapped_stoch_row_map,
774 me->get_p_map(i), sg_comm,
775 me->get_f_map()->NumMyElements()));
776 me_outargs.set_DfDp_sg(i,
777 SGDerivative(dfdp_sg,
778 dfdp.getMultiVectorOrientation()));
779 }
780 TEUCHOS_TEST_FOR_EXCEPTION(dfdp.getLinearOp() != Teuchos::null, std::logic_error,
781 "Error! Stokhos::SGModelEvaluator_Adaptive::evalModel " <<
782 "cannot handle operator form of df/dp!");
783 }
784 }
785
786 // Responses (g, dg/dx, dg/dp, ...)
787 for (int i=0; i<num_g_sg; i++) {
788 int ii = sg_g_index_map[i];
789
790 // g
791 Teuchos::RCP<Epetra_Vector> g = outArgs.get_g(i);
792 if (g != Teuchos::null) {
793 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> g_sg =
794 create_g_sg(sg_g_index_map[i], View, g.get());
795 me_outargs.set_g_sg(ii, g_sg);
796 }
797
798 // dg/dxdot
799 if (outArgs.supports(OUT_ARG_DgDx_dot, i).supports(DERIV_LINEAR_OP)) {
800 Derivative dgdx_dot = outArgs.get_DgDx_dot(i);
801 if (dgdx_dot.getLinearOp() != Teuchos::null) {
802 Teuchos::RCP<Stokhos::SGOperator> op =
803 Teuchos::rcp_dynamic_cast<Stokhos::SGOperator>(
804 dgdx_dot.getLinearOp(), true);
805 Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > sg_blocks =
806 op->getSGPolynomial();
807 if (me_outargs.supports(OUT_ARG_DgDx, ii).supports(DERIV_LINEAR_OP))
808 me_outargs.set_DgDx_dot_sg(ii, sg_blocks);
809 else {
810 for (unsigned int k=0; k<num_sg_blocks; k++) {
811 Teuchos::RCP<Epetra_MultiVector> mv =
812 Teuchos::rcp_dynamic_cast<Stokhos::EpetraMultiVectorOperator>(
813 sg_blocks->getCoeffPtr(k), true)->getMultiVector();
814 dgdx_dot_sg_blocks[i]->setCoeffPtr(k, mv);
815 }
816 if (me_outargs.supports(OUT_ARG_DgDx_dot_sg, ii).supports(DERIV_MV_BY_COL))
817 me_outargs.set_DgDx_dot_sg(ii, SGDerivative(dgdx_dot_sg_blocks[i],
818 DERIV_MV_BY_COL));
819 else
820 me_outargs.set_DgDx_dot_sg(ii, SGDerivative(dgdx_dot_sg_blocks[i],
821 DERIV_TRANS_MV_BY_ROW));
822 }
823 }
824 TEUCHOS_TEST_FOR_EXCEPTION(dgdx_dot.getLinearOp() == Teuchos::null &&
825 dgdx_dot.isEmpty() == false,
826 std::logic_error,
827 "Error! Stokhos::SGModelEvaluator_Adaptive::evalModel: " <<
828 "Operator form of dg/dxdot is required!");
829 }
830
831 // dg/dx
832 if (outArgs.supports(OUT_ARG_DgDx, i).supports(DERIV_LINEAR_OP)) {
833 Derivative dgdx = outArgs.get_DgDx(i);
834 if (dgdx.getLinearOp() != Teuchos::null) {
835 Teuchos::RCP<Stokhos::SGOperator> op =
836 Teuchos::rcp_dynamic_cast<Stokhos::SGOperator>(
837 dgdx.getLinearOp(), true);
838 Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > sg_blocks =
839 op->getSGPolynomial();
840 if (me_outargs.supports(OUT_ARG_DgDx, ii).supports(DERIV_LINEAR_OP))
841 me_outargs.set_DgDx_sg(ii, sg_blocks);
842 else {
843 for (unsigned int k=0; k<num_sg_blocks; k++) {
844 Teuchos::RCP<Epetra_MultiVector> mv =
845 Teuchos::rcp_dynamic_cast<Stokhos::EpetraMultiVectorOperator>(
846 sg_blocks->getCoeffPtr(k), true)->getMultiVector();
847 dgdx_sg_blocks[i]->setCoeffPtr(k, mv);
848 }
849 if (me_outargs.supports(OUT_ARG_DgDx_sg, ii).supports(DERIV_MV_BY_COL))
850 me_outargs.set_DgDx_sg(ii, SGDerivative(dgdx_sg_blocks[i],
851 DERIV_MV_BY_COL));
852 else
853 me_outargs.set_DgDx_sg(ii, SGDerivative(dgdx_sg_blocks[i],
854 DERIV_TRANS_MV_BY_ROW));
855 }
856 }
857 TEUCHOS_TEST_FOR_EXCEPTION(dgdx.getLinearOp() == Teuchos::null &&
858 dgdx.isEmpty() == false,
859 std::logic_error,
860 "Error! Stokhos::SGModelEvaluator_Adaptive::evalModel: " <<
861 "Operator form of dg/dxdot is required!");
862 }
863
864 // dg/dp
865 // Rembember, no derivatives w.r.t. sg parameters
866 for (int j=0; j<num_p; j++) {
867 if (!outArgs.supports(OUT_ARG_DgDp, i, j).none()) {
868 Derivative dgdp = outArgs.get_DgDp(i,j);
869 if (dgdp.getMultiVector() != Teuchos::null) {
870 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> dgdp_sg;
871 if (dgdp.getMultiVectorOrientation() == DERIV_MV_BY_COL)
872 dgdp_sg =
874 sg_basis, overlapped_stoch_row_map,
875 me->get_g_map(ii), sg_g_map[i], sg_comm,
876 View, *(dgdp.getMultiVector())));
877 else if (dgdp.getMultiVectorOrientation() == DERIV_TRANS_MV_BY_ROW) {
878 Teuchos::RCP<const Epetra_BlockMap> product_map =
879 Teuchos::rcp(&(dgdp.getMultiVector()->Map()),false);
880 dgdp_sg =
882 sg_basis, overlapped_stoch_row_map,
883 me->get_p_map(j), product_map, sg_comm,
884 View, *(dgdp.getMultiVector())));
885 }
886 me_outargs.set_DgDp_sg(ii, j,
887 SGDerivative(dgdp_sg,
888 dgdp.getMultiVectorOrientation()));
889 }
890 TEUCHOS_TEST_FOR_EXCEPTION(dgdp.getLinearOp() != Teuchos::null,
891 std::logic_error,
892 "Error! Stokhos::SGModelEvaluator_Adaptive::evalModel " <<
893 "cannot handle operator form of dg/dp!");
894 }
895 }
896
897 }
898
899 // Compute the functions
900 me->evalModel(me_inargs, me_outargs);
901
902 // Copy block SG components for W
903 if ((W_out != Teuchos::null || (eval_W_with_f && f_out != Teuchos::null)) ) {
904
905 Teuchos::RCP<Epetra_Operator> W;
906 if (W_out != Teuchos::null)
907 W = W_out;
908 else
909 W = my_W;
910
911 Teuchos::RCP<Epetra_CrsMatrix> W_sg = Teuchos::rcp_dynamic_cast<Epetra_CrsMatrix>(W, true);
912 adaptMngr->setupOperator(*W_sg,*Cijk,*W_sg_blocks);
913 }
914
915 // f
916 if (f_out!=Teuchos::null){
917 if (!scaleOP)
918 for (int i=0; i<sg_basis->size(); i++)
919 (*f_sg_blocks)[i].Scale(sg_basis->norm_squared(i));
920
921 adaptMngr->copyToAdaptiveVector(*f_sg_blocks,*f_out);
922 }
923
924 // df/dp
925 for (int i=0; i<num_p; i++) {
926 if (!outArgs.supports(OUT_ARG_DfDp, i).none()) {
927 Derivative dfdp = outArgs.get_DfDp(i);
928 SGDerivative dfdp_sg = me_outargs.get_DfDp_sg(i);
929 if (dfdp.getMultiVector() != Teuchos::null) {
930
931 dfdp.getMultiVector()->Export(
932 *(dfdp_sg.getMultiVector()->getBlockMultiVector()),
933 *adapted_overlapped_f_exporter, Insert);
934 }
935 }
936 }
937}
938
939void
941 const Stokhos::EpetraVectorOrthogPoly& x_sg_in)
942{
943 *sg_x_init = x_sg_in;
944}
945
946Teuchos::RCP<const Stokhos::EpetraVectorOrthogPoly>
948{
949 return sg_x_init;
950}
951
952void
954 int i, const Stokhos::EpetraVectorOrthogPoly& p_sg_in)
955{
956 *sg_p_init[i] = p_sg_in;
957}
958
959Teuchos::RCP<const Stokhos::EpetraVectorOrthogPoly>
961{
962 return sg_p_init[l];
963}
964
965Teuchos::Array<int>
967{
968 return sg_p_index_map;
969}
970
971Teuchos::Array<int>
973{
974 return sg_g_index_map;
975}
976
977Teuchos::Array< Teuchos::RCP<const Epetra_Map> >
979{
980 Teuchos::Array< Teuchos::RCP<const Epetra_Map> > base_maps(num_g);
981 for (int i=0; i<num_g; i++)
982 base_maps[i] = me->get_g_map(i);
983 return base_maps;
984 }
985
986Teuchos::RCP<const Epetra_BlockMap>
988{
989 return overlapped_stoch_row_map;
990}
991
992Teuchos::RCP<const Epetra_BlockMap>
994{
995 return adapted_overlapped_x_map;
996}
997
998Teuchos::RCP<const Epetra_Import>
1000{
1001 return adapted_overlapped_x_importer;
1002}
1003
1004Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
1006 const Epetra_Vector* v) const
1007{
1008 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> sg_x;
1009 if (v == NULL)
1010 sg_x = Teuchos::rcp(new Stokhos::EpetraVectorOrthogPoly(
1011 sg_basis, stoch_row_map, x_map, sg_x_map, sg_comm));
1012 else
1013 sg_x = Teuchos::rcp(new Stokhos::EpetraVectorOrthogPoly(
1014 sg_basis, stoch_row_map, x_map, sg_x_map, sg_comm,
1015 CV, *v));
1016 return sg_x;
1017}
1018
1019Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
1021 const Epetra_Vector* v) const
1022{
1023 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> sg_x;
1024 if (v == NULL)
1025 sg_x = Teuchos::rcp(new Stokhos::EpetraVectorOrthogPoly(
1026 sg_basis, overlapped_stoch_row_map, x_map,
1027 sg_overlapped_x_map, sg_comm));
1028 else
1029 sg_x = Teuchos::rcp(new Stokhos::EpetraVectorOrthogPoly(
1030 sg_basis, overlapped_stoch_row_map, x_map,
1031 sg_overlapped_x_map, sg_comm, CV, *v));
1032 return sg_x;
1033}
1034
1035Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly>
1037 const Epetra_MultiVector* v) const
1038{
1039 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_x;
1040 if (v == NULL)
1041 sg_x = Teuchos::rcp(new Stokhos::EpetraMultiVectorOrthogPoly(
1042 sg_basis, stoch_row_map, x_map, sg_x_map, sg_comm,
1043 num_vecs));
1044 else
1045 sg_x = Teuchos::rcp(new Stokhos::EpetraMultiVectorOrthogPoly(
1046 sg_basis, stoch_row_map, x_map, sg_x_map, sg_comm,
1047 CV, *v));
1048 return sg_x;
1049}
1050
1051Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly>
1053 int num_vecs,
1055 const Epetra_MultiVector* v) const
1056{
1057 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_x;
1058 if (v == NULL)
1059 sg_x = Teuchos::rcp(new Stokhos::EpetraMultiVectorOrthogPoly(
1060 sg_basis, overlapped_stoch_row_map, x_map,
1061 sg_overlapped_x_map, sg_comm, num_vecs));
1062 else
1063 sg_x = Teuchos::rcp(new Stokhos::EpetraMultiVectorOrthogPoly(
1064 sg_basis, overlapped_stoch_row_map, x_map,
1065 sg_overlapped_x_map, sg_comm, CV, *v));
1066 return sg_x;
1067}
1068
1069Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
1071 const Epetra_Vector* v) const
1072{
1073 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> sg_p;
1074 Teuchos::Array<int>::const_iterator it = std::find(sg_p_index_map.begin(),
1075 sg_p_index_map.end(),
1076 l);
1077 TEUCHOS_TEST_FOR_EXCEPTION(it == sg_p_index_map.end(), std::logic_error,
1078 "Error! Invalid p map index " << l);
1079 int ll = it - sg_p_index_map.begin();
1080 if (v == NULL)
1081 sg_p = Teuchos::rcp(new Stokhos::EpetraVectorOrthogPoly(
1082 sg_basis, overlapped_stoch_p_map, me->get_p_map(l),
1083 sg_p_map[ll], sg_comm));
1084 else
1085 sg_p = Teuchos::rcp(new Stokhos::EpetraVectorOrthogPoly(
1086 sg_basis, overlapped_stoch_p_map, me->get_p_map(l),
1087 sg_p_map[ll], sg_comm, CV, *v));
1088 return sg_p;
1089}
1090
1091Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
1093 const Epetra_Vector* v) const
1094{
1095 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> sg_f;
1096 if (v == NULL)
1097 sg_f = Teuchos::rcp(new Stokhos::EpetraVectorOrthogPoly(
1098 sg_basis, stoch_row_map, f_map, sg_f_map, sg_comm));
1099 else
1100 sg_f = Teuchos::rcp(new Stokhos::EpetraVectorOrthogPoly(
1101 sg_basis, stoch_row_map, f_map, sg_f_map, sg_comm,
1102 CV, *v));
1103 return sg_f;
1104}
1105
1106Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
1108 const Epetra_Vector* v) const
1109{
1110 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> sg_f;
1111 if (v == NULL)
1112 sg_f = Teuchos::rcp(new Stokhos::EpetraVectorOrthogPoly(
1113 sg_basis, overlapped_stoch_row_map, f_map,
1114 sg_overlapped_f_map, sg_comm));
1115 else
1116 sg_f = Teuchos::rcp(new Stokhos::EpetraVectorOrthogPoly(
1117 sg_basis, overlapped_stoch_row_map, f_map,
1118 sg_overlapped_f_map, sg_comm, CV, *v));
1119 return sg_f;
1120}
1121
1122Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly>
1124 int num_vecs,
1126 const Epetra_MultiVector* v) const
1127{
1128 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_f;
1129 if (v == NULL)
1130 sg_f = Teuchos::rcp(new Stokhos::EpetraMultiVectorOrthogPoly(
1131 sg_basis, stoch_row_map, f_map, sg_f_map, sg_comm,
1132 num_vecs));
1133 else
1134 sg_f = Teuchos::rcp(new Stokhos::EpetraMultiVectorOrthogPoly(
1135 sg_basis, stoch_row_map, f_map, sg_f_map, sg_comm,
1136 CV, *v));
1137 return sg_f;
1138}
1139
1140Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly>
1142 int num_vecs,
1144 const Epetra_MultiVector* v) const
1145{
1146 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_f;
1147 if (v == NULL)
1148 sg_f = Teuchos::rcp(new Stokhos::EpetraMultiVectorOrthogPoly(
1149 sg_basis, overlapped_stoch_row_map, f_map,
1150 sg_overlapped_f_map, sg_comm, num_vecs));
1151 else
1152 sg_f = Teuchos::rcp(new Stokhos::EpetraMultiVectorOrthogPoly(
1153 sg_basis, overlapped_stoch_row_map, f_map,
1154 sg_overlapped_f_map, sg_comm, CV, *v));
1155 return sg_f;
1156}
1157
1158Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly>
1160 const Epetra_Vector* v) const
1161{
1162 Teuchos::RCP<Stokhos::EpetraVectorOrthogPoly> sg_g;
1163 Teuchos::Array<int>::const_iterator it = std::find(sg_g_index_map.begin(),
1164 sg_g_index_map.end(),
1165 l);
1166 TEUCHOS_TEST_FOR_EXCEPTION(it == sg_g_index_map.end(), std::logic_error,
1167 "Error! Invalid g map index " << l);
1168 int ll = it - sg_g_index_map.begin();
1169 if (v == NULL)
1170 sg_g = Teuchos::rcp(new Stokhos::EpetraVectorOrthogPoly(
1171 sg_basis, overlapped_stoch_row_map,
1172 me->get_g_map(l),
1173 sg_g_map[ll], sg_comm));
1174 else
1175 sg_g = Teuchos::rcp(new Stokhos::EpetraVectorOrthogPoly(
1176 sg_basis, overlapped_stoch_row_map,
1177 me->get_g_map(l),
1178 sg_g_map[ll], sg_comm, CV, *v));
1179 return sg_g;
1180}
1181
1182Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly>
1185 const Epetra_MultiVector* v) const
1186{
1187 Teuchos::RCP<Stokhos::EpetraMultiVectorOrthogPoly> sg_g;
1188 Teuchos::Array<int>::const_iterator it = std::find(sg_g_index_map.begin(),
1189 sg_g_index_map.end(),
1190 l);
1191 TEUCHOS_TEST_FOR_EXCEPTION(it == sg_g_index_map.end(), std::logic_error,
1192 "Error! Invalid g map index " << l);
1193 int ll = it - sg_g_index_map.begin();
1194 if (v == NULL)
1195 sg_g = Teuchos::rcp(new Stokhos::EpetraMultiVectorOrthogPoly(
1196 sg_basis, overlapped_stoch_row_map,
1197 me->get_g_map(l),
1198 sg_g_map[ll], sg_comm, num_vecs));
1199 else
1200 sg_g = Teuchos::rcp(new Stokhos::EpetraMultiVectorOrthogPoly(
1201 sg_basis, overlapped_stoch_row_map,
1202 me->get_g_map(l),
1203 sg_g_map[ll], sg_comm, CV, *v));
1204 return sg_g;
1205}
Insert
Epetra_DataAccess
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,...
Abstract base class for multivariate orthogonal polynomials.
Abstract base class for orthogonal polynomial-based expansions.
Abstract base class for multivariate orthogonal polynomials generated from tensor products of univari...
Abstract base class for quadrature methods.
Teuchos::RCP< const Epetra_Vector > get_x_init() const
Return initial solution.
Teuchos::RCP< const Epetra_Import > get_x_sg_importer() const
Return x sg importer.
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.
unsigned int num_sg_blocks
Number of stochastic blocks.
void evalModel(const InArgs &inArgs, const OutArgs &outArgs) const
Evaluate model on InArgs.
Teuchos::RCP< const Epetra_Vector > get_p_init(int l) const
Return initial parameters.
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::Array< Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > > sg_p_init
SG initial p.
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::RCP< const Epetra_Map > sg_f_map
Block SG residual map.
Teuchos::Array< Teuchos::RCP< Teuchos::Array< std::string > > > sg_p_names
SG coefficient parameter names.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > f_sg_blocks
f stochastic Galerkin components
Teuchos::RCP< const Stokhos::EpetraVectorOrthogPoly > get_x_sg_init() const
Return initial SG x.
Teuchos::Array< Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > > dgdx_dot_sg_blocks
dg/dxdot stochastic Galerkin components
Teuchos::RCP< const Epetra_BlockMap > stoch_row_map
Map for stochastic blocks.
Teuchos::Array< Teuchos::RCP< Stokhos::EpetraMultiVectorOrthogPoly > > dgdx_sg_blocks
dg/dx stochastic Galerkin components
Teuchos::Array< int > sg_p_index_map
Index map between block-p and p_sg maps.
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< EpetraExt::ModelEvaluator > me
Underlying model evaluator.
std::vector< Teuchos::RCP< const Stokhos::ProductBasis< int, double > > > sg_row_dof_basis
Teuchos::RCP< const Epetra_Map > get_x_map() const
Return solution vector map.
Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > W_sg_blocks
W stochastic Galerkin components.
Teuchos::RCP< const Epetra_BlockMap > get_overlap_stochastic_map() const
Return overlap stochastic map.
Teuchos::RCP< const Epetra_BlockMap > overlapped_stoch_row_map
Overlapped map for stochastic blocks (local map)
Teuchos::RCP< const Epetra_Map > f_map
Underlying residual map.
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.
int num_g_sg
Number of stochastic response vectors.
Teuchos::Array< Teuchos::RCP< const Epetra_Map > > sg_g_map
Block SG response map.
void set_p_sg_init(int i, const Stokhos::EpetraVectorOrthogPoly &p_sg_in)
Set initial parameter polynomial.
Teuchos::RCP< const Epetra_Map > get_g_map(int l) const
Return response 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< 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::Array< Teuchos::RCP< const Epetra_Map > > sg_p_map
Block SG parameter map.
Teuchos::RCP< const Teuchos::Array< std::string > > get_p_names(int l) const
Return array of parameter names.
Teuchos::RCP< Epetra_Import > adapted_overlapped_x_importer
Importer from SG to SG-overlapped maps.
unsigned int num_p_blocks
Number of p stochastic blocks (may be smaller than num_sg_blocks)
Teuchos::Array< Teuchos::RCP< const Epetra_Map > > get_g_sg_base_maps() const
Get base maps of SG responses.
Teuchos::RCP< const Epetra_Map > adapted_overlapped_f_map
Adapted block SG overlapped residual map.
unsigned int num_W_blocks
Number of W stochastic blocks (may be smaller than num_sg_blocks)
Teuchos::RCP< const Epetra_Map > sg_overlapped_f_map
Block SG overlapped residual map.
Teuchos::RCP< const Stokhos::ParallelData > sg_parallel_data
Parallel SG data.
Teuchos::RCP< Epetra_Vector > my_x
x pointer for evaluating preconditioner
Teuchos::Array< int > get_g_sg_map_indices() const
Get indices of SG responses.
Teuchos::RCP< const Epetra_Map > sg_overlapped_x_map
Block SG overlapped unknown map.
Teuchos::RCP< const Epetra_Map > adapted_f_map
Adapted block SG residual map.
Teuchos::RCP< Stokhos::AdaptivityManager > adaptMngr
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::Array< int > sg_g_index_map
Index map between block-g and g_sg maps.
Teuchos::RCP< const Stokhos::OrthogPolyBasis< int, double > > sg_basis
Stochastic Galerkin basis.
Teuchos::RCP< const Epetra_Map > x_map
Underlying 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 Stokhos::EpetraVectorOrthogPoly > get_p_sg_init(int l) const
Return initial SG parameters.
Teuchos::RCP< Stokhos::EpetraVectorOrthogPoly > x_dot_sg_blocks
x_dot stochastic Galerkin components
Teuchos::RCP< const Epetra_Map > get_f_map() const
Return residual vector map.
int num_p_sg
Number of stochastic parameter vectors.
Teuchos::RCP< Teuchos::ParameterList > params
Algorithmic parameters.
Teuchos::RCP< Epetra_Operator > create_W() const
Create W = alpha*M + beta*J matrix.
Teuchos::RCP< const Epetra_BlockMap > get_x_sg_overlap_map() const
Return x sg overlap map.
int num_g
Number of response vectors of underlying model evaluator.
Teuchos::RCP< const Epetra_Map > get_p_map(int l) const
Return parameter vector map.
Teuchos::Array< int > get_p_sg_map_indices() const
Get indices of SG parameters.
Teuchos::RCP< const Epetra_Map > adapted_x_map
Adapted lock SG unknown map.
Teuchos::RCP< const Epetra_Map > sg_x_map
Block SG unknown map.
void set_x_sg_init(const Stokhos::EpetraVectorOrthogPoly &x_sg_in)
Set initial solution polynomial.
Teuchos::RCP< const EpetraExt::MultiComm > sg_comm
Parallel SG communicator.
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.
Teuchos::RCP< const Epetra_Map > adapted_overlapped_x_map
Adapated block SG overlapped unknown map.
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< const Epetra_BlockMap > overlapped_stoch_p_map
Overlapped map for p stochastic blocks (local map)
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.
int num_p
Number of parameter vectors of underlying model evaluator.
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::EpetraVectorOrthogPoly > sg_x_init
SG initial x.
bool supports_x
Whether we support x (and thus f and W)
Teuchos::RCP< Epetra_Export > adapted_overlapped_f_exporter
Exporter from SG-overlapped to SG maps.
Teuchos::RCP< const Stokhos::EpetraSparse3Tensor > serialCijk
Serial Epetra Cijk for dgdx*.
Teuchos::RCP< const Stokhos::EpetraSparse3Tensor > epetraCijk
Epetra Cijk.
ScalarType g(const Teuchos::Array< ScalarType > &x, const ScalarType &y)