Series: Finite Element Method

Reference Domains and Quadrature

Computing the element stiffness matrices and load vectors from the previous post requires evaluating integrals over each element $T_e$. In a general mesh the elements can differ in position, size, and orientation, so it is impractical to design shape functions and integration rules separately for each one. Instead, every element of a given type is treated as the image of a single fixed reference domain $\hat{T}$ under a mapping $F_e : \hat{T} \to T_e$. Shape functions and quadrature rules are defined once on $\hat{T}$ and then pulled back to the physical element through $F_e$.

The reference domain for each element type is a standard, convenient geometric shape. A line element maps to the interval $[0, 1]$, a triangle to the right triangle with vertices $(0,0)$, $(1,0)$, $(0,1)$, and so on.

Integrals over $\hat{T}$ are evaluated numerically by quadrature: a rule of order $p$ consists of $Q$ sample points $\hat{x}_1, \ldots, \hat{x}_Q \in \hat{T}$ and associated weights $w_1, \ldots, w_Q$ such that

$$\int_{\hat{T}} f(\hat{x})\, \mathrm{d}\hat{x} \approx \sum_{q=1}^{Q} w_q\, f(\hat{x}_q),$$

with the approximation being exact whenever $f$ is a polynomial of degree at most $p$. Given a required polynomial degree, one selects a rule with enough points to integrate it exactly.

In the library, each reference domain is represented by the ReferenceDomain class in refdomains/refdomain.py, with a concrete subclass for each element type. Its quadrature(order) method returns the pair (points, weights) for the requested degree, where points has shape (Q, topological_dimension) and weights has shape (Q,).

Each element type has a dedicated post describing its reference domain and quadrature rules:

• 1D
  • The Line Domain
• 2D
  • The Triangle Domain
  • The Quadrilateral Domain
• 3D
  • The Tetrahedron Domain
  • The Hexahedron Domain