Flow123d  master-f44eb46
Public Member Functions | Protected Member Functions | Protected Attributes | Friends | List of all members
FiniteElement< dim > Class Template Reference

Abstract class for the description of a general finite element on a reference simplex in dim dimensions. More...

#include <finite_element.hh>

Inheritance diagram for FiniteElement< dim >:
Inheritance graph
[legend]
Collaboration diagram for FiniteElement< dim >:
Collaboration graph
[legend]

Public Member Functions

 FiniteElement ()
 Constructor. More...
 
unsigned int n_dofs () const
 Returns the number of degrees of freedom needed by the finite element. More...
 
double shape_value (const unsigned int i, const arma::vec::fixed< dim > &p, const unsigned int comp=0) const
 Calculates the value of the comp-th component of the i-th shape function at the point p on the reference element. More...
 
arma::vec::fixed< dim > shape_grad (const unsigned int i, const arma::vec::fixed< dim > &p, const unsigned int comp=0) const
 Calculates the comp-th component of the gradient of the i-th shape function at the point p on the reference element. More...
 
unsigned int n_components () const
 Returns numer of components of the basis function.
More...
 
const Dofdof (unsigned int i) const
 Returns i -th degree of freedom. More...
 
unsigned int n_space_components (unsigned int spacedim)
 Number of components of FE in a mapped space with dimension spacedim. More...
 
virtual std::vector< arma::vec::fixed< dim+1 > > dof_points () const
 
virtual ~FiniteElement ()
 Destructor. More...
 

Protected Member Functions

void init (bool primitive=true, FEType type=FEScalar)
 Clears all internal structures. More...
 
void setup_components ()
 Initialize vectors with information about components of basis functions. More...
 
virtual UpdateFlags update_each (UpdateFlags flags)
 Decides which additional quantities have to be computed for each cell. More...
 
virtual void compute_node_matrix ()
 Initializes the node_matrix for computing the coefficients of the shape functions in the raw basis of functions_space_. This is done by evaluating the dofs_ for the basis function and by inverting the resulting matrix. More...
 
bool is_primitive () const
 Indicates whether the basis functions have one or more nonzero components (scalar FE spaces are always primitive). More...
 
unsigned int system_to_component_index (unsigned sys_idx) const
 Returns the component index for vector valued finite elements. More...
 
const std::vector< bool > & get_nonzero_components (unsigned int sys_idx) const
 Returns the mask of nonzero components for given basis function. More...
 

Protected Attributes

FEType type_
 Type of FiniteElement. More...
 
bool is_primitive_
 Primitive FE is using componentwise shape functions, i.e. only one component is nonzero for each shape function. More...
 
std::vector< unsigned int > component_indices_
 Indices of nonzero components of shape functions (for primitive FE). More...
 
std::vector< std::vector< bool > > nonzero_components_
 Footprints of nonzero components of shape functions. More...
 
arma::mat node_matrix
 Matrix that determines the coefficients of the raw basis functions from the values at the support points. More...
 
std::shared_ptr< FunctionSpacefunction_space_
 Function space defining the FE. More...
 
std::vector< Dofdofs_
 Set of degrees of freedom (functionals) defining the FE. More...
 

Friends

class FESystem< dim >
 
class FEValues< 3 >
 
class FE_P_disc< dim >
 
class SubDOFHandlerMultiDim
 

Detailed Description

template<unsigned int dim>
class FiniteElement< dim >

Abstract class for the description of a general finite element on a reference simplex in dim dimensions.

The finite element is determined by a function_space_ and a set of dofs_. Further it must be set whether the finite element is_primitive_, which means that even if the functions in function_space_ are vector-valued, the basis functions have only one nonzero component (this is typical for tensor-product FE, e.g. vector-valued polynomial FE, but not for Raviart-Thomas FE).

Description of dofs:

The reference cell consists of lower dimensional entities (points, lines, triangles). Each dof is associated to one of these entities. If the entity is shared by 2 or more neighbouring cells in the mesh then this dof is shared by the finite elements on all of these cells. If a dof is associated to the cell itself then it is not shared with neighbouring cells.

Shape functions:

Sometimes it is convenient to describe the function space using a basis (called the raw basis) that is different from the set of shape functions for the finite element (the actual basis). For this reason we define the support points which play the role of nodal functionals associated to the particular dofs. To convert between the two bases one can use the node_matrix, which is constructed by the method compute_node_matrix(). In the case of non-Lagrangean finite elements the dofs are not associated to nodal functionals but e.g. to derivatives or averages. For that reason we distinguish between the unit support points which are uniquely associated to the dofs and the generalized support points that are auxiliary for the calculation of the dof functionals.

Definition at line 250 of file finite_element.hh.

Constructor & Destructor Documentation

◆ FiniteElement()

template<unsigned int dim>
FiniteElement< dim >::FiniteElement

Constructor.

Definition at line 74 of file finite_element.cc.

◆ ~FiniteElement()

template<unsigned int dim>
virtual FiniteElement< dim >::~FiniteElement ( )
inlinevirtual

Destructor.

Definition at line 309 of file finite_element.hh.

Member Function Documentation

◆ compute_node_matrix()

template<unsigned int dim>
void FiniteElement< dim >::compute_node_matrix
inlineprotectedvirtual

Initializes the node_matrix for computing the coefficients of the shape functions in the raw basis of functions_space_. This is done by evaluating the dofs_ for the basis function and by inverting the resulting matrix.

Reimplemented in FESystem< dim >.

Definition at line 98 of file finite_element.cc.

Here is the caller graph for this function:

◆ dof()

template<unsigned int dim>
const Dof& FiniteElement< dim >::dof ( unsigned int  i) const
inline

Returns i -th degree of freedom.

Definition at line 296 of file finite_element.hh.

◆ dof_points()

template<unsigned int dim>
vector< arma::vec::fixed< dim+1 > > FiniteElement< dim >::dof_points
virtual

Get barycentric coordinates of the points on the reference element associated with the dofs. Used in BDDC for unknown reason.

Reimplemented in FESystem< dim >.

Definition at line 208 of file finite_element.cc.

◆ get_nonzero_components()

template<unsigned int dim>
const std::vector<bool>& FiniteElement< dim >::get_nonzero_components ( unsigned int  sys_idx) const
inlineprotected

Returns the mask of nonzero components for given basis function.

Parameters
sys_idxIndex of basis function.

Definition at line 358 of file finite_element.hh.

◆ init()

template<unsigned int dim>
void FiniteElement< dim >::init ( bool  primitive = true,
FEType  type = FEScalar 
)
protected

Clears all internal structures.

Definition at line 81 of file finite_element.cc.

Here is the caller graph for this function:

◆ is_primitive()

template<unsigned int dim>
bool FiniteElement< dim >::is_primitive ( ) const
inlineprotected

Indicates whether the basis functions have one or more nonzero components (scalar FE spaces are always primitive).

Definition at line 344 of file finite_element.hh.

◆ n_components()

template<unsigned int dim>
unsigned int FiniteElement< dim >::n_components ( ) const
inline

Returns numer of components of the basis function.

Definition at line 292 of file finite_element.hh.

Here is the caller graph for this function:

◆ n_dofs()

template<unsigned int dim>
unsigned int FiniteElement< dim >::n_dofs ( ) const
inline

Returns the number of degrees of freedom needed by the finite element.

Definition at line 262 of file finite_element.hh.

Here is the caller graph for this function:

◆ n_space_components()

template<unsigned int dim>
unsigned int FiniteElement< dim >::n_space_components ( unsigned int  spacedim)

Number of components of FE in a mapped space with dimension spacedim.

Definition at line 178 of file finite_element.cc.

◆ setup_components()

template<unsigned int dim>
void FiniteElement< dim >::setup_components
protected

Initialize vectors with information about components of basis functions.

Definition at line 90 of file finite_element.cc.

Here is the caller graph for this function:

◆ shape_grad()

template<unsigned int dim>
arma::vec::fixed< dim > FiniteElement< dim >::shape_grad ( const unsigned int  i,
const arma::vec::fixed< dim > &  p,
const unsigned int  comp = 0 
) const

Calculates the comp-th component of the gradient of the i-th shape function at the point p on the reference element.

Parameters
iNumber of the shape function.
pPoint of evaluation.
compNumber of vector component.

Definition at line 127 of file finite_element.cc.

Here is the caller graph for this function:

◆ shape_value()

template<unsigned int dim>
double FiniteElement< dim >::shape_value ( const unsigned int  i,
const arma::vec::fixed< dim > &  p,
const unsigned int  comp = 0 
) const

Calculates the value of the comp-th component of the i-th shape function at the point p on the reference element.

Parameters
iNumber of the shape function.
pPoint of evaluation.
compNumber of vector component.

Definition at line 112 of file finite_element.cc.

Here is the caller graph for this function:

◆ system_to_component_index()

template<unsigned int dim>
unsigned int FiniteElement< dim >::system_to_component_index ( unsigned  sys_idx) const
inlineprotected

Returns the component index for vector valued finite elements.

Parameters
sys_idxIndex of shape function.

Definition at line 351 of file finite_element.hh.

◆ update_each()

template<unsigned int dim>
UpdateFlags FiniteElement< dim >::update_each ( UpdateFlags  flags)
inlineprotectedvirtual

Decides which additional quantities have to be computed for each cell.

Parameters
flagsComputed update flags.

Reimplemented in FESystem< dim >.

Definition at line 146 of file finite_element.cc.

Friends And Related Function Documentation

◆ FE_P_disc< dim >

template<unsigned int dim>
friend class FE_P_disc< dim >
friend

Definition at line 388 of file finite_element.hh.

◆ FESystem< dim >

template<unsigned int dim>
friend class FESystem< dim >
friend

Definition at line 388 of file finite_element.hh.

◆ FEValues< 3 >

template<unsigned int dim>
friend class FEValues< 3 >
friend

Definition at line 388 of file finite_element.hh.

◆ SubDOFHandlerMultiDim

template<unsigned int dim>
friend class SubDOFHandlerMultiDim
friend

Definition at line 394 of file finite_element.hh.

Member Data Documentation

◆ component_indices_

template<unsigned int dim>
std::vector<unsigned int> FiniteElement< dim >::component_indices_
protected

Indices of nonzero components of shape functions (for primitive FE).

Definition at line 373 of file finite_element.hh.

◆ dofs_

template<unsigned int dim>
std::vector<Dof> FiniteElement< dim >::dofs_
protected

Set of degrees of freedom (functionals) defining the FE.

Definition at line 388 of file finite_element.hh.

◆ function_space_

template<unsigned int dim>
std::shared_ptr<FunctionSpace> FiniteElement< dim >::function_space_
protected

Function space defining the FE.

Definition at line 385 of file finite_element.hh.

◆ is_primitive_

template<unsigned int dim>
bool FiniteElement< dim >::is_primitive_
protected

Primitive FE is using componentwise shape functions, i.e. only one component is nonzero for each shape function.

Definition at line 370 of file finite_element.hh.

◆ node_matrix

template<unsigned int dim>
arma::mat FiniteElement< dim >::node_matrix
protected

Matrix that determines the coefficients of the raw basis functions from the values at the support points.

Definition at line 382 of file finite_element.hh.

◆ nonzero_components_

template<unsigned int dim>
std::vector<std::vector<bool> > FiniteElement< dim >::nonzero_components_
protected

Footprints of nonzero components of shape functions.

Definition at line 376 of file finite_element.hh.

◆ type_

template<unsigned int dim>
FEType FiniteElement< dim >::type_
protected

Type of FiniteElement.

Definition at line 364 of file finite_element.hh.


The documentation for this class was generated from the following files: