A posteriori error analysis and adaptivity of a space-time finite element method for the wave equation in second order formulation

1. Problem Statement

The paper addresses the numerical approximation of the linear wave equation in its second-order formulation:
$$u'' - \Delta u = f \quad \text{in } J \times D$$
subject to homogeneous Dirichlet boundary conditions and initial conditions $u(0)=u_0, u'(0)=u_1$.

Challenges:

• Spatio-temporal features: Solutions often exhibit highly localized features (traveling fronts, pulses) that require high resolution in specific regions of space and time.
• Computational Cost: Standard methods using quasi-uniform grids lead to excessive degrees of freedom (DoFs) to capture these features.
• CFL Restrictions: Explicit time-stepping methods are limited by Courant-Friedrichs-Lewy (CFL) stability conditions, which restrict time-step sizes.
• Dynamic Meshes: While implicit space-time methods allow for local mesh refinement, handling dynamic mesh modifications (changing the spatial mesh between time steps) introduces discontinuities in the temporal discretization, complicating error analysis.
• Error Norms: Standard Galerkin approximations for the wave equation are known to suffer from suboptimal convergence in the $L^\infty(H^1)$ norm due to initial data approximation errors. Therefore, the authors focus on the $L^\infty(L^2)$ norm, which is more robust for wave problems.

2. Methodology

2.1 The Numerical Scheme

The authors propose a Space-Time Finite Element Method (ST-FEM) with the following characteristics:

• Spatial Discretization: Standard $H^1$-conforming finite elements (continuous piecewise polynomials of degree $p$) on simplicial or box-type meshes. The spatial mesh can change dynamically between time steps.
• Temporal Discretization:
  • Continuous piecewise polynomials of degree $q$ in time.
  • An upwind discontinuous Galerkin (DG) type approximation for the second temporal derivative.
  • This results in a scheme that is continuous in time for static meshes but becomes discontinuous in time when the spatial mesh changes (due to the projection of the solution from one mesh to another).
• Formulation: The method seeks $U \in V_{p,q}$ satisfying a weak form involving a jump term $[U']_{n-1}$ at time nodes to handle the discontinuity caused by mesh changes.

2.2 A Posteriori Error Analysis Framework

The core theoretical contribution is the derivation of fully computable a posteriori error bounds in the $L^\infty(L^2)$ norm. The proof relies on a combination of novel and existing reconstruction techniques:

1. Time Reconstruction ($\hat{V}$):

   • Based on the work of Makridakis and Nochetto for first-order problems, extended here to second-order wave equations.
   • It maps the discrete solution (piecewise polynomials) to a $C^1$ continuous-in-time function.
   • Crucially, it handles the temporal discontinuity introduced by mesh changes by averaging Hermite interpolation basis expansions to preserve the first temporal derivative while lifting the discontinuity.

2. Hermite-type Reconstruction ($I_c$):

   • A non-standard reconstruction operator that modifies coefficients of Hermite basis functions to ensure global $C^0$ continuity in time.
   • Used to control the error propagation caused by mesh modifications.

3. Elliptic Reconstruction ($\hat{U}$):

   • A spatial reconstruction that maps the discrete solution to a function in $H^2(D) \cap H^1_0(D)$.
   • Defined via a Poisson equation involving the discrete Laplacian and residual terms. This allows the use of elliptic regularity estimates to bound spatial errors.

4. Baker's Test Function:

   • Inspired by classical wave equation analysis, a test function $v(t) = \int_t^\xi (u-U) ds$ is used. This technique is essential for deriving $L^\infty(L^2)$ bounds rather than just energy norms.

2.3 Error Estimators

The total error bound is decomposed into five computable estimators:

• $\eta_{init}$: Initial data approximation error.
• $\eta_f$: Source term oscillations (data approximation).
• $\eta_{time}$: Temporal discretization error (related to time-step size and polynomial degree).
• $\eta_{space}$: Spatial discretization error (related to mesh size $h$ and polynomial degree $p$).
• $\eta_{mesh}$: Mesh-change error. This is a novel term specifically designed to quantify the error introduced when the spatial mesh is modified between time steps. It depends on the projection of the solution from the old mesh to the new one.

3. Key Contributions

1. Rigorous Error Bounds for Dynamic Meshes: The paper provides the first rigorous a posteriori error analysis for an arbitrary-order space-time Galerkin method that allows for dynamic mesh modification. Previous works were largely restricted to static meshes or specific low-order schemes.
2. Explicit Constants: The derived error bounds include explicit constants depending on the spatial and temporal polynomial orders ($p, q$) and mesh regularity, making the bounds fully computable and practical.
3. Novel Reconstruction Operators: The development of a specific time reconstruction that handles mesh-induced discontinuities while preserving temporal derivative continuity is a significant mathematical innovation.
4. Adaptive Algorithm: A practical space-time adaptive algorithm is proposed, utilizing the derived estimators to refine/coarsen both the spatial mesh and the time step dynamically.
5. Energy Dissipation Analysis: The authors analyze the energy dissipation properties of the method, showing that while the scheme is formally dissipative (due to the DG nature and mesh changes), the dissipation decreases exponentially as the temporal order $q$ increases.

4. Numerical Results

The authors validated their theory through extensive numerical experiments using the Gridap.jl library in Julia:

• Convergence: For smooth solutions, the error estimators ($\eta_{space}, \eta_{time}$) converge optimally with respect to $h$ and $\tau$. The effectivity indices (ratio of estimator to true error) remain bounded and do not grow significantly with mesh refinement.
• Mesh Modification Impact:
  • Experiments with "switching" meshes confirmed that mesh changes introduce additional error.
  • The $\eta_{mesh}$ estimator successfully captures this error.
  • The choice of projection operator ($P^n_p$) when changing meshes matters: $H^1$ projection yields smaller mesh-change errors compared to $L^2$ projection or nodal interpolation.
• Adaptive Performance:
  • The adaptive algorithm successfully concentrated DoFs around singularities (e.g., discontinuous initial data) and traveling waves.
  • It demonstrated long-time stability, with error growth much slower than the linear dependence seen in uniform refinement methods.
  • The algorithm effectively balanced spatial and temporal errors, adjusting time steps and mesh density automatically.
• Energy Conservation:
  • On static meshes, energy dissipation is negligible for high-order methods ($q \ge 4$).
  • With dynamic mesh changes, dissipation increases but remains manageable, especially for higher temporal orders.

5. Significance and Conclusion

This work bridges a critical gap in the numerical analysis of wave equations by enabling robust, adaptive space-time simulations with dynamic meshing.

• Practicality: The method allows for efficient simulation of wave propagation over large domains and long times by refining the mesh only where and when necessary, avoiding the "curse of dimensionality" associated with uniform grids.
• Theoretical Advancement: It extends the scope of a posteriori error analysis to arbitrary-order space-time methods with mesh changes, providing a theoretical foundation for adaptive algorithms in hyperbolic problems.
• Future Directions: The authors note that the analysis can be extended to variable coefficients, damped wave equations, and local time-stepping strategies. The inclusion of local time-stepping in the adaptive framework is identified as a key area for future research.

In summary, the paper presents a mathematically rigorous and numerically efficient framework for solving the wave equation, offering a powerful tool for problems where localized features and long-time integration are paramount.