Rheolef  7.1
an efficient C++ finite element environment
basis

finite element method

Description

The class constructor takes a string that characterizes the finite element method, e.g. "P0", "P1" or "P2". The letters characterize the finite element family and the number is the polynomial degree in this family. The finite element also depends upon the reference element, e.g. triangle, square, tetrahedron, etc. See reference_element for more details. For instance, on the square reference element, the "P1" string designates the common Q1 four-nodes finite element.

Available finite elements

Pk

The classical Lagrange finite element family, where k is any polynomial degree greater or equal to 1.

Pkd

The discontinuous Galerkin finite element family, where k is any polynomial degree greater or equal to 0. For convenience, P0d could also be written as P0.

bubble

The product of baycentric coordinates. Only simplicial elements are supported: edge, triangle and tetrahedron.

RTk
RTkd

The vector-valued Raviart-Thomas family, where k is any polynomial degree greater or equal to 0. The RTkd piecewise discontinuous version is fully implemented for the triangle. Its RTk continuous counterpart is still under development.

Bk

The Bernstein finite element family, where k is any polynomial degree greater or equal to 1. This basis was initially introduced by Bernstein (Comm. Soc. Math. Kharkov, 2th series, 1912) and more recently used in the context of finite elements. This feature is still under development.

Sk

The spectral finite element family, as proposed by Sherwin and Karniadakis (Oxford University Press, 2nd ed., 2005). Here, k is any polynomial degree greater or equal to 1. This feature is still under development.

Continuous feature

Some finite element families could support either continuous or discontinuous junction at inter-element boundaries. For instance, the Pk Lagrange finite element family, with k >= 1, has a continuous interelements junction: it defines a piecewise polynomial and globally continuous finite element functional space. Conversely, the Pkd Lagrange finite element family, with k >= 0, has a discontinuous interelements junction: it defines a piecewise polynomial and globally discontinuous finite element functional space.

Options

The basis class supports some options associated to each finite element method. These options are appended to the string bewteen square braces, and are separated by a coma, e.g. "P5[feteke,bernstein]". See basis for some examples.

Lagrange nodes option

equispaced

Nodes are equispaced. This is the default.

fekete

Node are non-equispaced: for high-order polynomial degree, e.g. greater or equal to 3, their placement is fully optimized, as proposed by Taylor, Wingate and Vincent, 2000, SIAM J. Numer. Anal. With this choice, the interpolation error is dramatically decreased for high order polynomials. This feature is still restricted to the triangle reference element and to polynomial degree lower or equal to 19. Otherwise, it fall back to the following warburton node set.

warburton

Node are non-equispaced: for high-order polynomial degree, e.g. greater or equal to 3, their placement is optimized thought an heuristic solution, as proposed by Warburton, 2006, J. Eng. Math. With this choice, the interpolation error is dramatically decreased for high order polynomials. This feature applies to any reference element and polynomial degree.

Raw polynomial basis

The raw (or initial) polynomial basis is used for building the dual basis, associated to degrees of freedom, via the Vandermonde matrix and its inverse. Changing the raw basis do not change the definition of the FEM basis, but only the way it is constructed. There are three possible raw basis:

monomial

The monomial basis is the simplest choice. While it is suitable for low polynomial degree (less than five), for higher polynomial degree, the Vandermonde matrix becomes ill-conditioned and its inverse leads to errors in double precision.

dubiner

The Dubiner basis (see Dubiner, 1991 J. Sci. Comput.) leads to much better condition number. This is the default.

bernstein

The Bernstein basis could also be an alternative raw basis.

Internals: evaluate and interpolate

The basis class defines member functions that evaluates all the polynomial basis functions of a finite element and their derivatives at a point or a set of point.

The basis polynomial functions are evaluated by the evaluate member function. This member function is called during the assembly process of the integrate function.

The interpolation nodes used by the interpolate function are available by the member function hat_node. Conversely, the member function compute_dofs compute the degrees of freedom.

Implementation

This documentation has been generated from file fem/lib/basis.h

The basis class is simply an alias to the basis_basic class

typedef basis_basic<Float> basis;

The basis_basic class provides an interface to a data container:

template<class T>
class basis_basic : public smart_pointer_nocopy<basis_rep<T> >,
public persistent_table<basis_basic<T> > {
public:
// typedefs:
typedef basis_rep<T> rep;
typedef smart_pointer_nocopy<rep> base;
typedef typename rep::size_type size_type;
typedef typename rep::value_type value_type;
typedef typename rep::valued_type valued_type;
// allocators:
basis_basic (std::string name = "");
void reset (std::string& name);
// accessors:
bool is_initialized() const { return base::operator->() != 0; }
size_type degree() const;
std::string family_name() const;
std::string name() const;
size_type ndof (reference_element hat_K) const;
size_type nnod (reference_element hat_K) const;
bool is_continuous() const;
bool is_discontinuous() const;
bool is_nodal() const;
bool is_hierarchical() const;
size_type size() const;
const basis_basic<T>& operator[] (size_type i_comp) const;
bool have_index_parameter() const;
const basis_option& option() const;
const std::string& valued() const;
const piola_fem<T>& get_piola_fem() const;
};
rheolef::basis_basic::have_index_parameter
bool have_index_parameter() const
Definition: basis.h:781
rheolef::basis_basic::reset
void reset(std::string &name)
Definition: basis_rep.cc:78
rheolef::basis_basic::get_piola_fem
const piola_fem< T > & get_piola_fem() const
Definition: basis.h:795
rheolef::basis_basic::is_initialized
bool is_initialized() const
Definition: basis.h:551
rheolef::basis_basic::value_type
rep::value_type value_type
Definition: basis.h:540
rheolef::basis_basic::base
smart_pointer_nocopy< rep > base
Definition: basis.h:538
rheolef::basis
basis_basic< Float > basis
Definition: basis.h:644
rheolef::basis_basic::is_nodal
bool is_nodal() const
Definition: basis.h:731
rheolef::basis_basic::operator[]
const basis_basic< T > & operator[](size_type i_comp) const
Definition: basis.h:774
rheolef::basis_basic::have_compact_support_inside_element
bool have_compact_support_inside_element() const
Definition: basis.h:745
rheolef::basis_basic::is_continuous
bool is_continuous() const
Definition: basis.h:717
rheolef::basis_rep::value_type
T value_type
Definition: basis.h:215
rheolef::basis_rep::size_type
reference_element::size_type size_type
Definition: basis.h:214
rheolef::basis_basic::family_index
size_type family_index() const
Definition: basis.h:675
rheolef::basis_basic::valued
const std::string & valued() const
Definition: basis.h:802
rheolef::basis_basic::basis_basic
basis_basic(std::string name="")
Definition: basis.h:652
rheolef::basis_basic::valued_tag
valued_type valued_tag() const
Definition: basis.h:788
rheolef::basis_basic::name
std::string name() const
Definition: basis.h:682
rheolef::basis_basic::size
size_type size() const
Definition: basis.h:752
rheolef::basis_basic::valued_type
rep::valued_type valued_type
Definition: basis.h:541
rheolef::basis_basic::nnod
size_type nnod(reference_element hat_K) const
Definition: basis.h:703
rheolef::basis_rep::valued_type
space_constant::valued_type valued_type
Definition: basis.h:216
rheolef::basis_basic::reset_family_index
void reset_family_index(size_type k)
Definition: basis_rep.cc:93
rheolef::basis_basic::size_type
rep::size_type size_type
Definition: basis.h:539
rheolef::basis_basic::is_hierarchical
bool is_hierarchical() const
Definition: basis.h:759
rheolef::basis_basic::ndof
size_type ndof(reference_element hat_K) const
Definition: basis.h:696
rheolef::smart_pointer_base< T, details::no_copy< T > >::operator->
const T * operator->() const
Definition: smart_pointer.h:273
rheolef::basis_basic::family_name
std::string family_name() const
Definition: basis.h:668
rheolef::basis_basic::rep
basis_rep< T > rep
Definition: basis.h:537
rheolef::basis_basic::is_discontinuous
bool is_discontinuous() const
Definition: basis.h:724
rheolef::basis_basic::have_continuous_feature
bool have_continuous_feature() const
Definition: basis.h:738
size_type
field::size_type size_type
Definition: branch.cc:425
rheolef::basis_basic::option
const basis_option & option() const
Definition: basis.h:710
rheolef::basis_basic::degree
size_type degree() const
Definition: basis.h:689