Error analysis for a Finite Element Discretization of a corotational harmonic map heat flow problem

1. Problem Statement

The paper addresses the Harmonic Map Heat Flow (HMHF) problem, which describes the $L^2$ gradient flow of the Dirichlet energy for vector fields mapping a domain $\Omega$ to the unit sphere $S^2$. The general equation is:
$$\partial_t u = \Delta u + |\nabla u|^2 u$$
subject to the pointwise unit length constraint $|u|=1$.

The authors focus specifically on a corotational case where the domain $\Omega$ is the unit disk in $\mathbb{R}^2$, and the solution possesses rotational symmetry. Using polar coordinates $(r, \psi)$, the vector field $u$ reduces to a scalar function $u(t, r)$ via the ansatz:
$$u(r, \psi, t) = (\cos \psi \sin u, \sin \psi \sin u, \cos u)^T$$
This reduction transforms the original PDE into a one-dimensional symmetric problem on the interval $I=[0,1]$:
$$\partial_t u = \partial_{rr} u + \frac{1}{r}\partial_r u - \frac{\sin(2u)}{2r^2}$$
with boundary conditions $u(t,0)=u(t,1)=0$. This specific case is mathematically significant because it serves as a testbed for analyzing finite-time singularities (blow-up) in the general HMHF problem.

2. Methodology

The authors propose and analyze a numerical scheme combining spatial and temporal discretization:

• Spatial Discretization: A standard $H^1$-conforming finite element method (FEM) is used. The domain is discretized using a uniform grid with mesh size $h$, and piecewise polynomials of degree $k$ ($S^k_h$) are employed.
• Temporal Discretization: A semi-implicit Euler method is applied. The linear diffusion terms are treated implicitly, while the nonlinear term is treated semi-implicitly (using the solution from the previous time step).
  • The discrete scheme at time step $j+1$ seeks $u_h^{j+1} \in S^k_{h,0}$ such that:
    $$\left( \frac{u_h^{j+1} - u_h^j}{\tau}, v_h \right)_{0,r} + (u_h^{j+1}, v_h)_{H^1_r} = \left( \frac{f(u_h^j)}{r^2}, v_h \right)_{0,r}$$
    where $f(u) = u - \frac{1}{2}\sin(2u)$ and the inner products are weighted by $r$ to account for the polar coordinate transformation.
• Key Feature: The semi-implicit treatment ensures that the resulting system at each time step is linear, making it computationally efficient to solve.

3. Key Contributions and Theoretical Analysis

The primary contribution is a rigorous error analysis proving optimal convergence rates under the assumption of smooth solutions. The analysis relies on several novel technical ingredients:

• Hidden Convexity and Energy Splitting:
  The Dirichlet energy $E(u)$ for this problem is not convex. However, the authors identify a "hidden convexity" by splitting the energy:
  $$E(u) = \frac{1}{2}\|u\|_{H^1_r}^2 - \int_0^1 \frac{F(u(r))}{r^2} r dr$$
  where $F(v) = \frac{1}{2}(v^2 - \sin^2 v)$. Crucially, the function $F$ is convex ($F'' \geq 0$). This property is essential for deriving stability estimates and controlling the linearization error.

• Discrete Energy Dissipation:
  The authors prove that the discrete scheme satisfies a discrete energy decay property:
  $$E(u_h^{j+1}) + \text{dissipation terms} \leq E(u_h^j)$$
  This mimics the continuous energy dissipation and ensures the numerical solution remains bounded.

• $L^\infty$ Bounds and Stability:
  A critical challenge in HMHF analysis is controlling the $L^\infty$ norm of the solution to prevent blow-up in the numerical scheme.

  • The authors prove that if the initial energy is sufficiently small, the $L^\infty$ norm of the discrete solution is bounded by $\pi$.
  • They establish a bound for the term $\|v/r\|_{L^\infty}$ for finite element functions, which is necessary to handle the singularity at $r=0$ in the nonlinear term.

• Error Bounds:
  Using the Galerkin projection and the convexity arguments mentioned above, the authors derive optimal error estimates. For a solution with regularity $u \in H^{k+1}$, the error in the $H^1_r$ norm satisfies:
  $$\|u(t_j) - u_h^j\|_{H^1_r} \leq C(\tau + h^k)$$
  This confirms that the method achieves optimal order convergence in both time ($O(\tau)$) and space ($O(h^k)$).

4. Results

• Theoretical: The paper provides a complete proof of stability and optimal convergence rates for the proposed semi-implicit FEM scheme.
• Numerical: The authors present numerical experiments using the Netgen/NGSolve software package.
  • Convergence Rates: Tables 5.1 and 5.2 demonstrate that for linear ($k=1$) and quadratic ($k=2$) elements, the observed convergence rates in the $H^1_r$ norm match the theoretical predictions ($O(h)$ and $O(h^2)$ respectively).
  • Superconvergence: In the $L^2_r$ norm, the numerical results show a convergence order one degree higher than in the $H^1_r$ norm (e.g., $O(h^2)$ for linear elements), which is a common phenomenon in FEM but not covered by the current theoretical analysis.
  • Energy Dissipation: Figure 5.1 confirms that the numerical solution strictly dissipates energy over time, validating Theorem 3.1.

5. Significance

• Bridging Theory and Computation: While the HMHF problem is well-studied mathematically, rigorous error analysis for its numerical discretization, particularly in the corotational case with singular coefficients ($1/r^2$), is scarce. This paper fills that gap.
• Handling Singularities: The analysis successfully addresses the singularity at the origin ($r=0$) inherent in polar coordinates, a common difficulty in such problems.
• Foundation for Blow-up Studies: Since the corotational case is the standard model for studying finite-time blow-up in HMHF, having a stable, convergent, and linear scheme is vital for future numerical investigations into singularity formation.
• Limitations: The authors note that their analysis relies on the specific convex structure of the corotational case. Extending this analysis to the general equivariant case or the Landau-Lifshitz-Gilbert (LLG) equation is currently obstructed by the lack of a similar convex energy splitting.

In summary, this paper provides a robust mathematical foundation for simulating corotational harmonic map heat flows, proving that a simple semi-implicit Euler finite element scheme is stable, energy-dissipative, and optimally convergent.