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 »