Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
Public Member Functions | Static Public Member Functions | Public Attributes | Private Member Functions | List of all members
internal::MatrixFreeFunctions::ShapeInfo< Number > Struct Template Reference

#include <deal.II/matrix_free/shape_info.h>

Inheritance diagram for internal::MatrixFreeFunctions::ShapeInfo< Number >:
[legend]

Public Member Functions

 ShapeInfo ()
 
template<int dim, int spacedim, int dim_q>
 ShapeInfo (const Quadrature< dim_q > &quad, const FiniteElement< dim, spacedim > &fe, const unsigned int base_element=0)
 
template<int dim, int spacedim, int dim_q>
void reinit (const Quadrature< dim_q > &quad, const FiniteElement< dim, spacedim > &fe_dim, const unsigned int base_element=0)
 
const UnivariateShapeData< Number > & get_shape_data (const unsigned int dimension=0, const unsigned int component=0) const
 
std::size_t memory_consumption () const
 

Static Public Member Functions

template<int dim, int spacedim>
static bool is_supported (const FiniteElement< dim, spacedim > &fe)
 
static Table< 2, unsigned intcompute_orientation_table (const unsigned int n_points_per_dim)
 

Public Attributes

ElementType element_type
 
std::vector< unsigned intlexicographic_numbering
 
std::vector< UnivariateShapeData< Number > > data
 
::Table< 2, UnivariateShapeData< Number > * > data_access
 
unsigned int n_dimensions
 
unsigned int n_components
 
unsigned int n_q_points
 
unsigned int dofs_per_component_on_cell
 
unsigned int n_q_points_face
 
std::vector< unsigned intn_q_points_faces
 
unsigned int dofs_per_component_on_face
 
::Table< 2, unsigned intface_to_cell_index_nodal
 
::Table< 2, unsigned intface_to_cell_index_hermite
 
::Table< 2, unsigned intface_orientations_dofs
 
::Table< 2, unsigned intface_orientations_quad
 

Private Member Functions

bool check_1d_shapes_symmetric (UnivariateShapeData< Number > &univariate_shape_data)
 
bool check_1d_shapes_collocation (const UnivariateShapeData< Number > &univariate_shape_data) const
 

Detailed Description

template<typename Number>
struct internal::MatrixFreeFunctions::ShapeInfo< Number >

This struct stores a tensor (Kronecker) product view of the finite element and quadrature formula used for evaluation. It is based on a single or a collection of UnivariateShapeData object(s) that describe one-dimensional ingredients, plus some additional information about how these are combined and how indices are laid out in the multi-dimensional case such as the hierarchical -> lexicographic ordering of FE_Q.

Definition at line 349 of file shape_info.h.

Constructor & Destructor Documentation

◆ ShapeInfo() [1/2]

template<typename Number >
internal::MatrixFreeFunctions::ShapeInfo< Number >::ShapeInfo ( )

Empty constructor. Does nothing.

◆ ShapeInfo() [2/2]

template<typename Number >
template<int dim, int spacedim, int dim_q>
internal::MatrixFreeFunctions::ShapeInfo< Number >::ShapeInfo ( const Quadrature< dim_q > &  quad,
const FiniteElement< dim, spacedim > &  fe,
const unsigned int  base_element = 0 
)

Constructor that initializes the data fields using the reinit method.

Member Function Documentation

◆ reinit()

template<typename Number >
template<int dim, int spacedim, int dim_q>
void internal::MatrixFreeFunctions::ShapeInfo< Number >::reinit ( const Quadrature< dim_q > &  quad,
const FiniteElement< dim, spacedim > &  fe_dim,
const unsigned int  base_element = 0 
)

Initializes the data fields. Takes a one-dimensional quadrature formula and a finite element as arguments and evaluates the shape functions, gradients and Hessians on the one-dimensional unit cell. This function assumes that the finite element is derived from a one- dimensional element by a tensor product and that the zeroth shape function in zero evaluates to one.

◆ is_supported()

template<typename Number >
template<int dim, int spacedim>
static bool internal::MatrixFreeFunctions::ShapeInfo< Number >::is_supported ( const FiniteElement< dim, spacedim > &  fe)
static

Return whether an element is supported by MatrixFree.

The following scalar elements are supported:

In the case of vectorial elements, FE_RaviartThomasNodal and FESystem with base elements from the scalar elements listed above are supported.

◆ compute_orientation_table()

template<typename Number >
static Table< 2, unsigned int > internal::MatrixFreeFunctions::ShapeInfo< Number >::compute_orientation_table ( const unsigned int  n_points_per_dim)
static

Compute a table with numbers of re-orientation for all versions of face flips, orientation, and rotation (relating only to 3d elements).

◆ get_shape_data()

template<typename Number >
const UnivariateShapeData< Number > & internal::MatrixFreeFunctions::ShapeInfo< Number >::get_shape_data ( const unsigned int  dimension = 0,
const unsigned int  component = 0 
) const
inline

Return data of univariate shape functions which defines the dimension dimension of tensor product shape functions regarding vector component component of the underlying finite element.

Definition at line 608 of file shape_info.h.

◆ memory_consumption()

template<typename Number >
std::size_t internal::MatrixFreeFunctions::ShapeInfo< Number >::memory_consumption ( ) const

Return the memory consumption of this class in bytes.

◆ check_1d_shapes_symmetric()

template<typename Number >
bool internal::MatrixFreeFunctions::ShapeInfo< Number >::check_1d_shapes_symmetric ( UnivariateShapeData< Number > &  univariate_shape_data)
private

Check whether we have symmetries in the shape values. In that case, also fill the shape_???_eo fields.

◆ check_1d_shapes_collocation()

template<typename Number >
bool internal::MatrixFreeFunctions::ShapeInfo< Number >::check_1d_shapes_collocation ( const UnivariateShapeData< Number > &  univariate_shape_data) const
private

Check whether symmetric 1d basis functions are such that the shape values form a diagonal matrix, i.e., the nodal points are collocated with the quadrature points. This allows for specialized algorithms that save some operations in the evaluation.

Member Data Documentation

◆ element_type

template<typename Number >
ElementType internal::MatrixFreeFunctions::ShapeInfo< Number >::element_type

Encodes the type of element detected at construction. FEEvaluation will select the most efficient algorithm based on the given element type.

Definition at line 356 of file shape_info.h.

◆ lexicographic_numbering

template<typename Number >
std::vector<unsigned int> internal::MatrixFreeFunctions::ShapeInfo< Number >::lexicographic_numbering

Renumbering from deal.II's numbering of cell degrees of freedom to lexicographic numbering used inside the FEEvaluation schemes of the underlying element in the DoFHandler. For vector-valued elements, the renumbering starts with a lexicographic numbering of the first component, then everything of the second component, and so on.

Definition at line 435 of file shape_info.h.

◆ data

template<typename Number >
std::vector<UnivariateShapeData<Number> > internal::MatrixFreeFunctions::ShapeInfo< Number >::data

Stores data of univariate shape functions defining the underlying tensor product finite element.

Definition at line 441 of file shape_info.h.

◆ data_access

template<typename Number >
::Table<2, UnivariateShapeData<Number> *> internal::MatrixFreeFunctions::ShapeInfo< Number >::data_access

Grants access to univariate shape function data of given dimension and vector component. Rows identify dimensions and columns identify vector components.

Definition at line 448 of file shape_info.h.

◆ n_dimensions

template<typename Number >
unsigned int internal::MatrixFreeFunctions::ShapeInfo< Number >::n_dimensions

Stores the number of space dimensions.

Definition at line 453 of file shape_info.h.

◆ n_components

template<typename Number >
unsigned int internal::MatrixFreeFunctions::ShapeInfo< Number >::n_components

Stores the number of vector components of the underlying vector-valued finite element.

Definition at line 459 of file shape_info.h.

◆ n_q_points

template<typename Number >
unsigned int internal::MatrixFreeFunctions::ShapeInfo< Number >::n_q_points

Stores the number of quadrature points in dim dimensions for a cell.

Definition at line 465 of file shape_info.h.

◆ dofs_per_component_on_cell

template<typename Number >
unsigned int internal::MatrixFreeFunctions::ShapeInfo< Number >::dofs_per_component_on_cell

Stores the number of DoFs per cell of the scalar element in dim dimensions.

Definition at line 471 of file shape_info.h.

◆ n_q_points_face

template<typename Number >
unsigned int internal::MatrixFreeFunctions::ShapeInfo< Number >::n_q_points_face

Stores the number of quadrature points per face in dim dimensions.

Definition at line 476 of file shape_info.h.

◆ n_q_points_faces

template<typename Number >
std::vector<unsigned int> internal::MatrixFreeFunctions::ShapeInfo< Number >::n_q_points_faces

Stores the number of quadrature points of a face in dim dimensions for simplex, wedge and pyramid reference cells.

Definition at line 482 of file shape_info.h.

◆ dofs_per_component_on_face

template<typename Number >
unsigned int internal::MatrixFreeFunctions::ShapeInfo< Number >::dofs_per_component_on_face

Stores the number of DoFs per face in dim dimensions.

Definition at line 487 of file shape_info.h.

◆ face_to_cell_index_nodal

template<typename Number >
::Table<2, unsigned int> internal::MatrixFreeFunctions::ShapeInfo< Number >::face_to_cell_index_nodal

For nodal basis functions with nodes located at the boundary of the unit cell, face integrals that involve only the values of the shape functions (approximations of first derivatives in DG) do not need to load all degrees of freedom of the cell but rather only the degrees of freedom located on the face. While it would also be possible to compute these indices on the fly, we choose to simplify the code and store the indirect addressing in this class.

The first table index runs through the faces of a cell, and the second runs through the nodal degrees of freedom of the face, using dofs_per_face entries.

The indices stored in this member variable are as follows. Consider for example a 2d element of degree 3 with the following degrees of freedom in lexicographic numbering:

12 13 14 15
8 9 10 11
4 5 6 7
0 1 2 3

The first row stores the indices on the face with index 0, i.e., the numbers 0, 4, 8, 12, the second row holds the indices 3, 7, 11, 15 for face 1, the third row holds the indices 0, 1, 2, 3 for face 2, and the last (fourth) row holds the indices 12, 13, 14, 15. Similarly, the indices are stored in 3d. (Note that the y faces in 3d use indices reversed in terms of the lexicographic numbers due to the orientation of the coordinate system.)

Note
This object is only filled in case nodal_at_cell_boundaries evaluates to true.

Definition at line 524 of file shape_info.h.

◆ face_to_cell_index_hermite

template<typename Number >
::Table<2, unsigned int> internal::MatrixFreeFunctions::ShapeInfo< Number >::face_to_cell_index_hermite

The face_to_cell_index_nodal provides a shortcut for the evaluation of values on the faces. For Hermite-type basis functions, one can go one step further and also use shortcuts to get derivatives more cheaply where only two layers of degrees of freedom contribute to the derivative on the face. In the lexicographic ordering, the respective indices is in the next "layer" of degrees of freedom as compared to the nodal values. This array stores the indirect addressing of both the values and the gradient.

The first table index runs through the faces of a cell, and the second runs through the pairs of the nodal degrees of freedom of the face and the derivatives, using 2*dofs_per_face entries.

The indices stored in this member variable are as follows. Consider for example a 2d element of degree 3 with the following degrees of freedom in lexicographic numbering:

20 21 22 23 24
15 16 17 18 19
10 11 12 13 14
5 6 7 8 9
0 1 2 3 4

The first row stores the indices for values and gradients on the face with index 0, i.e., the numbers 0, 1, 5, 6, 10, 11, 15, 16, 20, 21, the second row holds the indices 4, 3, 9, 8, 14, 13, 19, 18, 24, 23 for face 1, the third row holds the indices 0, 5, 1, 6, 2, 7, 3, 8, 4, 9 for face 2, and the last (fourth) row holds the indices 20, 15, 21, 16, 22, 17, 23, 18, 24, 19. Similarly, the indices are stored in 3d. (Note that the y faces in 3d use indices reversed in terms of the lexicographic numbers due to the orientation of the coordinate system.)

Note
This object is only filled in case element_type evaluates to tensor_symmetric_hermite.

Definition at line 564 of file shape_info.h.

◆ face_orientations_dofs

template<typename Number >
::Table<2, unsigned int> internal::MatrixFreeFunctions::ShapeInfo< Number >::face_orientations_dofs

For unknowns located on faces, the basis functions are not in the correct order if a face is not in the standard orientation to a given element. This data structure is used to re-order the basis functions to represent the correct order.

Definition at line 572 of file shape_info.h.

◆ face_orientations_quad

template<typename Number >
::Table<2, unsigned int> internal::MatrixFreeFunctions::ShapeInfo< Number >::face_orientations_quad

For interpretation of values at quadrature points, the order of points is not correct if a face is not in the standard orientation to a given element. This data structure is used to re-order the quadrature points to represent the correct order.

Definition at line 580 of file shape_info.h.


The documentation for this struct was generated from the following file: