General concepts » Ciarlet's definition of a finite element

The mathematical foundation of finite elements.

Introduction

The finite element method (FEM) is a powerful technique for numerically solving partial differential equations. At its heart lies the concept of a finite element, which provides a systematic way to construct piecewise-polynomial approximations on a mesh.

The most widely used formal definition of a finite element was introduced by Philippe G. Ciarlet. This definition provides a clean mathematical framework that separates three concerns:

  1. Geometry: on which shape is the element defined?
  2. Approximation: which polynomial space is used?
  3. Degrees of freedom: how are functions in the space uniquely determined?

This separation is precisely mirrored in Rodin's architecture: the Geometry module handles shapes, the Variational module defines polynomial spaces and DOFs, and the Assembly module puts them together.

The Ciarlet Triple

A finite element in the sense of Ciarlet is defined as a triple $ (K, P, \Sigma) $ where:

  • $ K $ is a geometric domain — a compact, connected subset of $ \mathbb{R}^d $ with nonempty interior (e.g. a triangle, tetrahedron, quadrilateral). This is also called the reference element. In Rodin, reference element geometry is defined by the Polytope::Type enumeration (Triangle, Tetrahedron, Quadrilateral, etc.).

  • $ P $ is a finite-dimensional function space defined on $ K $ — typically a space of polynomials of a certain degree. For instance, $ P = \mathbb{P}_k(K) $ , the space of polynomials of degree at most $ k $ , whose dimension is:

    \[ \dim \mathbb{P}_k(\mathbb{R}^d) = \binom{k + d}{d} \]

    For $ d = 2 $ : $ \dim \mathbb{P}_1 = 3 $ , $ \dim \mathbb{P}_2 = 6 $ . For $ d = 3 $ : $ \dim \mathbb{P}_1 = 4 $ , $ \dim \mathbb{P}_2 = 10 $ .

  • $ \Sigma = \{ \sigma_1, \sigma_2, \ldots, \sigma_n \} $ is a set of degrees of freedom (DOFs) — linear functionals $ \sigma_i : P \to \mathbb{R} $ — such that any function $ p \in P $ is uniquely determined by the values $ \sigma_1(p), \ldots, \sigma_n(p) $ .

The key requirement is that the DOFs must form a basis for the dual space $ P' $ . This means $ \dim(P) = n $ , and the map

\[ p \mapsto (\sigma_1(p), \ldots, \sigma_n(p)) \]

is an isomorphism from $ P $ to $ \mathbb{R}^n $ . Equivalently, if $ \sigma_i(p) = 0 $ for all $ i $ , then $ p = 0 $ . This property is called unisolvence.

Basis Functions

Given a Ciarlet triple $ (K, P, \Sigma) $ , there exists a unique set of basis functions $ \{ \phi_1, \ldots, \phi_n \} \subset P $ such that:

\[ \sigma_i(\phi_j) = \delta_{ij} \]

where $ \delta_{ij} $ is the Kronecker delta. This is called the nodal basis of the finite element. These basis functions are the shape functions that are used to represent any function in the finite element space.

Any function $ p \in P $ can then be written as:

\[ p = \sum_{i=1}^{n} \sigma_i(p) \, \phi_i \]

The values $ \sigma_i(p) $ are the degrees of freedom (DOF values) of $ p $ . In Rodin, these DOF values are stored as the coefficient vector of a GridFunction.

Common Finite Elements

P1 (Piecewise Linear) Elements

The simplest and most widely used finite element is the P1 element on a triangle:

  • $ K $ = a triangle with vertices $ a_1, a_2, a_3 $
  • $ P = \mathbb{P}_1(K) $ = polynomials of degree $ \leq 1 $ (i.e. $ p(x, y) = \alpha + \beta x + \gamma y $ ), dimension 3
  • $ \Sigma = \{ \sigma_1, \sigma_2, \sigma_3 \} $ where $ \sigma_i(p) = p(a_i) $ (point evaluations at vertices)

The basis functions are the barycentric coordinates $ \lambda_1, \lambda_2, \lambda_3 $ , satisfying $ \lambda_i(a_j) = \delta_{ij} $ . Each basis function is 1 at one vertex and 0 at the other two, varying linearly between them.

In Rodin, a P1 finite element space is created with:

P1 Vh(mesh);
// Number of DOFs = number of vertices
std::cout << Vh.getSize() << std::endl;

The P1 class constructs one DOF per mesh vertex, producing a space of continuous piecewise linear functions that is $ H^1 $ -conforming.

P0 (Piecewise Constant) Elements

An even simpler element is the P0 element, used for discontinuous quantities:

  • $ K $ = a polytope (triangle, tetrahedron, etc.)
  • $ P = \mathbb{P}_0(K) $ = constant functions, dimension 1
  • $ \Sigma = \{ \sigma_1 \} $ where $ \sigma_1(p) $ is the value of the constant, often defined as the average over $ K $

In Rodin:

P0 Wh(mesh);
// Number of DOFs = number of cells
std::cout << Wh.getSize() << std::endl;

The P0 class constructs one DOF per mesh cell. Functions in P0 are discontinuous across cell boundaries, making this space a subspace of $ L^2(\Omega) $ (not $ H^1 $ ).

Higher-Order Elements (H1<k>)

For polynomial degree $ k \geq 2 $ , additional DOFs are placed on edges and (for $ k \geq 3 $ ) on cell interiors. For example, the P2 element on a triangle has:

  • $ K $ = a triangle
  • $ P = \mathbb{P}_2(K) $ = quadratic polynomials, dimension 6
  • $ \Sigma $ = point evaluations at the 3 vertices and 3 edge midpoints

In Rodin, higher-order elements are provided by the H1<k> class:

// Quadratic (P2) finite element space
H1 Vh(std::integral_constant<size_t, 2>{}, mesh);

The polynomial order $ k $ is a compile-time template parameter.

From Local to Global

The Ciarlet definition describes a local finite element on a single reference cell $ \hat{K} $ . To construct a global finite element space on an entire mesh $ \mathcal{T}_h $ , we:

  1. Define a reference element $ (\hat{K}, \hat{P}, \hat{\Sigma}) $
  2. Map it to each physical cell $ K \in \mathcal{T}_h $ using a geometric transformation $ F_K : \hat{K} \to K $ . In Rodin, the Jacobian of this transformation is accessible through Point::getJacobian().
  3. Assemble local DOFs into global DOFs, enforcing inter-element continuity conditions at shared entities (vertices, edges, faces)

The result is a finite-dimensional subspace $ V_h $ of the function space in which we seek our approximate solution.

The assembly of local DOFs into global DOFs is handled automatically when you create a finite element space. For example, P1 Vh(mesh) constructs the global DOF numbering by identifying vertex-based DOFs on shared mesh vertices.

Conformity

Different types of inter-element continuity give rise to different types of conforming finite element spaces:

ConformityContinuity requirementTypical useRodin class
$ H^1 $ -conformingContinuous across element boundariesScalar PDEs (Poisson, heat)P1, H1
$ H(\mathrm{div}) $ -conformingNormal component continuousMixed methods, fluid flow—
$ H(\mathrm{curl}) $ -conformingTangential component continuousElectromagnetics—
$ L^2 $ (no conformity)No continuity requiredDiscontinuous Galerkin, indicatorsP0

Rodin currently provides $ H^1 $ -conforming elements (P1, H1) and $ L^2 $ elements (P0). The P1 and H1 spaces produce global functions that are continuous across element boundaries: values at shared vertices (and shared edge/face DOFs for higher order) are identified.

Connection to Rodin

In Rodin, the Ciarlet definition is realized through the interaction of several classes:

Ciarlet componentMathematical objectRodin class
$ K $ Reference element geometryPolytope + Polytope::Type
$ P $ Local polynomial spaceDefined internally by P0, P1, H1
$ \Sigma $ Degrees of freedomDOF numbering managed by the finite element space class
$ V_h $ Global finite element spaceThe space object itself (e.g. P1 Vh(mesh))
$ u_h \in V_h $ Discrete functionGridFunction
$ u $ (unknown)Trial functionTrialFunction
$ v $ (weight)Test functionTestFunction
// The mesh provides the cells K and the geometric transformations F_K
Mesh mesh;
mesh = mesh.UniformGrid(Polytope::Type::Triangle, {16, 16});

// P1 defines the local element (K, P1, point evaluations at vertices)
// and assembles them into a global space Vh
P1 Vh(mesh);

// Trial and test functions live in this global space
TrialFunction u(Vh);
TestFunction  v(Vh);

See Also

« Connectivity | General concepts | Finite element spaces »