Rodin/Variational/QuadratureRule.h file

Quadrature rule classes for numerical integration.

This file defines the QuadratureRule classes which implement numerical quadrature (integration) rules for evaluating integrals over mesh elements. Quadrature rules are fundamental to finite element assembly.

Mathematical Foundation

A quadrature rule approximates an integral as a weighted sum:

\[ \int_K f(x) \, dx \approx \sum_{i=1}^{n_q} w_i f(x_i) \]

where:

  • $ x_i $ are quadrature points
  • $ w_i $ are quadrature weights
  • $ n_q $ is the number of quadrature points

Accuracy

A quadrature rule is exact for polynomials up to degree $ p $ if:

\[ \int_K q(x) \, dx = \sum_{i=1}^{n_q} w_i q(x_i) \quad \forall q \in \mathbb{P}_p \]

Common Quadrature Rules

  • Gauss-Legendre: Optimal for polynomial integrands
  • Gauss-Lobatto: Includes element vertices
  • Vertex quadrature: Simple but low accuracy
  • Midpoint rule: Single point at element center

Applications in FEM

  • Evaluating stiffness matrix entries: $ A_{ij} = \int_K \nabla \phi_i \cdot \nabla \phi_j \, dx $
  • Evaluating load vector entries: $ b_i = \int_K f \phi_i \, dx $
  • Computing functionals and error norms

Usage Example

// Quadrature rule automatically selected based on element type
auto qr = QuadratureRule(integrand);
qr.setPolytope(element);
Real integral_value = qr.compute();

Namespaces

namespace Rodin
The Rodin library for Shape and Topology Optimization.
namespace Rodin::Variational
Module which provides the necessary tools for constructing variational problems.

Classes

template<class FunctionDerived>
class Rodin::Variational::QuadratureRule<FunctionBase<FunctionDerived>>
Quadrature rule for integrating functions on mesh polytopes.
template<class FES, class Data>
class Rodin::Variational::QuadratureRule<GridFunction<FES, Data>>
Integration of a GridFunction object.
template<class LHSDerived, class TrialFES, class RHSDerived, class TestFES>
class Rodin::Variational::QuadratureRule<Dot<ShapeFunctionBase<LHSDerived, TrialFES, TrialSpace>, ShapeFunctionBase<RHSDerived, TestFES, TestSpace>>>
Approximation of the integral of the the dot product between a trial shape function and a test shape function.
template<class NestedDerived, class FES>
class Rodin::Variational::QuadratureRule<ShapeFunctionBase<NestedDerived, FES, TestSpace>>
Approximation of the integral of a test shape function.