A Fully Discrete Surface Finite Element Method for the Navier--Stokes equations on Evolving Surfaces with prescribed Normal Velocity

We analyze two fully time-discrete numerical schemes for the incompressible Navier-Stokes equations posed on evolving surfaces in $\mathbb{R}^3$ with prescribed normal velocity using the evolving surface finite element method (ESFEM). We employ generalized Taylor-Hood finite elements $\mathrm{\mathbf{P}}{k_u}$– $\mathrm{P}{k_{pr}}$– $\mathrm{P}{k_λ}$, $k_u=k{pr}+1 \geq 2$, $k_λ\geq 1$, for the spatial discretization, where the normal velocity constraint is enforced weakly via a Lagrange multiplier $λ$, and a backward Euler discretization for the time-stepping procedure. Depending on the approximation order of $λ$ and weak formulation of the Navier-Stokes equations, we present stability and error analysis for two different discrete schemes, whose difference lies in the geometric information needed. We establish optimal velocity $L^{2}{a_h}$-norm error bounds ($a_h$ an energy norm) for both schemes when $k_λ=k_u$, but only for the more information intensive one when $k_λ=k_u-1$, using iso-parametric and super-parametric discretizations, respectively, with the help of a newly derived surface Ritz-Stokes projection. Similarly, stability and optimal convergence for the pressures is established in an $L^2{L^2}\times L^2_{H_h^{-1}}$-norm ($H_h^{-1}$ a discrete dual space) when $k_λ=k_u$, using a novel Leray time-projection to ensure weakly divergence conformity for our discrete velocity solution at two different time-steps (surfaces). Assuming further regularity conditions for the more information intensive scheme, along with an almost weak divergence conformity result at two different time-steps, we establish optimal $L^2_{L^2}\times L^2_{L^2}$-norm pressure error bounds when $k_λ=k_u-1$, using super-parametric approximation. Simulations verifying our results are provided, along with a comparison test against a penalty approach.

💡 Research Summary

This paper presents a rigorous numerical analysis of the incompressible Navier–Stokes equations posed on evolving closed surfaces in three‑dimensional space, under the assumption that the normal velocity of the surface is prescribed a priori. The authors adopt a Lagrangian, parametric evolving‑surface finite element method (ESFEM) and develop two fully discrete schemes based on backward Euler time stepping. Spatial discretisation employs a generalized Taylor–Hood element triple (\mathbf{P}{k_u})–(P{k_{pr}})–(P_{k_\lambda}) with (k_u = k_{pr}+1 \ge 2) and (k_\lambda \ge 1). The normal‑velocity constraint (u\cdot n_\Gamma = V_\Gamma) is enforced weakly by an additional Lagrange multiplier (\lambda).

The two schemes differ only in the polynomial degree chosen for (\lambda) and consequently in the amount of geometric information required.

• Scheme eNSW_dh (information‑rich) uses (k_\lambda = k_u). In this case the Lagrange multiplier is approximated with the same order as the velocity, allowing the use of iso‑parametric surface meshes.
• Scheme eNSW_ch (information‑light) works for both (k_\lambda = k_u) and (k_\lambda = k_u-1). When (k_\lambda = k_u-1) the geometric approximation error becomes dominant, and the authors resort to super‑parametric meshes to retain optimal convergence.

A central theoretical contribution is the introduction of a surface Ritz–Stokes projection (R_h). This operator projects the continuous solution onto the discrete velocity–pressure space while handling the fact that the material derivative does not commute with the projection. The authors derive sub‑optimal bounds for (|\partial^\circ_h (u - R_h u)|{L^2(\Gamma_h)}) but show that these bounds are sufficient to obtain optimal energy‑norm ((L^2{a_h})) error estimates for the velocity.

Pressure stability is achieved by a novel discrete Leray time‑projection. For each time step (n) the previous velocity (u_h^{n-1}) is projected onto the current surface, yielding (\widehat{u}h^{n-1}) that is weakly divergence‑free with respect to the new geometry. This construction, together with an inverse discrete Stokes operator and a discrete inf‑sup condition for (\lambda), yields stability and optimal error estimates for the pressure in the mixed norm (L^2(L^2)\times L^2(H_h^{-1})) when (k\lambda = k_u). For the case (k_\lambda = k_u-1) the authors prove optimal (L^2(L^2)\times L^2(L^2)) error bounds under additional regularity assumptions, again relying on the super‑parametric discretisation to control geometric errors.

The paper provides a detailed analysis of the geometric perturbation of bilinear forms, new transport formulae for moving surfaces, and a careful treatment of the material derivative on discrete surfaces. The authors also discuss the role of the mean curvature approximation (H_h) and define the quantities (r_u = \min{k_u, k_g-1}) and (b_{r_u} = \min{k_u, k_g}) to quantify the influence of geometry on convergence rates.

Numerical experiments confirm the theoretical predictions. Convergence studies on static spheres, deforming tori, and more complex geometries demonstrate that both schemes achieve the predicted rates for velocity, pressure, and the Lagrange multiplier. Comparisons with a penalty‑based approach for the normal‑velocity constraint show that the Lagrange‑multiplier formulation yields smaller errors for (\lambda) and better pressure stability, especially when the Leray time‑projection is employed.

In conclusion, the work delivers the first complete a‑priori error analysis for Navier–Stokes equations on evolving surfaces using parametric finite elements. The combination of a weakly enforced normal‑velocity constraint, the surface Ritz–Stokes projection, and the discrete Leray time‑projection constitutes a powerful framework that can be extended to more complex surface fluid models, adaptive mesh refinement, and higher‑order time integration schemes.