Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
Stokhos_SmolyakPseudoSpectralOperator.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_SMOLYAK_PSEUDO_SPECTRAL_OPERATOR_HPP
43 #define STOKHOS_SMOLYAK_PSEUDO_SPECTRAL_OPERATOR_HPP
44 
46 #include "Stokhos_SmolyakBasis.hpp"
48 #include "Teuchos_SerialDenseMatrix.hpp"
49 #include "Teuchos_BLAS.hpp"
50 
51 namespace Stokhos {
52 
57  template <typename ordinal_t,
58  typename value_t,
59  typename point_compare_type =
62  public PseudoSpectralOperator<ordinal_t,value_t,point_compare_type> {
63  public:
64 
65  typedef ordinal_t ordinal_type;
66  typedef value_t value_type;
68  typedef typename base_type::point_type point_type;
71  typedef typename base_type::iterator iterator;
75 
78  typedef Teuchos::Array< Teuchos::RCP<operator_type> > operator_set_type;
79 
81  template <typename coeff_compare_type>
84  bool use_smolyak_apply = true,
85  bool use_pst = true,
86  const point_compare_type& point_compare = point_compare_type());
87 
90 
92  ordinal_type point_size() const { return points.size(); }
93 
95  ordinal_type coeff_size() const { return coeff_sz; }
96 
98  iterator begin() { return point_map.begin(); }
99 
101  iterator end() { return point_map.end(); }
102 
104  const_iterator begin() const { return point_map.begin(); }
105 
107  const_iterator end() const { return point_map.end(); }
108 
110  set_iterator set_begin() { return points.begin(); }
111 
113  set_iterator set_end() { return points.end(); }
114 
116  const_set_iterator set_begin() const { return points.begin(); }
117 
119  const_set_iterator set_end() const { return points.end(); }
120 
123  const_set_iterator it = points.find(point);
124  TEUCHOS_TEST_FOR_EXCEPTION(
125  it == points.end(), std::logic_error, "Invalid term " << point);
126  return it->second.second;
127  }
128 
130  const point_type& point(ordinal_type n) const { return point_map[n]; }
131 
133 
140  void transformQP2PCE(
141  const value_type& alpha,
142  const Teuchos::SerialDenseMatrix<ordinal_type,value_type>& input,
143  Teuchos::SerialDenseMatrix<ordinal_type,value_type>& result,
144  const value_type& beta,
145  bool trans = false) const;
146 
148 
155  virtual void transformPCE2QP(
156  const value_type& alpha,
157  const Teuchos::SerialDenseMatrix<ordinal_type,value_type>& input,
158  Teuchos::SerialDenseMatrix<ordinal_type,value_type>& result,
159  const value_type& beta,
160  bool trans = false) const;
161 
162  protected:
163 
165  void apply_direct(
166  const Teuchos::SerialDenseMatrix<ordinal_type,value_type>& A,
167  const value_type& alpha,
168  const Teuchos::SerialDenseMatrix<ordinal_type,value_type>& input,
169  Teuchos::SerialDenseMatrix<ordinal_type,value_type>& result,
170  const value_type& beta,
171  bool trans) const;
172 
178  const value_type& alpha,
179  const Teuchos::SerialDenseMatrix<ordinal_type,value_type>& input,
180  Teuchos::SerialDenseMatrix<ordinal_type,value_type>& result,
181  const value_type& beta,
182  bool trans) const;
183 
189  const value_type& alpha,
190  const Teuchos::SerialDenseMatrix<ordinal_type,value_type>& input,
191  Teuchos::SerialDenseMatrix<ordinal_type,value_type>& result,
192  const value_type& beta,
193  bool trans) const;
194 
195  void gather(
196  const Teuchos::Array<ordinal_type>& map,
197  const Teuchos::SerialDenseMatrix<ordinal_type,value_type>& input,
198  bool trans,
199  Teuchos::SerialDenseMatrix<ordinal_type,value_type>& result) const;
200 
201  void scatter(
202  const Teuchos::Array<ordinal_type>& map,
203  const Teuchos::SerialDenseMatrix<ordinal_type,value_type>& input,
204  bool trans,
205  Teuchos::SerialDenseMatrix<ordinal_type,value_type>& result) const;
206 
207  protected:
208 
211 
214 
217 
220 
222  Teuchos::Array<value_type> smolyak_coeffs;
223 
226 
228  Teuchos::Array< Teuchos::Array<ordinal_type> > gather_maps;
229 
231  Teuchos::Array< Teuchos::Array<ordinal_type> > scatter_maps;
232 
234  Teuchos::SerialDenseMatrix<ordinal_type,value_type> qp2pce;
235 
237  Teuchos::SerialDenseMatrix<ordinal_type,value_type> pce2qp;
238 
240  Teuchos::BLAS<ordinal_type,value_type> blas;
241 
242  };
243 
244 }
245 
247 
248 #endif
Teuchos::SerialDenseMatrix< ordinal_type, value_type > pce2qp
Matrix mapping coefficients to points.
std::map< point_type, std::pair< value_type, ordinal_type >, point_compare_type > point_set_type
const_set_iterator set_end() const
Iterator to end of point set.
point_map_type point_map
Map index to sparse grid term.
Teuchos::SerialDenseMatrix< ordinal_type, value_type > qp2pce
Matrix mapping points to coefficients.
set_iterator set_end()
Iterator to end of point set.
Container storing a term in a generalized tensor product.
const_iterator end() const
Iterator to end of point set.
An operator interface for building pseudo-spectral approximations.
iterator begin()
Iterator to begining of point set.
virtual void transformPCE2QP(const value_type &alpha, const Teuchos::SerialDenseMatrix< ordinal_type, value_type > &input, Teuchos::SerialDenseMatrix< ordinal_type, value_type > &result, const value_type &beta, bool trans=false) const
Transform PCE coefficients to quadrature values.
ordinal_type coeff_size() const
Number of coefficients.
Teuchos::Array< Teuchos::RCP< operator_type > > operator_set_type
void apply_direct(const Teuchos::SerialDenseMatrix< ordinal_type, value_type > &A, const value_type &alpha, const Teuchos::SerialDenseMatrix< ordinal_type, value_type > &input, Teuchos::SerialDenseMatrix< ordinal_type, value_type > &result, const value_type &beta, bool trans) const
Apply transformation operator using direct method.
Teuchos::Array< Teuchos::Array< ordinal_type > > gather_maps
Gather maps for each operator for Smolyak apply.
ordinal_type index(const point_type &point) const
Get point index for given point.
Teuchos::Array< Teuchos::Array< ordinal_type > > scatter_maps
Scatter maps for each operator for Smolyak apply.
const_iterator begin() const
Iterator to begining of point set.
Top-level namespace for Stokhos classes and functions.
An operator for building pseudo-spectral coefficients using a sparse Smolyak construction.
Teuchos::BLAS< ordinal_type, value_type > blas
BLAS wrappers.
void transformQP2PCE(const value_type &alpha, const Teuchos::SerialDenseMatrix< ordinal_type, value_type > &input, Teuchos::SerialDenseMatrix< ordinal_type, value_type > &result, const value_type &beta, bool trans=false) const
Transform values at quadrature points to PCE coefficients.
Teuchos::Array< value_type > smolyak_coeffs
Smolyak coefficients.
Multivariate orthogonal polynomial basis generated from a Smolyak sparse grid.
const point_type & point(ordinal_type n) const
Get point for given index.
void transformPCE2QP_smolyak(const value_type &alpha, const Teuchos::SerialDenseMatrix< ordinal_type, value_type > &input, Teuchos::SerialDenseMatrix< ordinal_type, value_type > &result, const value_type &beta, bool trans) const
Transform PCE coefficients to values at quadrature points using Smolyak formula.
void scatter(const Teuchos::Array< ordinal_type > &map, const Teuchos::SerialDenseMatrix< ordinal_type, value_type > &input, bool trans, Teuchos::SerialDenseMatrix< ordinal_type, value_type > &result) const
const_set_iterator set_begin() const
Iterator to begining of point set.
SmolyakPseudoSpectralOperator(const SmolyakBasis< ordinal_type, value_type, coeff_compare_type > &smolyak_basis, bool use_smolyak_apply=true, bool use_pst=true, const point_compare_type &point_compare=point_compare_type())
Constructor.
An operator for building pseudo-spectral coefficients using tensor-product quadrature.
TensorProductPseudoSpectralOperator< ordinal_type, value_type > operator_type
void gather(const Teuchos::Array< ordinal_type > &map, const Teuchos::SerialDenseMatrix< ordinal_type, value_type > &input, bool trans, Teuchos::SerialDenseMatrix< ordinal_type, value_type > &result) const
set_iterator set_begin()
Iterator to begining of point set.
point_map_type::const_iterator const_iterator
LexographicLess< TensorProductElement< ordinal_type, value_type >, FloatingPointLess< value_type > > type
PseudoSpectralOperator< ordinal_type, value_type, point_compare_type > base_type
point_set_type::const_iterator const_set_iterator
void transformQP2PCE_smolyak(const value_type &alpha, const Teuchos::SerialDenseMatrix< ordinal_type, value_type > &input, Teuchos::SerialDenseMatrix< ordinal_type, value_type > &result, const value_type &beta, bool trans) const
Transform values at quadrature points to PCE coefficients using Smolyak formula.
operator_set_type operators
Tensor pseudospectral operators.