Stokhos Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
Stokhos_StieltjesGramSchmidtBuilderImp.hpp
Go to the documentation of this file.
1// $Id: Stokhos_SGModelEvaluator.cpp,v 1.10 2009/10/06 16:51:22 agsalin Exp $
2// $Source: /space/CVS/Trilinos/packages/stokhos/src/Stokhos_SGModelEvaluator.cpp,v $
3// @HEADER
4// ***********************************************************************
5//
6// Stokhos Package
7// Copyright (2009) Sandia Corporation
8//
9// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
10// license for use of this work by or on behalf of the U.S. Government.
11//
12// Redistribution and use in source and binary forms, with or without
13// modification, are permitted provided that the following conditions are
14// met:
15//
16// 1. Redistributions of source code must retain the above copyright
17// notice, this list of conditions and the following disclaimer.
18//
19// 2. Redistributions in binary form must reproduce the above copyright
20// notice, this list of conditions and the following disclaimer in the
21// documentation and/or other materials provided with the distribution.
22//
23// 3. Neither the name of the Corporation nor the names of the
24// contributors may be used to endorse or promote products derived from
25// this software without specific prior written permission.
26//
27// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
28// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
29// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
30// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
31// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
32// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
33// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
34// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
35// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
36// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
37// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
38//
39// Questions? Contact Eric T. Phipps (etphipp@sandia.gov).
40//
41// ***********************************************************************
42// @HEADER
43
47
48template <typename ordinal_type, typename value_type>
51 const Teuchos::RCP<const Stokhos::Quadrature<ordinal_type, value_type> >& quad_,
52 const Teuchos::Array< Stokhos::OrthogPolyApprox<ordinal_type, value_type> >& pces,
53 ordinal_type new_order, bool use_pce_qp, bool normalize) :
54 quad(quad_)
55{
56 // Create array to store new coordinate bases
57 ordinal_type new_dim = pces.size();
58 Teuchos::Array< Teuchos::RCP<const OneDOrthogPolyBasis<ordinal_type, value_type> > > new_coordinate_bases(new_dim);
59
60 // Create Stieltjes basis for each pce
61 for (ordinal_type k=0; k<new_dim; k++) {
62 new_coordinate_bases[k] = Teuchos::rcp(
63 new StieltjesPCEBasis<ordinal_type,value_type>(
64 new_order, Teuchos::rcp(&(pces[k]),false), quad, use_pce_qp,
65 normalize)
66 );
67 }
68
69 // Create tensor product basis from coordinate bases
70 tensor_basis = Teuchos::rcp(
71 new CompletePolynomialBasis<ordinal_type,value_type>(new_coordinate_bases)
72 );
73
74 // Use Gram-Schmidt to orthogonalize tensor product bases
75 const Teuchos::Array<value_type>& weights = quad->getQuadWeights();
76 const Teuchos::Array< Teuchos::Array<value_type> >& points =
77 quad->getQuadPoints();
78 ordinal_type nqp = points.size();
79 Teuchos::RCP< Teuchos::Array< Teuchos::Array<value_type> > > new_points =
80 Teuchos::rcp(new Teuchos::Array< Teuchos::Array<value_type> >(nqp));
81 Teuchos::RCP< Teuchos::Array<value_type> > new_weights =
82 Teuchos::rcp(new Teuchos::Array<value_type>(weights));
83 for (ordinal_type i=0; i<nqp; i++)
84 (*new_points)[i].resize(new_dim);
85 for (ordinal_type k=0; k<new_dim; k++) {
86 Teuchos::Array<value_type> st_points;
87 Teuchos::Array<value_type> st_weights;
88 Teuchos::Array< Teuchos::Array<value_type> > st_values;
89 new_coordinate_bases[k]->getQuadPoints(new_order+1, st_points, st_weights,
90 st_values);
91
92 for (ordinal_type i=0; i<nqp; i++)
93 (*new_points)[i][k] = st_points[i];
94 }
95 gs_basis = Teuchos::rcp(
96 new GramSchmidtBasis<ordinal_type,value_type>(tensor_basis,
97 *new_points,
98 *new_weights,
99 1e-15)
100 );
101
102 // Create new quadrature object
103 Teuchos::RCP<const OrthogPolyBasis<ordinal_type,value_type> > new_basis =
104 gs_basis;
105 gs_quad = Teuchos::rcp(
106 new UserDefinedQuadrature<ordinal_type,value_type>(new_basis,
107 new_points,
108 new_weights)
109 );
110
111}
112
113template <typename ordinal_type, typename value_type>
114Teuchos::RCP<const Stokhos::OrthogPolyBasis<ordinal_type, value_type> >
116getReducedBasis() const
117{
118 return gs_basis;
119}
120
121template <typename ordinal_type, typename value_type>
122Teuchos::RCP<Stokhos::Quadrature<ordinal_type, value_type> >
125{
126 return gs_quad;
127}
128
129template <typename ordinal_type, typename value_type>
130void
133 const Teuchos::Array< Stokhos::OrthogPolyApprox<ordinal_type, value_type> >& pces,
134 Teuchos::Array< Stokhos::OrthogPolyApprox<ordinal_type, value_type> >& new_pces)
135{
136 // Map pce coefficients to tensor basis to Gram-Schmidt basis
137 ordinal_type dim = pces.size();
138 if (new_pces.size() != pces.size())
139 new_pces.resize(dim);
140 for (ordinal_type k=0; k<dim; k++) {
141 OrthogPolyApprox<ordinal_type,value_type> p_tensor(tensor_basis);
142 p_tensor.term(k, 0) = pces[k].mean();
143 p_tensor.term(k, 1) = 1.0;
144 new_pces[k].reset(gs_basis);
145 gs_basis->transformCoeffs(p_tensor.coeff(), new_pces[k].coeff());
146 }
147}
Class to store coefficients of a projection onto an orthogonal polynomial basis.
pointer coeff()
Return coefficient array.
reference term(ordinal_type dimension, ordinal_type order)
Get coefficient term for given dimension and order.
Abstract base class for quadrature methods.
void computeReducedPCEs(const Teuchos::Array< OrthogPolyApprox< ordinal_type, value_type > > &pces, Teuchos::Array< OrthogPolyApprox< ordinal_type, value_type > > &new_pces)
Get reduced PCEs.
StieltjesGramSchmidtBuilder(const Teuchos::RCP< const Quadrature< ordinal_type, value_type > > &quad, const Teuchos::Array< OrthogPolyApprox< ordinal_type, value_type > > &pces, ordinal_type new_order, bool use_pce_qp, bool normalize)
Constructor.
Teuchos::RCP< const OrthogPolyBasis< ordinal_type, value_type > > getReducedBasis() const
Get reduced basis.
Teuchos::RCP< Quadrature< ordinal_type, value_type > > getReducedQuadrature() const
Get reduced quadrature.