🔢 mathematics

Quadrature rules for the mass and stiffness matrices of finite elements for the wave equation on the 2-, 3- and 4-simplex

This paper presents new quadrature rules for mass and stiffness matrices on 2-, 3-, and 4-dimensional simplices to facilitate efficient explicit time stepping and reduced computational cost in finite element discretizations of the wave equation.

Original authors: W. A. Mulder

Published 2026-07-07

📖 5 min read🧠 Deep dive

CC BY 4.0

Original authors: W. A. Mulder

Original paper licensed under CC BY 4.0 (http://creativecommons.org/licenses/by/4.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

Imagine you are trying to simulate how sound waves travel through a complex, bumpy terrain—like a canyon or a building with many rooms. To do this on a computer, scientists break the space down into tiny, simple shapes called simplexes. In 2D, these are triangles; in 3D, they are tetrahedrons (pyramids with triangular bases); and in this paper, the author even tackles 4D shapes (4-simplices), which are the mathematical cousins of those pyramids in a fourth dimension.

The paper by W. A. Mulder is essentially a recipe book for finding the most efficient way to calculate two specific things on these shapes: the Mass Matrix and the Stiffness Matrix.

Here is the breakdown using simple analogies:

1. The Problem: The "Heavy" Calculation

When simulating waves, the computer has to solve a giant system of equations.

The Mass Matrix is like a list of how much "stuff" (mass) is in each tiny shape.
The Stiffness Matrix is like a list of how stiff or rigid the connections are between those shapes.

Usually, calculating these matrices is like trying to weigh a pile of sand grain-by-grain with a super-precise scale. It's accurate, but it takes forever. If you want to simulate waves moving in real-time (explicit time stepping), you need to do these calculations over and over again, very quickly.

2. The Solution: "Mass Lumping" (The Quick-Weigh Trick)

To speed things up, the author uses a trick called Mass Lumping.

The Analogy: Imagine you have a complex, wobbly jelly. Instead of weighing every single molecule, you pretend all the weight is concentrated at specific "anchor points" (nodes) on the jelly.
The Catch: If you just pick random anchor points, your simulation becomes inaccurate (the jelly might wobble the wrong way). To fix this, you need to find the perfect spots to place these anchors and assign them the perfect weights so that the total "feel" of the jelly remains accurate, even though you've simplified the math.

The paper's main goal is to find these perfect anchor points and weights for shapes in 2, 3, and 4 dimensions.

3. The New Discoveries: Finding the Perfect Anchor Points

The author spent a lot of time solving complex mathematical puzzles (polynomial systems) to find these points. Here is what they found:

For 4D Shapes (The 4-Simplex): This is the big news. The author created new rules for 4-dimensional shapes. They found a way to place 20 anchor points on a degree-2 shape (a shape with some curvature) that works perfectly. They also found rules for degree-3 shapes with 80 points.
- Note: The author mentions that while they found a 4D rule with fewer points than some previous attempts, it has a "speed limit" (CFL number) that is lower, meaning the computer simulation has to run slower to stay stable. So, it's a trade-off: fewer points to calculate, but you have to take smaller time steps.
For 3D Shapes (Tetrahedrons): The author found a new version of a 3D pyramid shape that uses 14 points instead of the usual 15.
- The Twist: Just like the 4D case, this new 14-point version is slightly "slower" (has a lower CFL number) than the old 15-point version. It saves a tiny bit of math, but you have to run the simulation more slowly. The author concludes the old 15-point version is still better for practical use.
For 2D Shapes (Triangles): The author didn't invent new triangles but refined the rules for how to calculate the "Stiffness" part. They found ways to use fewer points to calculate stiffness without losing accuracy.

4. The "Stiffness" Shortcut

The paper also tackles the Stiffness Matrix.

The Analogy: Calculating stiffness exactly is like measuring the tension in every single string of a guitar. The author found "numerical quadrature" rules, which are like measuring the tension at just a few key spots and guessing the rest.
The Benefit: These shortcuts require fewer "measurement spots" (nodes) than the exact method. This makes the computer run much faster. The author provides new lists of these shortcut spots for triangles, tetrahedrons, and 4D shapes.

5. The Bottom Line

The paper is a technical manual for mathematicians and engineers who build wave simulation software.

What it achieved: It provided new, verified lists of "anchor points" and weights for 2D, 3D, and 4D shapes.
The Reality Check: While they found ways to use fewer points (which sounds great), the author is honest that sometimes using fewer points forces the simulation to run slower to stay stable.
The Challenge: Finding these points is incredibly hard. It involves solving massive systems of equations. For higher dimensions and more complex shapes, the author admits they couldn't check every single possibility, so there might be even better, simpler rules out there waiting to be discovered.

In short, this paper is about optimizing the math behind wave simulations. It gives us new, slightly different ways to break down space into pieces so computers can calculate wave movements faster, though sometimes with a small trade-off in speed limits.

Technical Summary: Quadrature Rules for Mass and Stiffness Matrices on Simplex Elements

Problem Statement
The paper addresses the computational challenges associated with the finite-element discretization of the wave equation, specifically focusing on the efficiency of explicit time-stepping schemes. While mass lumping (replacing the mass matrix with a diagonal matrix) enables explicit time stepping by avoiding the inversion of large sparse matrices, standard polynomial simplicial elements suffer from a loss of spatial accuracy when lumped. To restore accuracy, higher-degree polynomials are often required in the interior of element faces. Furthermore, the evaluation of stiffness matrices, which involves repeated matrix-vector multiplications, constitutes a significant computational bottleneck. Exact integration of these matrices is costly, and while numerical quadrature can reduce operations, it must be carefully constructed to maintain stability (positive weights) and accuracy (preserving the order of the scheme and avoiding spurious modes).

Methodology
The author employs a systematic approach to construct quadrature rules for both mass and stiffness matrices on 2-, 3-, and 4-dimensional simplices (triangles, tetrahedra, and 4-simplices).

Polynomial Spaces and Node Patterns: The construction utilizes hierarchical polynomial spaces defined by bubble functions on $m$ -faces. The polynomial space $U(p_0, p_1, \dots, p_d)$ ensures face conformity. A necessary condition for unisolvence dictates the node patterns, which are organized into symmetric equivalence classes based on barycentric coordinates.
Mass Matrix Lumping: For the mass matrix, the goal is to find quadrature weights and node positions that are exact for a specific set of polynomials ( $P_{p_1-2} \otimes U$ ). This leads to a system of polynomial equations. The author seeks real, positive solutions where weights correspond to the diagonal entries of the lumped mass matrix.
Stiffness Matrix Quadrature: For stiffness matrices, the paper proposes dedicated numerical quadrature rules that differ from generic rules exact for polynomials up to a given degree. The criteria require the rules to be exact for polynomials in the set $P_{p_1-1} \otimes DU$ (where $DU$ represents partial derivatives of the basis functions) and to be "spurious-free" (ensuring that a zero gradient at all nodes implies a constant solution).
Optimization and Selection: When polynomial systems yield multiple or infinite solutions, selection criteria are applied:
- For mass matrices: Maximizing the Courant-Friedrichs-Lewy (CFL) number to allow larger time steps.
- For stiffness matrices: Minimizing the error $\varepsilon$ between the exact stiffness matrix and the numerically integrated approximation, measured via the Frobenius norm.
Computational Tools: The systems were generated and solved using Mathematica® with extended precision, though results are presented with 24 digits.

Key Contributions and Results

4-Simplex (4D) Extensions: The paper successfully extends mass lumping to 4-dimensional space.
- A degree-2 element with 20 nodes is presented for the 4-simplex.
- A degree-3 element with 80 nodes for the mass matrix and 46 nodes for the stiffness matrix is constructed.
- New dedicated quadrature rules for the stiffness matrix are provided for these 4D elements, often requiring fewer nodes than generic rules.
3-Simplex (Tetrahedron) Improvements:
- A new degree-2 tetrahedral element with 14 nodes is introduced. However, the author notes that while it has fewer nodes, its CFL number is significantly smaller than that of the known 15-node element, making it less attractive for practical applications despite lower computational cost.
- Several new 14-node quadrature rules for the stiffness matrix are listed, with one optimized rule ( $s_{1,2}$ ) showing a significantly smaller error compared to existing rules.
2-Simplex (Triangle) Refinements:
- For degree-2 and degree-3 triangles, the paper lists known elements and provides new dedicated stiffness matrix quadrature rules.
- For degree-4 triangles, new rules with 15 and 16 nodes are presented. One 16-node rule is exact, while optimized 15-node rules offer reduced error compared to exact evaluation.
Stiffness Matrix Efficiency: The study demonstrates that dedicated numerical quadrature for stiffness matrices can be more efficient than exact evaluation by requiring fewer integration nodes. In several cases (e.g., the degree-2 triangle and degree-2 tetrahedron), the number of nodes required for the stiffness matrix is lower than that for the mass matrix.

Significance and Claims
The paper claims to provide a necessary extension of mass-lumped finite elements to four dimensions, a capability previously limited to lower dimensions. It emphasizes that the construction of these elements involves solving complex polynomial systems where the existence of positive solutions is not guaranteed, particularly for higher degrees and dimensions.

The author modestly notes that while new rules were found, the exhaustive treatment of all possible cases was impossible due to computational complexity. Consequently, the paper does not claim to have found the absolute optimal rules for all cases but rather presents a collection of viable candidates. The significance lies in providing specific, tested quadrature rules that balance computational efficiency (fewer nodes) with the stability and accuracy requirements of wave equation simulations. The work also highlights the trade-off between node count and the CFL number, suggesting that a reduction in nodes does not always yield a more efficient element if the time-step restriction becomes too severe.