Originally, this concept was intermixed with the idea of the vector component. Since the introduction of non-primitive elements, they became different. Take for instance the solution of the mixed Laplacian system with FE_RaviartThomas. There, the first dim
components are the directional derivatives. Since the shape functions are linear combinations of those, they constitute only a single block. The primal function u would be in the second block, but in the dim+1
st component.
In most cases, when you subdivide a matrix or vector into blocks, you do so by creating one block for each vector component. However, this is not always so, and the DoFRenumbering::component_wise function allows to group several vector components into the same block (see, for example, the step-22 or step-31 tutorial programs, as opposed to step-20).
dim
components are the derivatives in each coordinate direction and the last component is the primal function u.Originally, components were not distinguished from blocks, but since the introduction of non-primitive elements, they have to be distinguished. See FiniteElementData::n_components() and the documentation of FiniteElement
However, it turns out that a significant number of 3d meshes cannot satisfy this convention. This is due to the fact that the face convention for one cell already implies something for the neighbor, since they share a common face and fixing it for the first cell also fixes the normal vectors of the opposite faces of both cells. It is easy to construct cases of loops of cells for which this leads to cases where we cannot find orientations for all faces that are consistent with this convention.
For this reason, above convention is only what we call the standard orientation. deal.II actually allows faces in 3d to have either the standard direction, or its opposite, in which case the lines that make up a cell would have reverted orders, and the normal vector would have the opposite direction. You can ask a cell whether a given face has standard orientation by calling cell->face_orientation(face_no)
: if the result is true
, then the face has standard orientation, otherwise its normal vector is pointing the other direction. There are not very many places in application programs where you need this information actually, but a few places in the library make use of this. Note that in 2d, the result is always true
.
The only places in the library where face orientations play a significant role are in the Triangulation and its accessors, and in the QProjector class and its users.
Note that there is no simple relation between shape functions and generalized support points as for regular support points. Instead, FiniteElement defines a couple of interpolation functions doing the actual interpolation.
If a finite element is Lagrangian, generalized support points and support points coincide.
The full reference for this paper is as follows:
Article{BK07,
author = {Wolfgang Bangerth and Oliver Kayser-Herold},
title = {Data Structures and Requirements for hp Finite Element
Software},
journal = {ACM Trans. Math. Softw.},
year = 2009,
volume = 36,
number = 1,
pages = {4/1--4/31}
}
The numerical examples shown in that paper are generated with a slightly modified version of step-27. The main difference to that tutorial program is that various operations in the program were timed for the paper to compare different options and show that methods are really not all that expensive.
This splitting has several advantages, concerning analysis as well as implementation. For the analysis, it means that conformity with certain spaces (FiniteElementData::Conformity), e.g. continuity, is up to the node values. In deal.II, it helps simplifying the implementation of more complex elements like FE_RaviartThomas considerably.
Examples for node functionals are values in support points and moments with respect to Legendre polynomials. Let us give some examples:
Element | Function space | Node values |
---|---|---|
FE_Q, FE_DGQ | Qk | values in support points |
FE_DGP | Pk | moments with respect to Legendre polynomials |
FE_RaviartThomas (2d) | Qk+1,k x Qk,k+1 | moments on edges and in the interior |
FE_RaviartThomasNodal | Qk+1,k x Qk,k+1 | Gauss points on edges(faces) and anisotropic Gauss points in the interior |
Lagrangian elements fill the vector accessed by FiniteElementBase::get_unit_support_points(), such that the function FiniteElementBase::has_support_points() returns true
. Naturally, these support points are on the reference cell. Then, FEValues can be used (in conjunction with a Mapping) to access support points on the actual grid cells.
x=0
and x=1
(and similarly, in higher dimensions at the vertices of the unit square or cube). On the other hand, higher order Lagrangian elements have unit support points also in the interior of the unit line, square, or cube.