Stokhos Package Browser (Single Doxygen Collection)
Version of the Day
|
An operator for building pseudo-spectral coefficients using tensor-product quadrature. More...
#include <Stokhos_TensorProductPseudoSpectralOperator.hpp>
Public Member Functions | |
TensorProductPseudoSpectralOperator (const ProductBasis< ordinal_type, value_type > &product_basis, bool use_pst_=false, multiindex_type multiindex=multiindex_type(), const point_compare_type &point_compare=point_compare_type()) | |
Constructor. More... | |
virtual | ~TensorProductPseudoSpectralOperator () |
Destructor. More... | |
ordinal_type | point_size () const |
Number of points. More... | |
ordinal_type | coeff_size () const |
Number of coefficients. More... | |
iterator | begin () |
Iterator to begining of point set. More... | |
iterator | end () |
Iterator to end of point set. More... | |
const_iterator | begin () const |
Iterator to begining of point set. More... | |
const_iterator | end () const |
Iterator to end of point set. More... | |
set_iterator | set_begin () |
Iterator to begining of point set. More... | |
set_iterator | set_end () |
Iterator to end of point set. More... | |
const_set_iterator | set_begin () const |
Iterator to begining of point set. More... | |
const_set_iterator | set_end () const |
Iterator to end of point set. More... | |
ordinal_type | index (const point_type &point) const |
Get point index for given point. More... | |
const point_type & | point (ordinal_type n) const |
Get point for given index. More... | |
virtual 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. More... | |
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. More... | |
Public Member Functions inherited from Stokhos::PseudoSpectralOperator< ordinal_t, value_t, point_compare_type > | |
PseudoSpectralOperator () | |
Constructor. More... | |
virtual | ~PseudoSpectralOperator () |
Destructor. More... | |
Protected Member Functions | |
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. More... | |
void | apply_pst (const Teuchos::Array< Teuchos::SerialDenseMatrix< ordinal_type, value_type > > &Ak, 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, bool reorder_input, bool reorder_result) const |
Apply tranformation operator using PST method. More... | |
Protected Attributes | |
bool | use_pst |
Use partial-summation-transformation. More... | |
ordinal_type | coeff_sz |
Number of coefficients. More... | |
ordinal_type | dim |
Dimension. More... | |
point_set_type | points |
Quadrature points. More... | |
point_map_type | point_map |
Map index to point term. More... | |
bool | reorder |
Do we need to reorder coefficients for PST. More... | |
Teuchos::Array< ordinal_type > | perm |
Permutation array when reordering for PST. More... | |
Teuchos::SerialDenseMatrix < ordinal_type, value_type > | qp2pce |
Matrix mapping points to coefficients. More... | |
Teuchos::SerialDenseMatrix < ordinal_type, value_type > | pce2qp |
Matrix mapping coefficients to points. More... | |
Teuchos::Array < Teuchos::SerialDenseMatrix < ordinal_type, value_type > > | qp2pce_k |
Matrix mapping points to coefficients for each dimension for PST. More... | |
Teuchos::Array < Teuchos::SerialDenseMatrix < ordinal_type, value_type > > | pce2qp_k |
Matrix mapping coefficients to points for each dimension for PST. More... | |
Teuchos::BLAS< ordinal_type, value_type > | blas |
BLAS wrappers. More... | |
An operator for building pseudo-spectral coefficients using tensor-product quadrature.
Definition at line 29 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
typedef ordinal_t Stokhos::TensorProductPseudoSpectralOperator< ordinal_t, value_t, point_compare_type >::ordinal_type |
Definition at line 33 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
typedef value_t Stokhos::TensorProductPseudoSpectralOperator< ordinal_t, value_t, point_compare_type >::value_type |
Definition at line 34 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
typedef PseudoSpectralOperator<ordinal_type,value_type,point_compare_type> Stokhos::TensorProductPseudoSpectralOperator< ordinal_t, value_t, point_compare_type >::base_type |
Definition at line 35 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
typedef base_type::point_type Stokhos::TensorProductPseudoSpectralOperator< ordinal_t, value_t, point_compare_type >::point_type |
Definition at line 36 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
typedef base_type::point_set_type Stokhos::TensorProductPseudoSpectralOperator< ordinal_t, value_t, point_compare_type >::point_set_type |
Definition at line 37 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
typedef base_type::point_map_type Stokhos::TensorProductPseudoSpectralOperator< ordinal_t, value_t, point_compare_type >::point_map_type |
Definition at line 38 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
typedef base_type::iterator Stokhos::TensorProductPseudoSpectralOperator< ordinal_t, value_t, point_compare_type >::iterator |
Definition at line 39 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
typedef base_type::const_iterator Stokhos::TensorProductPseudoSpectralOperator< ordinal_t, value_t, point_compare_type >::const_iterator |
Definition at line 40 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
typedef base_type::set_iterator Stokhos::TensorProductPseudoSpectralOperator< ordinal_t, value_t, point_compare_type >::set_iterator |
Definition at line 41 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
typedef base_type::const_set_iterator Stokhos::TensorProductPseudoSpectralOperator< ordinal_t, value_t, point_compare_type >::const_set_iterator |
Definition at line 42 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
typedef MultiIndex<ordinal_type> Stokhos::TensorProductPseudoSpectralOperator< ordinal_t, value_t, point_compare_type >::multiindex_type |
Definition at line 44 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inline |
Constructor.
Definition at line 47 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlinevirtual |
Destructor.
Definition at line 185 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlinevirtual |
Number of points.
Implements Stokhos::PseudoSpectralOperator< ordinal_t, value_t, point_compare_type >.
Definition at line 188 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlinevirtual |
Number of coefficients.
Implements Stokhos::PseudoSpectralOperator< ordinal_t, value_t, point_compare_type >.
Definition at line 191 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlinevirtual |
Iterator to begining of point set.
Implements Stokhos::PseudoSpectralOperator< ordinal_t, value_t, point_compare_type >.
Definition at line 194 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlinevirtual |
Iterator to end of point set.
Implements Stokhos::PseudoSpectralOperator< ordinal_t, value_t, point_compare_type >.
Definition at line 197 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlinevirtual |
Iterator to begining of point set.
Implements Stokhos::PseudoSpectralOperator< ordinal_t, value_t, point_compare_type >.
Definition at line 200 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlinevirtual |
Iterator to end of point set.
Implements Stokhos::PseudoSpectralOperator< ordinal_t, value_t, point_compare_type >.
Definition at line 203 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlinevirtual |
Iterator to begining of point set.
Implements Stokhos::PseudoSpectralOperator< ordinal_t, value_t, point_compare_type >.
Definition at line 206 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlinevirtual |
Iterator to end of point set.
Implements Stokhos::PseudoSpectralOperator< ordinal_t, value_t, point_compare_type >.
Definition at line 209 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlinevirtual |
Iterator to begining of point set.
Implements Stokhos::PseudoSpectralOperator< ordinal_t, value_t, point_compare_type >.
Definition at line 212 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlinevirtual |
Iterator to end of point set.
Implements Stokhos::PseudoSpectralOperator< ordinal_t, value_t, point_compare_type >.
Definition at line 215 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlinevirtual |
Get point index for given point.
Implements Stokhos::PseudoSpectralOperator< ordinal_t, value_t, point_compare_type >.
Definition at line 218 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlinevirtual |
Get point for given index.
Implements Stokhos::PseudoSpectralOperator< ordinal_t, value_t, point_compare_type >.
Definition at line 226 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlinevirtual |
Transform values at quadrature points to PCE coefficients.
input
is a vector storing values of a function at the quadrature points, and result
will contain the resulting polynomial chaos coefficients. input
and result
can have multiple columns for vector-valued functions and set trans
to true if these (multi-) vectors are layed out in a transposed fashion.
Implements Stokhos::PseudoSpectralOperator< ordinal_t, value_t, point_compare_type >.
Definition at line 236 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlinevirtual |
Transform PCE coefficients to quadrature values.
input
is a vector storing polynomial chaos coefficients and result
will contain the resulting values at the quadrature points. input
and result
can have multiple columns for vector-valued functions and set trans
to true if these (multi-) vectors are layed out in a transposed fashion.
Implements Stokhos::PseudoSpectralOperator< ordinal_t, value_t, point_compare_type >.
Definition at line 256 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlineprotected |
Apply transformation operator using direct method.
Definition at line 280 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
inlineprotected |
Apply tranformation operator using PST method.
For k=1,...,d, let A_k be the m_k-by-n_k quadrature operator for dimension k and A = A_1 ... A_d be m--n where m = m_1...m_d and n = n_1...n_d. For any two matrices B and C and vector x, (B C)x = vec( C*X*B^T ) = vec( C*(B*X^T)^T ) where x = vec(X), and the vec() operator makes a vector out of a matrix by stacking columns. Applying this formula recursively to A yields the simple algorithm for computing y = A*x (x is n-by-1 and y is m-by-1): X = x; for k=1:d n = n / n_k; X = reshape(X, n, n_k); X = A_k*X'; n = n * m_k; end y = reshape(X, m, 1);
When x has p columns, it is somehwat more complicated because the standard transpose above isn't correct. Instead the transpose needs to be applied to each p block in X.
The strategy here for dealing with transposed input and result is to transpose them when copying to/from the temporary buffers used in the algorithm.
Definition at line 334 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
protected |
Use partial-summation-transformation.
Definition at line 466 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
protected |
Number of coefficients.
Definition at line 469 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
protected |
Dimension.
Definition at line 472 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
protected |
Quadrature points.
Definition at line 475 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
protected |
Map index to point term.
Definition at line 478 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
protected |
Do we need to reorder coefficients for PST.
Definition at line 481 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
protected |
Permutation array when reordering for PST.
Definition at line 484 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
protected |
Matrix mapping points to coefficients.
Definition at line 487 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
protected |
Matrix mapping coefficients to points.
Definition at line 490 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
protected |
Matrix mapping points to coefficients for each dimension for PST.
Definition at line 493 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
protected |
Matrix mapping coefficients to points for each dimension for PST.
Definition at line 496 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.
|
protected |
BLAS wrappers.
Definition at line 499 of file Stokhos_TensorProductPseudoSpectralOperator.hpp.