Rates of convergence of finite element approximations of second-order mean field games with nondifferentiable Hamiltonians

1. Problem Statement

The paper addresses the numerical approximation of Stationary Second-Order Mean Field Games (MFGs) where the Hamiltonian is nondifferentiable.

• Context: Standard MFG theory (Lasry-Lions) typically assumes the Hamiltonian $H(x, p)$ is differentiable with respect to the gradient variable $p$. This differentiability ensures unique optimal feedback controls. However, in many practical applications (e.g., minimal exit-time problems, bang-bang controls), the Hamiltonian is convex but nondifferentiable, leading to non-unique optimal controls.
• Mathematical Formulation: The problem is modeled as a system of coupled equations involving a Partial Differential Inclusion (PDI) rather than a standard PDE:
  $$\begin{aligned} -\nu\Delta u + H(x, \nabla u) &= F[m](x) \quad \text{in } \Omega, \\ -\nu\Delta m - G(x) &\in \text{div}(m \partial_p H(x, \nabla u)) \quad \text{in } \Omega, \end{aligned}$$
  where:
  • $u$ is the value function.
  • $m$ is the player density.
  • $\partial_p H$ denotes the subdifferential of the Hamiltonian (a set-valued map).
  • The inclusion in the second equation implies that the drift term is not uniquely determined but belongs to a set of possible drifts derived from optimal controls.
• Challenges:
  1. Lack of Differentiability: The term involving $\partial_p H$ is discontinuous and only defined implicitly via measurable selections.
  2. Low Regularity: Solutions $m$ may have very low regularity (potentially only in $L^2$ or $H^1$) even in 1D, making standard error analysis difficult.
  3. Non-coercivity: The system lacks standard coercivity properties due to the nonlinear coupling.
  4. Quantitative Control: Previous works established qualitative convergence (convergence without rates) for such problems, but quantitative error bounds were missing.

2. Methodology

The authors employ a regularization approach combined with stabilized Finite Element Methods (FEM) to derive a priori error bounds.

A. Regularization Strategy

To handle the nondifferentiable Hamiltonian, the authors approximate the original PDI system with a regularized PDE system using Moreau–Yosida regularization:

• The Hamiltonian $H$ is replaced by a smooth approximation $H_\lambda$ (where $\lambda \in (0, 1]$ is the regularization parameter).
• $H_\lambda$ is $C^{1,1}$ (continuously differentiable with Lipschitz gradient), allowing the use of standard PDE techniques.
• The regularized system is:
  $$\begin{aligned} -\nu\Delta u_\lambda + H_\lambda(x, \nabla u_\lambda) &= F[m_\lambda](x), \\ -\nu\Delta m_\lambda - \text{div}\left(m_\lambda \frac{\partial H_\lambda}{\partial p}(x, \nabla u_\lambda)\right) &= G(x). \end{aligned}$$
• Key Result: The error between the exact solution $(u, m)$ and the regularized solution $(u_\lambda, m_\lambda)$ is bounded by $O(\lambda^{1/2})$.

B. Finite Element Discretization

• Space: Piecewise linear, $H^1_0$-conforming finite element spaces ($V_k$) on shape-regular simplicial meshes $\mathcal{T}_k$ with mesh size $h_k$.
• Stabilization: A stabilized FEM is used to handle the convection-dominated nature of the KFP equation and the lack of coercivity. The method introduces a numerical diffusion tensor $D_k$ satisfying:
  1. Consistency: $D_k$ vanishes optimally as $h_k \to 0$.
  2. Discrete Maximum Principle (DMP): The stabilization ensures the discrete density remains non-negative, which is crucial for the monotonicity of the coupling operator $F$.
• Discrete System: The regularized PDE system is discretized to find $(u_{k,\lambda}, m_{k,\lambda}) \in V_k \times V_k$.

C. Error Analysis Framework

The total error is decomposed using the triangle inequality into three components:

1. Regularization Error: $\|(u, m) - (u_\lambda, m_\lambda)\|$.
2. Discretization Error of Regularized Problem: $\|(u_\lambda, m_\lambda) - (u_{k,\lambda}, m_{k,\lambda})\|$.
3. Regularization vs. Discrete PDI Error: $\|(u_{k,\lambda}, m_{k,\lambda}) - (u_k, m_k)\|$, where $(u_k, m_k)$ is the solution to the discretized PDI.

The authors establish bounds for each component, specifically leveraging:

• Semismoothness: A property of the Hamiltonian allowing linearization estimates even without differentiability.
• Elliptic Regularity: The regularity of the solution depends on the domain geometry (parameter $\gamma \in (0, 1]$), where $\gamma=1$ for convex domains.
• Duality Arguments: Used to bound the $L^2$ error of the density despite low regularity.

3. Key Contributions

1. First A Priori Error Bounds for Nondifferentiable MFGs: This is the first work to provide explicit convergence rates for FEM approximations of MFG systems with nondifferentiable Hamiltonians.
2. Handling Set-Valued Maps: The analysis rigorously deals with the subdifferential inclusion $\partial_p H$ and measurable selections, proving that the numerical scheme converges to the weak solution of the PDI.
3. Regularization as a Computational Tool: The paper demonstrates that solving the regularized discrete system (which is easier for solvers due to continuity) provides a valid approximation to the original PDI problem with controlled error.
4. Optimal Parameter Choice: The authors derive the optimal relationship between the mesh size $h_k$ and the regularization parameter $\lambda_k$ to minimize the total error.

4. Main Results

The paper proves the following convergence rate for the joint error in the value function (in $H^1$) and density (in $L^2$):

$$\|u - u_k\|_{H^1(\Omega)} + \|m - m_k\|_{L^2(\Omega)} \lesssim h_k^{\gamma/3}$$

• $\gamma$ Definition: The exponent $\gamma \in (0, 1]$ is determined by the elliptic regularity of the Poisson equation on the domain $\Omega$ with $L^2$ source terms.
  • If $\Omega$ is convex, $\gamma = 1$, yielding a rate of $O(h_k^{1/3})$.
  • If $\Omega$ is a non-convex polygon, $\gamma < 1$, resulting in a slower rate.
• Optimal Parameter Scaling: This rate is achieved by choosing the regularization parameter as a function of the mesh size: $\lambda_k \approx h_k^{2\gamma/3}$.
• Corollary for Regularized Solvers: The paper also provides a bound for the error between the exact solution and the regularized discrete solution $(u_{k,\lambda_k}, m_{k,\lambda_k})$, showing that this approach is viable for practical computations.

5. Significance

• Theoretical Advancement: It bridges a gap in the numerical analysis of MFGs by moving from qualitative convergence (existence of a limit) to quantitative convergence (error bounds) for the challenging case of nondifferentiable Hamiltonians.
• Practical Relevance: Many real-world control problems involve nonsmooth costs or constraints (e.g., $L^1$ costs, bang-bang controls). This work provides a theoretical justification for using regularized FEM schemes to solve these problems, offering practitioners a method to control the error introduced by both discretization and regularization.
• Methodological Innovation: The combination of Moreau–Yosida regularization with stabilized FEM and the use of semismoothness properties to handle the lack of differentiability offers a robust framework that could be extended to other classes of variational inequalities and differential inclusions in control theory.

In summary, the paper establishes that despite the mathematical difficulties posed by nondifferentiable Hamiltonians and set-valued drifts, finite element methods can approximate MFG solutions with a provable, albeit sub-optimal ($h^{1/3}$), convergence rate, provided the regularization parameter is tuned correctly relative to the mesh size.