44#ifndef STOKHOS_DERIVORTHOGPOLYEXPANSION_HPP
45#define STOKHOS_DERIVORTHOGPOLYEXPANSION_HPP
50#include "Teuchos_RCP.hpp"
51#include "Teuchos_Array.hpp"
52#include "Teuchos_SerialDenseMatrix.hpp"
53#include "Teuchos_SerialDenseVector.hpp"
54#include "Teuchos_LAPACK.hpp"
59 template <
typename ordinal_type,
typename value_type>
68 const Teuchos::RCP<
const Teuchos::SerialDenseMatrix<ordinal_type, value_type> >&
Bij,
76 ordinal_type
size()
const {
return sz; }
79 Teuchos::RCP<const OrthogPolyBasis<ordinal_type, value_type> >
83 virtual Teuchos::RCP<const Sparse3Tensor<ordinal_type, value_type> >
121 const value_type& b);
130 const value_type& b);
139 const value_type& b);
148 const value_type& b);
168 const value_type& b);
187 template <
typename OpT>
188 void quad(
const OpT& quad_func,
225 const value_type& b);
234 const value_type& b);
249 Teuchos::RCP< const Stokhos::DerivBasis<ordinal_type, value_type> >
basis;
252 Teuchos::RCP<const Teuchos::SerialDenseMatrix<ordinal_type, value_type> >
Bij;
255 Teuchos::RCP<const Stokhos::Sparse3Tensor<ordinal_type, value_type> >
Cijk;
258 Teuchos::RCP<const Stokhos::Dense3Tensor<ordinal_type, value_type> >
Dijk;
264 Teuchos::SerialDenseMatrix<ordinal_type, value_type>
A;
267 Teuchos::SerialDenseMatrix<ordinal_type, value_type>
B;
270 Teuchos::Array<ordinal_type>
piv;
278 ordinal_type
solve(ordinal_type s, ordinal_type nrhs);
300 return std::log(a+std::sqrt(a*a-value_type(1.0)));
306 return std::log(a+std::sqrt(a*a+value_type(1.0)));
312 return 0.5*std::log((value_type(1.0)+a)/(value_type(1.0)-a));
Data structure storing a dense 3-tensor C(i,j,k).
Abstract base class for multivariate orthogonal polynomials that support computing double and triple ...
Othogonal polynomial expansions based on derivative calculations.
Teuchos::SerialDenseMatrix< ordinal_type, value_type > B
RHS.
ordinal_type size() const
Get expansion size.
void minus(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const OrthogPolyApprox< ordinal_type, value_type, node_type > &b)
Teuchos::SerialDenseMatrix< ordinal_type, value_type > A
Matrix.
void times(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const OrthogPolyApprox< ordinal_type, value_type, node_type > &b)
Teuchos::RCP< const Stokhos::DerivBasis< ordinal_type, value_type > > basis
Basis.
DerivOrthogPolyExpansion(const DerivOrthogPolyExpansion &)
Teuchos::RCP< const Stokhos::Dense3Tensor< ordinal_type, value_type > > Dijk
Derivative Triple-product tensor.
void quad(const OpT &quad_func, OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const OrthogPolyApprox< ordinal_type, value_type, node_type > &b)
void max(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const OrthogPolyApprox< ordinal_type, value_type, node_type > &b)
void unaryMinus(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void plusEqual(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const value_type &x)
Teuchos::RCP< const Stokhos::Sparse3Tensor< ordinal_type, value_type > > Cijk
Triple-product tensor.
void sinhcosh(OrthogPolyApprox< ordinal_type, value_type, node_type > &s, OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void derivative(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void divide(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const OrthogPolyApprox< ordinal_type, value_type, node_type > &b)
void minusEqual(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const value_type &x)
void timesEqual(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const value_type &x)
void min(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const OrthogPolyApprox< ordinal_type, value_type, node_type > &b)
ordinal_type sz
Workspace size.
void pow(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const OrthogPolyApprox< ordinal_type, value_type, node_type > &b)
ordinal_type solve(ordinal_type s, ordinal_type nrhs)
Solve linear system.
DerivOrthogPolyExpansion(const Teuchos::RCP< const DerivBasis< ordinal_type, value_type > > &basis, const Teuchos::RCP< const Teuchos::SerialDenseMatrix< ordinal_type, value_type > > &Bij, const Teuchos::RCP< const Stokhos::Sparse3Tensor< ordinal_type, value_type > > &Cijk, const Teuchos::RCP< const Stokhos::Dense3Tensor< ordinal_type, value_type > > &Dijk)
Constructor.
void plus(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const OrthogPolyApprox< ordinal_type, value_type, node_type > &b)
void divideEqual(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const value_type &x)
Teuchos::RCP< const OrthogPolyBasis< ordinal_type, value_type > > getBasis() const
Get basis.
Teuchos::Array< ordinal_type > piv
Pivot array.
virtual Teuchos::RCP< const Sparse3Tensor< ordinal_type, value_type > > getTripleProduct() const
Get triple product.
Teuchos::RCP< const Teuchos::SerialDenseMatrix< ordinal_type, value_type > > Bij
Derivative double-product tensor.
Stokhos::StandardStorage< ordinal_type, value_type > node_type
Teuchos::LAPACK< ordinal_type, value_type > lapack
LAPACK wrappers.
DerivOrthogPolyExpansion & operator=(const DerivOrthogPolyExpansion &b)
virtual ~DerivOrthogPolyExpansion()
Destructor.
void sincos(OrthogPolyApprox< ordinal_type, value_type, node_type > &s, OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
Class to store coefficients of a projection onto an orthogonal polynomial basis.
Abstract base class for orthogonal polynomial-based expansions.
Data structure storing a sparse 3-tensor C(i,j,k) in a a compressed format.
Specialization for Sacado::UQ::PCE< Storage<...> >
Top-level namespace for Stokhos classes and functions.
value_type operator()(const value_type &a) const
value_type operator()(const value_type &a) const
value_type operator()(const value_type &a) const
value_type operator()(const value_type &a) const
value_type operator()(const value_type &a) const
value_type operator()(const value_type &a) const