Feel++  0.92.0
Public Types | Static Public Attributes
Feel::Dubiner< Dim, RealDim, Degree, NormalizationPolicy, T, StoragePolicy > Class Template Reference

Dubiner polynomial orthonormal basis. More...

#include <dubiner.hpp>

List of all members.

Public Types

typedef DubinerTraits< Dim,
RealDim, Degree,
NormalizationPolicy, T,
StoragePolicy > 
traits_type
typedef Discontinuous continuity_type
Typedefs
typedef Dubiner< Dim, RealDim,
Degree, NormalizationPolicy, T,
StoragePolicy > 
self_type
typedef self_type basis_type
typedef traits_type::value_type value_type
typedef traits_type::convex_type convex_type
typedef
traits_type::reference_convex_type 
reference_convex_type
typedef
traits_type::diff_pointset_type 
diff_pointset_type
typedef traits_type::storage_policy storage_policy
typedef traits_type::matrix_type matrix_type
typedef
traits_type::vector_matrix_type 
vector_matrix_type
typedef
traits_type::matrix_node_type 
matrix_node_type
typedef traits_type::points_type points_type
typedef traits_type::node_type node_type

Public Member Functions

Constructors, destructor
 Dubiner ()
 Dubiner (Dubiner const &d)
 ~Dubiner ()
Operator overloads
self_type const & operator= (self_type const &d)
matrix_type operator() (node_type const &pt) const
matrix_type operator() (points_type const &pts) const
Accessors
size_type size () const
uint16_type degree () const
self_type const & basis () const
bool isNormalized () const
std::string familyName () const

Static Public Attributes

static const uint16_type nDim = traits_type::nDim
static const uint16_type nRealDim = traits_type::nRealDim
static const uint16_type nOrder = traits_type::nOrder
static const uint16_type nConvexOrderDiff = traits_type::nConvexOrderDiff
static const bool is_normalized = traits_type::is_normalized
static const bool isTransformationEquivalent = true
static const bool isContinuous = false
static const bool is_product = true

Methods

matrix_type coeff () const
static matrix_type evaluate (points_type const &__pts)
template<typename AE >
static vector_matrix_type derivate (ublas::matrix_expression< AE > const &__pts)
static matrix_type const & d (uint16_type i)
 derivatives of Dubiner polynomials the derivatives are computed at the nodes of the lattice
static matrix_type const & derivate (uint16_type i)
 derivatives of Dubiner polynomials the derivatives are computed at the nodes of the lattice

Detailed Description

template<uint16_type Dim, uint16_type RealDim, uint16_type Degree, typename NormalizationPolicy = Normalized<true>, typename T = double, template< class > class StoragePolicy = StorageUBlas>
class Feel::Dubiner< Dim, RealDim, Degree, NormalizationPolicy, T, StoragePolicy >

Dubiner polynomial orthonormal basis.

This class represents the Dubiner polynomials up to degree Degree on a simplex in dimension Dim.

The dubiner polynomials in 1D, the segment $[-1;1]$ are defined using Jacobi polynomials as follows: $ \phi_i(x) = P_i^{0,0}(x) $ where $P_i^{0,0}(x)$ is the i-th Jacobi polynomial evaluated at $x \in [-1;1]$ with weights $(0,0)$.

Author:
Christophe Prud'homme
See also:
Robert C. Kirby Algorithm 839: FIAT - A New Paradigm for Computing Finite Element Basis Functions, ACM Trans. Math. Software Vol. 30 No. 4 pp 502-516
M. Dubiner. Spectral methods on triangles and other domains. J. Sci. Comput., 6:345–390, 1991.
G.E. Karniadakis and S.J. Sherwin, ''Spectral/hp Element Methods for CFD,'' Oxford University Press, March 1999.

Member Function Documentation

template<uint16_type Dim, uint16_type RealDim, uint16_type Degree, typename NormalizationPolicy = Normalized<true>, typename T = double, template< class > class StoragePolicy = StorageUBlas>
self_type const& Feel::Dubiner< Dim, RealDim, Degree, NormalizationPolicy, T, StoragePolicy >::basis ( ) const [inline]
Returns:
self as a basis
template<uint16_type Dim, uint16_type RealDim, uint16_type Degree, typename NormalizationPolicy = Normalized<true>, typename T = double, template< class > class StoragePolicy = StorageUBlas>
matrix_type Feel::Dubiner< Dim, RealDim, Degree, NormalizationPolicy, T, StoragePolicy >::coeff ( ) const [inline]

Dubiner polynomials is an orthonormal basis, the coefficients of the polynomials of the basis are the canonical vectors and represented by the identity matrix (lines are polynomials and columns are the polynomial basis )

This function is correct only if we use the Dubiner polynomials as a basis

template<uint16_type Dim, uint16_type RealDim, uint16_type Degree, typename NormalizationPolicy = Normalized<true>, typename T = double, template< class > class StoragePolicy = StorageUBlas>
static matrix_type const& Feel::Dubiner< Dim, RealDim, Degree, NormalizationPolicy, T, StoragePolicy >::d ( uint16_type  i) [inline, static]

derivatives of Dubiner polynomials the derivatives are computed at the nodes of the lattice

  • i index of the derivative (0 : x, 1 : y, 2 : z )
template<uint16_type Dim, uint16_type RealDim, uint16_type Degree, typename NormalizationPolicy = Normalized<true>, typename T = double, template< class > class StoragePolicy = StorageUBlas>
uint16_type Feel::Dubiner< Dim, RealDim, Degree, NormalizationPolicy, T, StoragePolicy >::degree ( ) const [inline]
Returns:
the maximum degree of the Dubiner polynomial to be constructed
template<uint16_type Dim, uint16_type RealDim, uint16_type Degree, typename NormalizationPolicy = Normalized<true>, typename T = double, template< class > class StoragePolicy = StorageUBlas>
static matrix_type const& Feel::Dubiner< Dim, RealDim, Degree, NormalizationPolicy, T, StoragePolicy >::derivate ( uint16_type  i) [inline, static]

derivatives of Dubiner polynomials the derivatives are computed at the nodes of the lattice

  • i index of the derivative (0 : x, 1 : y, 2 : z )
template<uint16_type Dim, uint16_type RealDim, uint16_type Degree, typename NormalizationPolicy = Normalized<true>, typename T = double, template< class > class StoragePolicy = StorageUBlas>
static matrix_type Feel::Dubiner< Dim, RealDim, Degree, NormalizationPolicy, T, StoragePolicy >::evaluate ( points_type const &  __pts) [inline, static]

evaluate the Dubiner polynomials at a set of points __pts

  • __x is a set of points

Referenced by Feel::Dubiner< Dim, Order, Normalized< false >, T, StorageUBlas >::evaluate().

template<uint16_type Dim, uint16_type RealDim, uint16_type Degree, typename NormalizationPolicy = Normalized<true>, typename T = double, template< class > class StoragePolicy = StorageUBlas>
std::string Feel::Dubiner< Dim, RealDim, Degree, NormalizationPolicy, T, StoragePolicy >::familyName ( ) const [inline]
Returns:
the familyName()
template<uint16_type Dim, uint16_type RealDim, uint16_type Degree, typename NormalizationPolicy = Normalized<true>, typename T = double, template< class > class StoragePolicy = StorageUBlas>
bool Feel::Dubiner< Dim, RealDim, Degree, NormalizationPolicy, T, StoragePolicy >::isNormalized ( ) const [inline]
Returns:
true if the Dubiner polynomials are normalized, false otherwise
template<uint16_type Dim, uint16_type RealDim, uint16_type Degree, typename NormalizationPolicy = Normalized<true>, typename T = double, template< class > class StoragePolicy = StorageUBlas>
size_type Feel::Dubiner< Dim, RealDim, Degree, NormalizationPolicy, T, StoragePolicy >::size ( ) const [inline]

Number of polynomials in set