Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Stokhos_DerivOrthogPolyExpansion.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Stokhos Package
4 //
5 // Copyright 2009 NTESS and the Stokhos contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 #ifndef STOKHOS_DERIVORTHOGPOLYEXPANSION_HPP
11 #define STOKHOS_DERIVORTHOGPOLYEXPANSION_HPP
12 
14 #include "Stokhos_DerivBasis.hpp"
15 
16 #include "Teuchos_RCP.hpp"
17 #include "Teuchos_Array.hpp"
20 #include "Teuchos_LAPACK.hpp"
21 
22 namespace Stokhos {
23 
25  template <typename ordinal_type, typename value_type>
26  class DerivOrthogPolyExpansion : public OrthogPolyExpansion<ordinal_type, value_type> {
27  public:
28 
30 
37 
40 
42  ordinal_type size() const { return sz; }
43 
46  getBasis() const {return basis; }
47 
50  getTripleProduct() const { return Cijk; }
51 
52  // Operations
53  void unaryMinus(
56 
58  const value_type& x);
60  const value_type& x);
62  const value_type& x);
64  const value_type& x);
65 
66  void plusEqual(
69  void minusEqual(
72  void timesEqual(
75  void divideEqual(
78 
83  const value_type& a,
87  const value_type& b);
92  const value_type& a,
96  const value_type& b);
101  const value_type& a,
105  const value_type& b);
110  const value_type& a,
114  const value_type& b);
115 
130  const value_type& a,
134  const value_type& b);
153  template <typename OpT>
154  void quad(const OpT& quad_func,
164 // void atan2(OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
165 // const OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
166 // const OrthogPolyApprox<ordinal_type, value_type, node_type>& b);
167 // void atan2(OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
168 // const T& a,
169 // const OrthogPolyApprox<ordinal_type, value_type, node_type>& b);
170 // void atan2(OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
171 // const OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
172 // const T& b);
187  const value_type& a,
191  const value_type& b);
196  const value_type& a,
200  const value_type& b);
203 
204  private:
205 
206  // Prohibit copying
208 
209  // Prohibit Assignment
211 
212  protected:
213 
216 
219 
222 
225 
228 
231 
234 
237 
240 
241  protected:
242 
245 
246  struct acos_quad_func {
247  value_type operator() (const value_type& a) const {
248  return std::acos(a);
249  }
250  };
251 
252  struct asin_quad_func {
253  value_type operator() (const value_type& a) const {
254  return std::asin(a);
255  }
256  };
257 
258  struct atan_quad_func {
259  value_type operator() (const value_type& a) const {
260  return std::atan(a);
261  }
262  };
263 
264  struct acosh_quad_func {
265  value_type operator() (const value_type & a) const {
266  return std::log(a+std::sqrt(a*a-value_type(1.0)));
267  }
268  };
269 
270  struct asinh_quad_func {
271  value_type operator() (const value_type& a) const {
272  return std::log(a+std::sqrt(a*a+value_type(1.0)));
273  }
274  };
275 
276  struct atanh_quad_func {
277  value_type operator() (const value_type& a) const {
278  return 0.5*std::log((value_type(1.0)+a)/(value_type(1.0)-a));
279  }
280  };
281 
282  }; // class DerivOrthogPolyExpansion
283 
284 } // namespace Stokhos
285 
287 
288 #endif // STOKHOS_DERIVORTHOGPOLYEXPANSION_HPP
KOKKOS_INLINE_FUNCTION PCE< Storage > sqrt(const PCE< Storage > &a)
void cbrt(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
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 acos(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
Abstract base class for multivariate orthogonal polynomials that support computing double and triple ...
Teuchos::SerialDenseMatrix< ordinal_type, value_type > A
Matrix.
virtual Teuchos::RCP< const Sparse3Tensor< ordinal_type, value_type > > getTripleProduct() const
Get triple product.
void minusEqual(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const value_type &x)
Teuchos::RCP< const OrthogPolyBasis< ordinal_type, value_type > > getBasis() const
Get basis.
Data structure storing a sparse 3-tensor C(i,j,k) in a a compressed format.
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)
ordinal_type size() const
Get expansion size.
void cos(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
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 timesEqual(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const value_type &x)
void asin(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void cosh(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
Teuchos::RCP< const Teuchos::SerialDenseMatrix< ordinal_type, value_type > > Bij
Derivative double-product tensor.
KOKKOS_INLINE_FUNCTION PCE< Storage > acos(const PCE< Storage > &a)
void unaryMinus(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
Abstract base class for orthogonal polynomial-based expansions.
Teuchos::RCP< const Stokhos::DerivBasis< ordinal_type, value_type > > basis
Basis.
void acosh(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 divideEqual(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const value_type &x)
void tan(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void log10(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 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)
DerivOrthogPolyExpansion & operator=(const DerivOrthogPolyExpansion &b)
void atan(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
Teuchos::RCP< const Stokhos::Dense3Tensor< ordinal_type, value_type > > Dijk
Derivative 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 sqrt(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
Othogonal polynomial expansions based on derivative calculations.
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::RCP< const Stokhos::Sparse3Tensor< ordinal_type, value_type > > Cijk
Triple-product tensor.
KOKKOS_INLINE_FUNCTION PCE< Storage > atan(const PCE< Storage > &a)
void tanh(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void asinh(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void abs(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.
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)
void fabs(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
Teuchos::LAPACK< ordinal_type, value_type > lapack
LAPACK wrappers.
Teuchos::Array< ordinal_type > piv
Pivot array.
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)
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 log(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
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)
Data structure storing a dense 3-tensor C(i,j,k).
void sin(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void atanh(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
Teuchos::SerialDenseMatrix< ordinal_type, value_type > B
RHS.
Stokhos::StandardStorage< ordinal_type, value_type > node_type
KOKKOS_INLINE_FUNCTION PCE< Storage > log(const PCE< Storage > &a)
void sinh(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)
KOKKOS_INLINE_FUNCTION PCE< Storage > asin(const PCE< Storage > &a)
ordinal_type solve(ordinal_type s, ordinal_type nrhs)
Solve linear system.
void exp(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)