# Numerical analysis of the large deviation regime of a kinetic equation with a nonlocal Hamilton-Jacobi limit

Hélène Hivert, Tino Laidin

<figure><img src="image_1.png" /></figure>

## To cite this version:

Hélène Hivert, Tino Laidin. Numerical analysis of the large deviation regime of a kinetic equation with a nonlocal Hamilton-Jacobi limit. 2025. hal-05251400

HAL Id: hal-05251400

https://hal.science/hal-05251400v1

Preprint submitted on 12 Sep 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
NUMERICAL ANALYSIS OF THE LARGE DEVIATION REGIME OF A KINETIC EQUATION WITH A NONLOCAL HAMILTON-JACOBI LIMIT

HÉLÈNE HIVERT AND TINO LAIDIN

ABSTRACT. We develop and study an asymptotic-preserving (AP) numerical scheme for a linear kinetic equation in a large deviation regime. After applying a Hopf-Cole transform to the distribution function, the system exhibits the behavior of rare events, which in the limit is governed by a non-standard, nonlocal Hamilton-Jacobi equation, as identified in [E. Bouin et al., J. Lond. Math. Soc., II. Ser., 2023].

The proposed scheme efficiently handles the stiffness introduced by scaling, with a computational cost that remains uniform with respect to the small parameter. It takes advantage of the conservation properties of the original kinetic model to overcome the numerical challenges posed by stiffness. The scheme satisfies a discrete maximum principle, preserves equilibrium states, and correctly captures the asymptotic limit, recovering the viscosity solution of the limit nonlocal Hamilton-Jacobi equation.

As the limit problem is non-standard, convergence results from the literature are not directly applicable. We introduce new analytical tools based on a discrete representation formula that links the numerical scheme with the continuous setting. This allows us to prove the convergence and establish key structural properties of the method. Numerical tests support the analysis and illustrate the robustness of the scheme and the original behavior of the limit system.

KEYWORDS: Asymptotic preserving scheme, Kinetic equation, BGK equation, nonlocal Hamilton-Jacobi equation.

2020 MATHEMATICS SUBJECT CLASSIFICATION: 65M12, 35B40, 82C40, 35F21, 35D40, (Primary) 35Q92, (secondary)

# 1. INTRODUCTION

In this paper, we introduce and analyze a new Asymptotic Preserving (AP) scheme for a linear kinetic equation under a scaling tailored to capture rare events in the dynamics. Specifically, we consider the linear BGK equation for the unknown distribution function $f(t, x, v)$:

$$ \partial_t f(t, x, v) + v \cdot \nabla_x f(t, x, v) = \rho(t, x) \mathcal{M}(v) - f(t, x, v), \quad (t, x, v) \in \mathbb{R}^+ \times \mathbb{R}^{d_x} \times \mathbb{R}^{d_v}, \quad (1.1) $$

where $d_x, d_v \in \mathbb{N}$ denote the spatial and velocity dimensions. The macroscopic density is given by $\rho(t, x) = \langle f(t, x, v) \rangle := \int_{\mathbb{R}^{d_v}} f(t, x, v) dv$, and the equilibrium distribution $\mathcal{M}$ is defined as

$$ \mathcal{M}(v) = \frac{1}{(2\pi)^{d_v/2}} \exp\left(-\frac{|v|^2}{2}\right). $$

The equation is supplemented with an initial data $f_{in}$, which will be specified later. The function $f(t, x, v)$ represents the probability density of finding an individual at time $t$, position $x$, and moving at velocity $v$. The relaxation operator $\rho\mathcal{M} - f$ models interactions that cause individuals to change velocity: they are removed according to the distribution $f$ and reintroduced with a Maxwellian profile $\rho\mathcal{M}$ while preserving mass. This corresponds to a stochastic process in which individuals randomly change velocity following a velocity jump process given by $\rho\mathcal{M} - f$.

In this work, we consider the rescaling, introduced in [12]:

$$ \left(\frac{t}{\varepsilon}, \frac{x}{\varepsilon^{3/2}}, \frac{v}{\varepsilon^{1/2}}\right). \qquad (1.2) $$

Under this scaling, equation (1.1) then reads

$$ \partial_t f^\varepsilon(t, x, v) + v \cdot \nabla_x f^\varepsilon(t, x, v) = \frac{1}{\varepsilon} (\rho^\varepsilon(t, x) \mathcal{M}^\varepsilon(v) - f^\varepsilon(t, x, v)), \quad (t, x, v) \in \mathbb{R}^+ \times \mathbb{R}^{d_x} \times \mathbb{R}^{d_v}, \quad (\mathcal{P}^\varepsilon) $$

where $\mathcal{M}^\varepsilon$ is the centered Gaussian of variance $\varepsilon$:

$$ \mathcal{M}^\varepsilon(v) = \frac{1}{(2\varepsilon\pi)^{d_v/2}} \exp\left(-\frac{|v|^2}{2\varepsilon}\right). $$
Equation (1.1) has been extensively studied in both the hydrodynamic [2, 36] and diffusive [9, 19] regimes, but its analysis under the scaling (1.2) is a recent research subject. One motivation for considering this scaling comes from the following kinetic reaction-transport equation, where individuals also undergo a reproduction process:

$$ \partial_t f + v \cdot \nabla_x f = \rho M - f + r\rho(M - f), \qquad (1.3) $$

where $r > 0$ is a reproduction rate. Equation (1.3) can be viewed as a kinetic analogue of the well-known Fisher-KPP equation. In [13] it was shown that (1.3) exhibits accelerated front propagation, with a velocity that grows like $t^{3/2}$. Due to this acceleration, the classical diffusive scaling is not appropriate to capture the motion of the fronts, which motivates the use of scaling (1.2).

The main challenge in analyzing the limit $\varepsilon \to 0$ in $(\mathcal{P}^\varepsilon)$ lies in the interplay between the stiffness of the right-hand side and the concentration in velocity induced by the resampling through the distribution $\rho^\varepsilon M^\varepsilon$. To accurately capture the asymptotic behavior as $\varepsilon \to 0$, a typical strategy is to study the limit of the Hopf-Cole transform of the unknown. It is defined as

$$ \varphi^{\varepsilon}(t, x, v) = -\varepsilon \ln(f^{\varepsilon}(t, x, v)). \qquad (1.4) $$

This logarithmic transform has, for example, been used in [20, 6, 5, 35, 11, 15] for the study of adaptive dynamics where concentration phenomena typically occur. Another way to interpret (1.4) is that it captures rare events occurring in the tails of the distribution $f^\varepsilon$, which becomes increasingly concentrated as $\varepsilon \to 0$. These rare events are precisely what drive the accelerating propagation fronts in (1.3) [13].

The analysis of $(\mathcal{P}^\varepsilon)$ follows a line of reasoning similar to that used for the classical Fisher-KPP equation [25, 4, 17], based on the so-called geometric optics approximation [27, 25]. In [12], a parallel is drawn with the heat equation under hyperbolic scaling, where the Hopf-Cole transform is shown to converge locally uniformly to the viscosity solution of a Hamilton-Jacobi equation (see, e.g., [17, 3]). In the case of $(\mathcal{P}^\varepsilon)$ the authors showed that $\varphi^\varepsilon$ converges, as $\varepsilon \to 0$, towards a function $\varphi(t, x, v)$ solution to:

$$ \left\{ \begin{array}{l} \max \left( \partial_t \varphi(t, x, v) + v \cdot \nabla_x \varphi(t, x, v) - 1, \varphi(t, x, v) - \min_{v \in \mathbb{R}^{d_v}} \{\varphi(t, x, v)\} - \frac{|v|^2}{2} \right) = 0, \\ \\ \partial_t \left( \min_{v \in \mathbb{R}^{d_v}} \{\varphi\}(t, x) \right) \le 0, \\ \\ \partial_t \left( \min_{v \in \mathbb{R}^{d_v}} \{\varphi\}(t, x) \right) = 0 \quad \text{if} \quad \underset{v \in \mathbb{R}^{d_v}}{\operatorname{argmin}} \{\varphi(t, x, v)\} = \{0\}, \quad (t, x, v) \in \mathbb{R}^+ \times \mathbb{R}^{d_x} \times \mathbb{R}^{d_v}. \end{array} \right. \qquad (\mathcal{P}^0) $$

For simplicity, the initial data of $(\mathcal{P}^\varepsilon)$ is now assumed to be of the form

$$ f_{\text{in}}^{\varepsilon}(x, v) = \exp\left(-\frac{\varphi_{\text{in}}(x, v)}{\varepsilon}\right) $$

so that $(\mathcal{P}^0)$ is supplemented with an initial data $\varphi_{\text{in}}$ such that $\varphi_{\text{in}} - |v|^2/2$ is bounded and satisfies Lipschitz regularity. The precise assumptions are stated later. In addition, we assume that $\varphi_{\text{in}}$ does not depend on $\varepsilon$, but this assumption can be relaxed by assuming local uniform convergence in $\varepsilon$.

They are several key features to $(\mathcal{P}^0)$. First, contrary to standard small-scale asymptotics of kinetic equations, the limit still depends on the velocity variable $v$. Secondly, the limit remains nonlocal in velocity because of the term $\min_{v \in \mathbb{R}^{d_v}} \{\varphi\}(t, x)$. Thirdly, although the system at $\varepsilon > 0$ consists of a single equation for the distribution, the limiting system becomes coupled through a condition involving the minimum in velocity of the unknown. Lastly, a remarkable feature of $(\mathcal{P}^\varepsilon)$ is that the limits $\varepsilon \to 0$ and $t \to \infty$ do not commute: for fixed $\varepsilon > 0$, the system relaxes to a constant equilibrium as $t \to \infty$ [22], whereas in the limit $\varepsilon = 0$, it retains a memory of the initial data even at large times. A schematic representation of this behavior is shown in Figure 1, where the spatial amplitude of the solution is plotted as a function of time. This quantity effectively illustrates the phenomenon: for $\varepsilon$ large, the amplitude converges to 0, whereas in our setting, as $\varepsilon \to 0$, it converges toward a no-zero constant.

The limit system $(\mathcal{P}^0)$ is referred to as a nonlocal Hamilton-Jacobi equation, both by analogy with the classical use of the Hopf-Cole transform for parabolic equations and because it admits a unique viscosity solution, characterized via the method of sub- and super-solutions. It is well known
<figure><img src="image_3.png" /><figcaption>FIGURE 1. Schematic depiction of the long-time spatial amplitude at a fixed ve-
locity of the solution to ($P^\varepsilon$) as $\varepsilon \to 0$.</figcaption></figure>

that solutions to Hamilton-Jacobi equations are, at best, expected to have Lipschitz regularity. This
low regularity presents additional challenges both for the analysis and from a numerical perspective.

It was shown in [12] that ($P^0$) admits a representation formula that characterizes the unique viscosity solution. The rigorous notion of viscosity solutions for ($P^0$) is also developed in that work, based on appropriate definitions of sub- and super-solutions. We refer the reader to [12] for a detailed discussion; the explicit form of the representation formula will be recalled later. For standard Hamilton-Jacobi equations, there is a duality between the Eulerian viewpoint [24, 25, 26], where the dynamics are described by a PDE as in ($P^0$), and the Lagrangian perspective, which focuses on trajectories, as in a representation formula. In the latter case, minimizing trajectories are typically smooth. In contrast, the trajectories associated with the system ($P^0$) are piecewise linear in position and piecewise constant in velocity [12]. The construction of AP schemes offers a promising perspective for extending the analysis of [12] to other systems, such as neutron transport [23], where the accurate treatment of rare events is equally critical.

The design of AP schemes for kinetic equations dates back to [30, 32]. The key idea is to build a scheme for the $\varepsilon > 0$ problem that, in the limit $\varepsilon \to 0$, converges towards a consistent discretization of the asymptotic model. Another fundamental feature of AP schemes is their ability to take time steps independent of the scaling parameter, avoiding a restrictive condition of the type $\Delta t \le \varepsilon \Delta x$, commonly encountered in classical schemes. For an overview of such methods in the kinetic context, we refer to [21, 31]. However, these approaches typically work directly on the distribution function and do not involve the Hopf-Cole transform (1.4).

In the context of biological systems, AP schemes based on logarithmic transformations have been proposed in [16, 28, 1]. To the best of our knowledge, the only instance of such a strategy applied to kinetic equations is found in [29], where (1.3) is treated in the case of bounded velocities to study front propagation. In this work the author considered a multiplicative micro-macro decomposition that allows to treat a macroscopic part, that does not depend on velocity, and a perturbative part separately. Micro-macro schemes are typically additive [8, 33, 34], but the logarithmic transform prompted the multiplicative one allowing the construction of an AP scheme for (1.3).

A major difference between the study in [29] and ours is the fact that they considered bounded velocities. Under this modeling assumption, it is sufficient to rescale (1.3) under a hyperbolic scaling $(t/\epsilon, x/\epsilon, v)$. The asymptotic limit $\epsilon \to 0$ of this problem was studied in [11, 10, 14] and its front propagation property showed in [13]. Because of the bounded velocities, the asymptotic limit does not depend on the velocity variable which is the typical behavior that prompted the use of micro-macro decomposition. For unbounded velocities, which is the novelty of our study compared to [29], it is not the case and such decomposition is not relevant anymore.

The standard approach to discretize Hamilton-Jacobi equations is based on finite differences, for which robust analytical tools have been developed to study convergence towards viscosity solutions [18, 37, 7]. These tools were adapted in [29], where the limiting system admits a classical Hamilton-Jacobi formulation, albeit with a Hamiltonian defined implicitly. However, in our setting, recasting ($P^0$) in Hamiltonian form remains an open problem, necessitating an alternative approach.
The goal of this paper is therefore twofold. First, we construct a numerical scheme capable of capturing the asymptotic regime $\epsilon \to 0$ in $(\mathcal{P}^\epsilon)$ as well as the key features mentioned above. Unlike in [29], the nature of the limit system leads us to develop a new class of AP schemes, based on the Duhamel representation of $(\mathcal{P}^\epsilon)$ and a kind of “meso-macro” framework. In this approach, a macroscopic component evolves along the function $\varphi^\epsilon$. Second, beyond the AP property, we aim to ensure that the limit scheme converges, when the discretization parameters tend to 0, towards the viscosity solution of $(\mathcal{P}^0)$. To this end, we make use of the representation formula of $(\mathcal{P}^0)$, along with a discrete analogue, to establish convergence to the correct solution.

At a deeper level, AP schemes are not just computational tools. They can also guide analytical understanding. While they were initially introduced to efficiently handle stiffness and capture asymptotic limits, their structure can serve as a bridge to the limiting problem itself. In our case, directly discretizing the variational formulation $(\mathcal{P}^0)$ is far from straightforward. By instead designing a scheme that preserves the key features of the scaled kinetic equation $(\mathcal{P}^\epsilon)$, we obtain, in the limit $\epsilon \to 0$, a discretization that naturally approximates the viscosity solution of $(\mathcal{P}^0)$. In addition, the derivation of a discrete representation formula for $(\mathcal{P}^0)$ not only reinforces the consistency of our scheme but also opens a path for extending this approach to other models with similar variational structures.

## Outline of the paper.

In Section 2, we begin by introducing a discretization $(S^\epsilon)$ of $(\mathcal{P}^\epsilon)$ and presenting the main results of the paper. These are the AP property of $(S^\epsilon)$ (Theorem 2.1) and the fact that the asymptotic limit $(S^0)$ captures the viscosity solution of $(\mathcal{P}^0)$ (Theorem 2.2). The theorems are proven in the following sections. Theorem 2.1 is proven in Section 3. It relies on several stability properties of $(S^\epsilon)$. Regarding Theorem 2.2, the main difficulty lies in the fact that classical techniques are not applicable to establish the convergence of $(S^0)$. As in [12], where non-standard features of $(\mathcal{P}^0)$ prompted the introduction of new notions of viscosity solutions, we reformulate $(S^0)$ in Section 4 to connect the discrete setting with its continuous counterpart. Specifically, we analyze stability properties of $(S^0)$ in Section 4.1, and derive a discrete representation formula in Section 4.3. The full convergence result of Theorem 2.2 is proven in Section 5. Finally, Section 6 presents numerical experiments that both confirm the theoretical results and investigate the long time behavior of $(\mathcal{P}^0)$.

## 2. CONSTRUCTION OF THE NUMERICAL SCHEME AND MAIN RESULTS

### 2.1. Construction of the numerical scheme.
We now propose the derivation of an asymptotic preserving scheme for $(\mathcal{P}^\epsilon)$. The main ideas of the discretization are as follows. First, the Duhamel formulation of $(\mathcal{P}^\epsilon)$ is considered. Then, we adopt a meso-macro approach by introducing the Hopf-Cole transform of the distribution $f^\epsilon$. This method allows for an explicit treatment of the stiff term. We proceed in three steps, namely the time, velocity, and space discretizations, which we detail below.

**Time discretization.** Let $N_t \in \mathbb{N}^*$. We set $(0, T] = \cup_{k=0}^{N_t}(t^k, t^{k+1}]$ with $t^n = n\Delta t$ and $t^{N_t} = T$. Our goal is to find a reformulation of $(\mathcal{P}^\epsilon)$ that allows for the derivation of a numerical scheme that is stable in the limit $\epsilon \to 0$. To this achieve this, we start by considering its Duhamel formulation at time $t^{n+1} = t^n + \Delta t$:

$$f^\epsilon(t^n + \Delta t, x, v) = f^\epsilon(t^n, x - \Delta t v, v) e^{-\Delta t/\epsilon} - \frac{1}{\epsilon c^\epsilon} e^{-\Delta t/\epsilon} \int_0^{\Delta t} e^{s/\epsilon} e^{-|v|^2/2\epsilon} \rho^\epsilon(t^n + s, x + (s - \Delta t)v) \, ds, \quad (2.1)$$

with

$$c^\epsilon = \left\langle \exp\left(-\frac{|v|^2}{2\epsilon}\right) \right\rangle. \quad (2.2)$$

In the spirit of [29], we introduce $\varphi^\epsilon$ and $\mu^\epsilon$ defined respectively as the Hopf-Cole transforms of $f^\epsilon$ and $\rho^\epsilon$. Specifically,

$$f^\epsilon(t, x, v) = e^{-\frac{\varphi^\epsilon(t, x, v)}{\epsilon}}, \quad \rho^\epsilon(t, x) = e^{-\frac{\mu^\epsilon(t, x)}{\epsilon}} \quad \text{and} \quad \left\langle e^{-\frac{\varphi^\epsilon(t, x, v)}{\epsilon}} \right\rangle = e^{-\frac{\mu^\epsilon(t, x)}{\epsilon}}. \quad (2.3)$$

The semi-discrete scheme in time arises from selecting a quadrature rule to approximate the integral in (2.1). Due to the stiffness in $\mu^\epsilon$, we adopt a right-point quadrature approximation. Let $\varphi^{\epsilon,n}(x, v)$ denote the approximation of $\varphi^\epsilon(t^n, x, v)$, and similarly, $\mu^{\epsilon,n}(x)$ approximates $\mu^\epsilon(t^n, x)$. The term
$e^{-\mu(t^n+s,x+(s-\Delta t v))/\epsilon}$ is approximated by $e^{-\mu^{n+1}(x)/\epsilon}$ allowing what remains in the integral to be computed explicitly. One then obtains

$$e^{-\varphi^{\epsilon,n+1}(x,v)/\epsilon} = e^{-\varphi^{\epsilon,n}(x-\Delta t v,v)/\epsilon} e^{-\Delta t/\epsilon} + (1 - e^{-\Delta t/\epsilon}) \frac{1}{c^{\epsilon}} e^{-\frac{|v|^2}{2\epsilon}} e^{-\mu^{\epsilon,n+1}(x)/\epsilon}. \quad (2.4)$$

Remark that (2.4) is well-posed, despite its implicit expression. Indeed, because of the conservation of mass of the original kinetic equation ($P^\epsilon$), $\mu^{\epsilon,n+1}$ can be obtained explicitly by integrating (2.4) in velocity and using (2.3), along with the definition (2.2) of $c_\epsilon$. The update of $\mu^{\epsilon,n+1}$ then reads:

$$e^{-\mu^{\epsilon,n+1}(x)/\epsilon} = \left\langle e^{-\varphi^{\epsilon,n}(x-\Delta t v,v)/\epsilon} \right\rangle. \quad (2.5)$$

Recalling that our aim is to build an AP scheme for ($P^\epsilon$), the asymptotic behavior of (2.5) can be understood thanks to the following reformulation:

$$\mu^{\epsilon,n+1}(x) = m_{\varphi}^{\epsilon,n}(x) - \epsilon \ln \left( \left\langle e^{-[\varphi^{\epsilon,n}(x-\Delta t v,v) - m_{\varphi}(x)]/\epsilon} \right\rangle \right), \quad (2.6)$$

where

$$m_{\varphi}^{\epsilon,n}(x) = \min_{v \in \mathbb{R}^{d_v}} \{\varphi^{\epsilon,n}(x - \Delta t v, v)\}. \quad (2.7)$$

Indeed, assuming that one can formally pass to the limit $\epsilon \to 0$ in (2.6), the log term vanishes and the limit of $\mu^{\epsilon,n+1}$ is obtained as the minimum (2.7) of the limit of $\varphi^{\epsilon,n}$. Similarly, we introduce the quantity

$$M^{\epsilon,n}(x,v) = \min \left\{ \varphi^{\epsilon,n}(x - \Delta t v, v) + \Delta t, \frac{|v|^2}{2} + \mu^{\epsilon,n+1}(x) \right\}, \quad (2.8)$$

allowing us to rewrite (2.4) as a semi-discrete in time numerical scheme. For all $x \in \mathbb{R}$, $v \in \mathbb{R}$ and $n \in [1, N_t]$, it reads:

$$\left\{ \begin{array}{l} \mu^{\epsilon,n+1}(x) = m_{\varphi}^{\epsilon,n}(x) - \epsilon \ln \left( \left\langle e^{-[\varphi^{\epsilon,n}(x-\Delta t v,v) - m_{\varphi}^{\epsilon,n}(x)]/\epsilon} \right\rangle \right), \\ \\ \varphi^{\epsilon,n+1}(x,v) = M^{\epsilon,n}(x,v) - \epsilon \ln \left( c^{\epsilon} e^{-[\varphi^{\epsilon,n}(x-\Delta t v,v) + \Delta t - M^{\epsilon,n}(x,v)]/\epsilon} \right. \\ \qquad \left. + (1 - e^{-\Delta t/\epsilon}) e^{-\left[\frac{|v|^2}{2} + \mu^{\epsilon,n+1}(x) - M^{\epsilon,n}(x,v)\right]/\epsilon} \right) + \epsilon \ln(c^{\epsilon}), \end{array} \right. \quad (S_{\Delta t}^{\epsilon})$$

initialized for all $x \in \mathbb{R}$ and all $v \in \mathbb{R}$, by

$$\mu^{\epsilon,0}(x) = \min_{v \in \mathbb{R}^{d_v}} \{\varphi_{\text{in}}(x - \Delta t - v, v)\} \quad \text{and} \quad \varphi^{\epsilon,0}(x, v) = \varphi_{\text{in}}(x, v).$$

**The Mesh.** From now on, we restrict the presentation to the $d_x = d_v = 1$ setting. In order to fit with the continuous setting, we consider unbounded discrete velocities in $\Delta v \mathbb{Z}$ for some $\Delta v > 0$ and set $v_j = j \Delta v$ for $j \in \mathbb{Z}$. We also consider unbounded positions in $\Delta x \mathbb{Z}$ for some $\Delta x > 0$ and set $x_i = i \Delta x$ for $i \in \mathbb{Z}$. The low dimensional setting is chosen to simplify the notations and the extension to higher dimensions will be discussed at the end of Section 5.

**Velocity discretization.** In what follows, we denote by $\varphi_j^{\epsilon,n}(x)$ an approximation of $\varphi^{\epsilon,n}(x, v_j)$. The semi-discretization of the first line of ($S_{\Delta t}^{\epsilon}$) is based on the following velocity quadrature. For a sequence $u = (u_j)_{j \in \mathbb{Z}}$, it is defined by

$$\langle u \rangle_{\Delta v} = \sum_{j \in \mathbb{Z}} u_j \Delta v,$$

and is used to approximate integrals in velocity. It also relies on an approximation of (2.7). For each position $x \in \mathbb{R}$, this is given by

$$m^{\epsilon,n}(x) = \min_{j \in \mathbb{Z}} \{\varphi_j^{\epsilon,n}(x - \Delta t v_j)\}, \quad (2.9)$$

where the minimum is taken over the discrete velocity set ($v_j)_{j \in \mathbb{Z}}$. Similarly, to discretize the second line of ($S_{\Delta t}^{\epsilon}$), we introduce $M_j^{\epsilon,n}(x)$ as an approximation of (2.8):

$$M_j^{\epsilon,n}(x) = \min \left\{ \varphi_j^{\epsilon,n}(x - \Delta t v_j) + \Delta t, \frac{v_j^2}{2} + \mu^{\epsilon,n+1}(x) \right\}. \quad (2.10)$$
The numerical scheme, semi-discrete in time and velocity for a fixed $\epsilon > 0$, is then given, for all $x \in \mathbb{R}$, $j \in \mathbb{Z}$, and $n \in [1, N_t]$, by:

$$ 
\left\{ 
\begin{array}{l}
\mu^{\epsilon,n+1}(x) = m^{\epsilon,n}(x) - \epsilon \ln \left( \left\langle e^{-[\varphi_j^{\epsilon,n}(x-\Delta t v_j) - m^{\epsilon,n}(x)]/\epsilon} \right\rangle_{\Delta v} \right), \\
\\
\varphi_j^{\epsilon,n+1}(x) = M_j^{\epsilon,n}(x) - \epsilon \ln \left( c_{\Delta v}^{\epsilon} e^{-[\varphi_j^{\epsilon,n}(x-\Delta t v_j) + \Delta t - M_j^{\epsilon,n}(x)]/\epsilon} \right. \\
\left. \qquad + \left(1 - e^{-\Delta t/\epsilon}\right) e^{-\left[\frac{v_j^2}{2} + \mu^{\epsilon,n+1}(x) - M_j^{\epsilon,n}(x)\right]/\epsilon} \right) + \epsilon \ln(c_{\Delta v}^{\epsilon}).
\end{array}
\right. 
\quad (S_{\Delta t, \Delta v}^{\epsilon}) 
$$

For all $x \in \mathbb{R}$ and all $j \in \mathbb{Z}$, it is initialized with,

$$ \mu^{\epsilon,0}(x) := \mu_{\text{in}}(x) = \min_{j \in \mathbb{Z}} \{\varphi_{\text{in}}(x - \Delta t v_j, v_j)\} \quad \text{and} \quad \varphi_j^{\epsilon,0}(x) = \varphi_{\text{in}}(x, v_j), \quad (2.11) $$

and the constant $c_{\Delta v}^{\epsilon}$ in the second line of $(S_{\Delta t, \Delta v}^{\epsilon})$ is defined by

$$ c_{\Delta v}^{\epsilon} = \sum_{j} e^{-\frac{v_j^2}{2\epsilon}} \Delta v. $$

Note that, since it is the approximation of the integral of a Gaussian, the use of a first-order discrete integration introduces an error of order $\Delta v$, even though the integral is not on a bounded domain:

$$ c_{\Delta v}^{\epsilon} \underset{\Delta v \to 0, \epsilon \text{ fixed}}{=} \sqrt{2\pi\epsilon} + \mathcal{O}(\Delta v). \quad (2.12) $$

**Position discretization.** Let $i \in \mathbb{Z}$. Denote $\varphi_{ij}^{\epsilon,n}$ and approximation of $\varphi_j^{\epsilon,n}(x_i)$. Remark that in $(S_{\Delta t, \Delta v}^{\epsilon})$, all terms are evaluated at grid points, except for $\varphi_j^{\epsilon,n}(x - \Delta t v_j)$. This term is approximated using a linear interpolation between the two grid points surrounding its argument. Denoting by $\bar{\varphi}_{ij,j}^{\epsilon,n}$ this approximation, it is defined as:

$$ \bar{\varphi}_{ij,j}^{\epsilon,n} = \alpha_j \varphi_{ij,j}^{\epsilon,n} + (1 - \alpha_j) \varphi_{ij-1,j}^{\epsilon,n}, \quad (2.13) $$

where

$$ i_j = i - \beta_j, \quad \beta_j = \left\lfloor v_j \frac{\Delta t}{\Delta x} \right\rfloor, \quad \alpha_j = v_j \frac{\Delta t}{\Delta x} - \beta_j, $$

and $\lfloor \cdot \rfloor$ denotes the floor function. In particular, this interpolation method guarantees that the transport preserves the $L^\infty$-norm and inequalities, just as it does in the continuous setting. Indeed, $\bar{\varphi}_{ij,j}^{\epsilon,n}$ is a convex combination of the two pre-transport values, ensuring that it is bounded by these initial quantities.

With the above notations, we define $\mu_i^{\epsilon,n}$ the approximation of $\mu^{\epsilon,n}(x_i)$, and $m_i^{\epsilon,n}$ the approximation of (2.9) given by

$$ m_i^{\epsilon,n} = \min_{j \in \mathbb{Z}} \{\bar{\varphi}_{ij,j}^{\epsilon,n}\}. \quad (2.14) $$

The quantity (2.10) is approximated by $M_{ij}^{\epsilon,n}$ defined as

$$ M_{ij}^{\epsilon,n} = \min \left\{ \bar{\varphi}_{ij,j}^{\epsilon,n} + \Delta t, \frac{v_j^2}{2} + \mu_i^{\epsilon,n+1} \right\}. \quad (2.15) $$

Finally, for a fixed $\epsilon > 0$, the fully discrete numerical scheme reads for all $(i,j) \in \mathbb{Z}^2$, and for all $n \in [1, N_t]$:

$$ 
\left\{ 
\begin{array}{l}
\mu_i^{\epsilon,n+1} = m_i^{\epsilon,n} - \epsilon \ln \left( \left\langle e^{-[\bar{\varphi}_{ij,j}^{\epsilon,n} - m_i^{\epsilon,n}]/\epsilon} \right\rangle_{\Delta v} \right), \\
\\
\varphi_{ij}^{\epsilon,n+1} = M_{ij}^{\epsilon,n} - \epsilon \ln \left( c_{\Delta v}^{\epsilon} e^{-[\bar{\varphi}_{ij,j}^{\epsilon,n} + \Delta t - M_{ij}^{\epsilon,n}]/\epsilon} \right. \\
\left. \qquad + \left(1 - e^{-\Delta t/\epsilon}\right) e^{-\left[\frac{v_j^2}{2} + \mu_i^{\epsilon,n+1} - M_{ij}^{\epsilon,n}\right]/\epsilon} \right) + \epsilon \ln(c_{\Delta v}^{\epsilon}).
\end{array}
\right. 
\quad (S^{\epsilon}) 
$$

This fully discrete scheme is initialized with

$$ \mu_i^{\epsilon,0} := \mu_{\text{in},i} = \min_{j \in \mathbb{Z}} \{\bar{\varphi}_{\text{in},i,j}\} \quad \text{and} \quad \varphi_{ij}^{\epsilon,0}(x) = \varphi_{\text{in}}(x_i, v_j), \quad (2.16) $$
and the notation $\bar{\varphi}_{\text{in},i,j}$ stand for the linear interpolation (2.13) of $(\varphi_{\text{in}}(x_i, v_j))_{i \in \mathbb{Z}}$ at points $x_i - \Delta t v_j$. Note that the introduction of (2.14) and (2.15) in the construction of the scheme not only provides intuition for the asymptotic behavior of the method, but also plays a crucial role in ensuring numerical stability with respect to $\varepsilon$. Indeed, in this formulation, the exponential always has a non-positive argument, which prevents blow-up in the limit $\varepsilon \to 0$.

**2.2. Main results.** The scheme $(S^\varepsilon)$ presented above is a structure-preserving discretization of $(\mathcal{P}^\varepsilon)$, for which the asymptotic limit can be rigorously justified. Moreover, its asymptotics naturally lead to a consistent and effective numerical scheme for the limit model $(\mathcal{P}^0)$. We emphasize that directly constructing a reliable numerical method for $(\mathcal{P}^0)$ is not straightforward, due to the non-classical nature of the problem. In addition, the convergence of the asymptotic scheme toward the viscosity solution of $(\mathcal{P}^0)$ as the discretization parameters vanish can also be established. These two results constitute the main contributions of this paper and are detailed below. Before stating them, we make the following assumption on the regularity of the initial data :

**Assumption 1.** First, there exists $M > 0$ such that for all $(x, v) \in \mathbb{R}^2$ and all $\varepsilon \in (0, 1]$

$$ \left| \varphi_{\text{in}}(x, v) - \frac{v^2}{2} \right| \leq M \qquad (\text{H1.a}) $$

$$ |\mu_{\text{in}} + \varepsilon \ln(c_{\Delta v}^{\varepsilon})| \leq M. \qquad (\text{H1.b}) $$

Second, for all $v \in \mathbb{R}$, the function $x \mapsto \varphi_{\text{in}}(x, v) - \frac{v^2}{2}$ is L-Lipschitz.

**Theorem 2.1 (Asymptotic Preserving property).** Let $\Delta t > 0$, $\Delta x > 0$, and $\Delta v > 0$ be fixed. For all $(i, j) \in \mathbb{Z}^2$, we define $\mu_i^{\varepsilon,0}$ and $\varphi_{ij}^{\varepsilon,0}$ according to (2.16). Assume moreover that $(\varphi_{\text{in}}, \mu_{\text{in}})$ satisfies Assumption 1. Then, in the limit $\varepsilon \to 0$, the solution $(\mu_i^{\varepsilon,n}, \varphi_{ij}^{\varepsilon,n})$ to $(S^\varepsilon)$ converges locally uniformly in $i \in \mathbb{Z}$ towards $(\mu_i^n, \varphi_{ij}^n)$ solution to

$$ \left\{ \begin{array}{l} \mu_i^{n+1} = \min_{j \in \mathbb{Z}} \{\bar{\varphi}_{ij,j}^n\}, \\ \\ \varphi_{ij}^{n+1} = \min \left\{ \bar{\varphi}_{ij,j}^n + \Delta t, \frac{v_j^2}{2} + \mu_i^{n+1} \right\}, \end{array} \right. \qquad (S^0) $$

with $\bar{\varphi}_{ij,j}^n$ defined as in (2.13), and the initial data is defined for all $j \in \mathbb{Z}$ and $i \in \mathbb{Z}$ by

$$ \mu_i^0 = \min_{j \in \mathbb{Z}} \{\bar{\varphi}_{\text{in},i,j}\} \quad \text{and} \quad \varphi_{ij}^0 = \varphi_{\text{in}}(x_i, v_j). \qquad (2.18) $$

In addition, the numerical cost of $(S^\varepsilon)$ is independent of $\varepsilon$.

As stated in the following theorem, the asymptotic limit $(S^0)$ of $(S^\varepsilon)$ is actually a good discretization of $(\mathcal{P}^0)$ in the sense that it captures the viscosity solution of $(\mathcal{P}^0)$ when the discretization parameters tend to 0.

**Theorem 2.2 (Convergence towards the viscosity solution).** Let $(n, i, j) \in [1, N_t] \times \mathbb{Z}^2$, and suppose that

$$ \Delta t \leq \frac{\Delta v^2}{2}. \qquad (2.19) $$

Assume moreover that the initial data $\varphi_{\text{in}}$ satisfies Assumption 1 as well as L-Lipschitz regularity in velocity. Then, the scheme $(S^0)-(2.18)$ captures the viscosity solution of $(\mathcal{P}^0)$ when $\Delta t$, $\Delta x$, and $\Delta v$ tend to 0 under condition (2.19) and $\frac{\Delta x}{\Delta t} \to 0$. More precisely, there exists a constant C such that

$$ \left| \varphi_{ij}^{N_t} - \varphi(T, x_i, v_j) \right| \leq C \left( \Delta t + \frac{\Delta x}{\Delta t} + \Delta v \right). \qquad (2.20) $$

The constant C depends on the final time T and on the velocity $v_j$.

The proof of these results both rely on stability properties of $(S^\varepsilon)$ and $(S^0)$, namely an equilibrium-preserving property and a maximum principle. They are established in Sections 3 and 4.1 respectively.

Regarding the proof of Theorem 2.2, it unfolds as follows. The convergence of a semi-discrete scheme in time and velocity $(S_{\Delta t, \Delta v}^0)$ towards $(\mathcal{P}^0)$ is first established. In particular, we show that the viscosity solution is indeed captured. This is obtained by introducing a discrete representation formula. This reformulation is equivalent to $(S_{\Delta t, \Delta v}^0)$ under the condition (2.19) whose nature is
<figure><img src="image_4.png" /><figcaption>FIGURE 2. Structure of the convergence proofs</figcaption></figure>

unique to this problem and will be detailed in Section 4.3. The full convergence of the solution to ($S^0$) towards that of ($P^0$) is then established by proving the convergence of ($S^0$) towards ($S^0_{\Delta t, \Delta v}$) through the propagation of interpolation errors. It is worth mentioning that, proving this result directly, without the intermediate step, would require a deep understanding of a fully discrete representation formula to connect with the continuous notion of viscosity solutions—a significantly more complex problem.

The various convergence results in this work are summarized in Figure 2. Although we do not emphasize this aspect, the convergence of ($S^\epsilon$) to ($P^\epsilon$) is expected for smooth solutions. However, a detailed analysis of this convergence lies outside the scope of this article and will instead be illustrated numerically in Section 6.

### 3. AP PROPERTY

We begin by focusing on the semi-discrete scheme ($S^\epsilon_{\Delta t, \Delta v}$), where the position variable remains continuous. This choice allows to highlight the main analytical tools involved in the study of the numerical scheme. Both in this section and in Section 5, the velocity and position variables play fundamentally different roles in the system’s dynamics. This distinction is already present in the original kinetic equation (1.1), where transport and relaxation affect the distribution $f^\epsilon$ in qualitatively different ways.

The extension of the analysis to the fully discrete case is conceptually straightforward, thanks to the interpolation procedure defined in (2.13). However, the notations involved are considerably heavier. For this reason, we do not detail the fully discrete setting here, but some hints are provided at the end of the section.

**Lemma 3.1 (Equilibrium Preserving).** The semi-discrete scheme ($S^\epsilon_{\Delta t, \Delta v}$) preserves the equilibrium ($v_j^2/2$)$_{j \in \mathbb{Z}}$: if $\varphi_{\text{in}}(x, v) = v^2/2$ for all $x \in \mathbb{R}$ and $v \in \mathbb{R}$, then

$$ \forall n \in [1, N_t], \forall j \in \mathbb{Z}, \forall x \in \mathbb{R}, \quad \varphi_j^{\epsilon,n}(x) = \frac{v_j^2}{2}. $$

*Proof.* Let $n \in [1, N_t]$ and assume that $\forall j \in \mathbb{Z}$ and $\forall x \in \mathbb{R}$, $\varphi_j^{\epsilon,n}(x) = \frac{v_j^2}{2}$. Then, $\forall x \in \mathbb{R}$, $\varphi_j^{\epsilon,n}(x - \Delta t v_j) = \frac{v_j^2}{2}$. Moreover, since $v_0 = 0$, definition (2.9) yields $m^{\epsilon,n}(x) = 0$, and

$$ \mu^{\epsilon,n+1}(x) = -\epsilon \ln \left( \left\langle e^{-\frac{v_j^2}{2}/\epsilon} \right\rangle_{\Delta v} \right) = -\epsilon \ln (c_{\Delta v}^{\epsilon}). \quad (3.1) $$

Equation (3.1) is then used to express $M_j^{\epsilon,n}(x)$, defined in (2.10), as

$$ M_j^{\epsilon,n}(x) = \frac{v_j^2}{2} + \min \{\Delta t, -\epsilon \ln(c_{\Delta v}^{\epsilon})\}. \quad (3.2) $$
To determine $\varphi_j^{\varepsilon,n+1}(x)$ and conclude, we distinguish two cases depending on the value of the minimum in (3.2). In both cases, inserting (3.1) and (3.2) into the scheme on $\varphi_j^{\varepsilon,n+1}$ in $(S_{\Delta t, \Delta v}^{\varepsilon})$ and simplifying the expressions yields

$$ \varphi_j^{\varepsilon,n+1}(x) = \frac{v_j^2}{2}, $$

which concludes the proof.

<figure><img src="image_5.png" /></figure>

The next Lemma states the monotonicity of $(S_{\Delta t, \Delta v}^{\varepsilon})$ as well as invariance and commutative properties.

**Lemma 3.2.** Let $(n, j) \in [1, N_t] \times \mathbb{Z}$ and $x \in \mathbb{R}$. Let $(\mu^{\varepsilon,n}(x))_{n \in [1, N_t], j \in \mathbb{Z}}$ and $(\varphi_j^{\varepsilon,n}(x))_{n \in [1, N_t], j \in \mathbb{Z}}$ be a solution to $(S_{\Delta t, \Delta v}^{\varepsilon})$ with initial data (2.11). Denote $(\psi^{\varepsilon,n}(x))_{n \in [1, N_t], j \in \mathbb{Z}}$ another solution to $(S_{\Delta t, \Delta v}^{\varepsilon})$ and

$$ \nu^{\varepsilon,n}(x) = -\varepsilon \ln \left( \langle \psi^{\varepsilon,n}(x - \Delta t v_j) \rangle_{\Delta v} \right). \quad (3.3) $$

Then, the following properties hold:

i) **Monotonicity.** Suppose that $\varphi_{\text{in}} \le \psi_{\text{in}}$. Then, for all $j \in \mathbb{Z}$, $x \in \mathbb{R}$ and $n \in [1, N_t]$,

$$ \mu^{\varepsilon,n}(x) \le \nu^{\varepsilon,n}(x), \quad (3.4) $$

$$ \varphi_j^{\varepsilon,n}(x) \le \psi_j^{\varepsilon,n}(x). \quad (3.5) $$

ii) **Invariance by translation.** Let $y \in \mathbb{R}$ and suppose that $\psi_{\text{in}} = \varphi_{\text{in}}(\cdot + y, \cdot)$. Then, for all $j \in \mathbb{Z}$, $x \in \mathbb{R}$ and $n \in [1, N_t]$,

$$ \nu^{\varepsilon,n}(x) = \mu^{\varepsilon,n}(x + y), \quad (3.6) $$

$$ \psi_j^{\varepsilon,n}(x) = \varphi_j^{\varepsilon,n}(x + y). \quad (3.7) $$

iii) **Commutation with constants.** Let $K \in \mathbb{R}$ and suppose that $\psi_{\text{in}} = \varphi_{\text{in}} + K$. Then, for all $j \in \mathbb{Z}$, $x \in \mathbb{R}$ and $n \in [1, N_t]$,

$$ \nu^{\varepsilon,n}(x) = \mu^{\varepsilon,n}(x) + K, \quad (3.8) $$

$$ \psi_j^{\varepsilon,n}(x) = \varphi_j^{\varepsilon,n}(x) + K. \quad (3.9) $$

*Proof.* For each property, we proceed by induction.

i) Let $n \in [0, N_t - 1]$ and suppose that $\varphi_j^{\varepsilon,n}(x) \le \psi_j^{\varepsilon,n}(x)$ for all $j \in \mathbb{Z}$ and $x \in \mathbb{R}$. Then, one has in particular that $\varphi_j^{\varepsilon,n}(x - \Delta t v_j) \le \psi_j^{\varepsilon,n}(x - \Delta t v_j)$. It implies

$$ \left\langle e^{-\left(\varphi_j^{\varepsilon,n}(x-\Delta t v_j)\right)/\varepsilon} \right\rangle_{\Delta v} \ge \left\langle e^{-\left(\psi_j^{\varepsilon,n}(x-\Delta t v_j)\right)/\varepsilon} \right\rangle_{\Delta v}. $$

Using definitions (2.5) and (3.3), it yields (3.4). Similar computations lead to (3.5) by using that the second line of $(S_{\Delta t, \Delta v}^{\varepsilon})$ can be reformulated as a sum of positive terms as follows:

$$ e^{-\varphi_j^{\varepsilon,n+1}(x)/\varepsilon} = e^{-\varphi_j^{\varepsilon,n}(x-\Delta t v_j)/\varepsilon} e^{-\Delta t/\varepsilon} + (1 - e^{-\Delta t/\varepsilon}) \frac{1}{c_{\Delta v}^{\varepsilon}} e^{-\frac{v_j^2}{2\varepsilon}} e^{-\mu^{\varepsilon,n+1}(x)/\varepsilon}. $$

ii) Let $y \in \mathbb{R}$ and define $\psi_j^{\varepsilon,n}(x) = \varphi_j^{\varepsilon,n}(x + y)$, for all $j \in \mathbb{Z}$ and $x \in \mathbb{R}$. Starting from definition (3.3), direct computations yield:

$$ e^{\nu^{\varepsilon,n+1}(x)/\varepsilon} = \left\langle e^{-\left(\varphi_j^{\varepsilon,n}(x+y-\Delta t v_j)\right)/\varepsilon} \right\rangle_{\Delta v} = e^{\mu^{\varepsilon,n+1}(x+y)/\varepsilon}, $$

from which (3.6) is deduced. Using (3.6), similar computations then lead to (3.7).

iii) The proof of (3.8)–(3.9) is analogous to the previous one and uses the algebraic properties of the exponentials.

<figure><img src="image_6.png" /></figure>

From Lemma 3.2 one can then deduce the exact propagation of Lipschitz bounds in position by the semi-discrete numerical scheme $(S_{\Delta t, \Delta v}^{\varepsilon})$ (see [18]), and:

**Corollary 3.3.** For all $n \in [1, N_t]$ and $\varepsilon > 0$, the solution to $(S_{\Delta t, \Delta v}^{\varepsilon})$ satisfies Assumption 1.
Another key aspect of the analysis is the identification of the minimizer $\argmin_{j \in \mathbb{Z}} \{\varphi_j^{\epsilon,n}(x - \Delta t v_j)\}$. Since $\varphi_{\text{in}}(x, v)$ grows like $v^2/2$ we expect that, its minima belong to a bounded interval. In other term, the minima cannot be achieved for a too large velocity. This is the content of the following lemma.

**Lemma 3.4.** Let $\Delta t, \Delta x, \Delta v > 0$ be fixed. Assume that the initial data $\varphi_{\text{in}}$ satisfies Assumption 1. Then, there exists a bounded set $\bar{V} \subset \mathbb{Z}$ such that for all $n \in [1, N_t]$, $\epsilon > 0$ and $x \in \mathbb{R}$,

$$ \argmin_{j \in \mathbb{Z}} \{\varphi_j^{\epsilon,n}(x - \Delta t v_j)\} \subset \Delta v \bar{V}. \qquad (3.10) $$

*Proof.* We proceed in two steps. We first localize $\argmin_{j \in \mathbb{Z}} \{\varphi_j^{\epsilon,n}(x)\}$ and then consider the effect of the transport. Let $(n, j) \in [1, N_t] \times \mathbb{Z}$, $\epsilon > 0$ and $x \in \mathbb{R}$. Thanks to Corollary 3.3,

$$ -M + \frac{v_j^2}{2} \le \varphi_j^{\epsilon,n}(x) \le M + \frac{v_j^2}{2}, \qquad (3.11) $$

which, for $j = 0$, gives in particular $|\varphi_0^{\epsilon,n}(x)| \le M$ and thus

$$ \min_{j \in \mathbb{Z}} \{\varphi_j^{\epsilon,n}(x)\} \le M, \qquad (3.12) $$

where we recall that $M$ depends neither on $n$, $\epsilon$ nor $x$. However, for $\Delta v > 0$ fixed there exists $J_1 \in \mathbb{N}$, such that for all $|k| > J_1$

$$ -M + \frac{v_k^2}{2} \ge M + 1. \qquad (3.13) $$

Combining (3.12) and (3.13) then gives that for all $|k| > J_1$

$$ \min_{j \in \mathbb{Z}} \{\varphi_j^{\epsilon,n}(x)\} + 1 < \varphi_k^{\epsilon,n}(x), $$

which means that the minimum is achieved for some $|j| < J_1$. The next step consists in dealing with the effect of the transport on the localization of the minima. Thanks to Corollary 3.3, $x \mapsto \varphi_j^{\epsilon,n}(x) - v^2/2$ is $L$-Lipschitz. Combined with (3.11), it yields

$$ -M + \frac{v_j^2}{2} - L\Delta t|v_j| \le \varphi_j^{\epsilon,n}(x - \Delta t v_j) \le L\Delta t|v_j| + M + \frac{v_j^2}{2}. \qquad (3.14) $$

Starting from (3.14), the same arguments as before can be used. The upper bound is in particular true at $j = 0$ giving $\min_{j \in \mathbb{Z}} \{\varphi_j^{\epsilon,n}(x - \Delta t v_j)\} \le M$. Considering the lower bound, there exists $J_2 \in \mathbb{N}$ such that for all $|k| > J_2$,

$$ -M + \frac{v_j^2}{2} - L\Delta t|v_j| > M + 1. $$

Concluding in same way as the transport-free minimization problem, we obtain (3.10) by setting $\bar{V} = [-\max\{J_1, J_2\}, \max\{J_1, J_2\}]$. $\square$

**Remark 1.** The extension of the stability results to the fully discrete scheme $(S^\epsilon)$ are actually straightforward given our choice of interpolation. Regarding Lemma 3.1, since the property is independent of the position variable, it remains unaffected by the spatial discretization. The maximum principle and the monotonicity of the scheme follow directly from the monotonicity of the linear interpolation (2.13). Finally, the commutative properties of $(S^\epsilon)$ with respect to constants and translations also carry over, as these properties are preserved by the linear interpolation. In the case of translations, it should be noted that they must correspond to shifts by grid points, i.e., of the form $y = i\Delta x$ for some $i \in \mathbb{Z}$. This is a natural assumption, given that the unknowns are localized on the grid.

**3.2.** **Limit $\epsilon \to 0$ of $(S^\epsilon)$**. This section is dedicated to the proof of the asymptotic preserving property of $(S^\epsilon)$ stated in Theorem 2.1. It is achieved in two steps. First, the logarithmic terms in $(S^\epsilon_{\Delta t, \Delta v})$ are shown to converge towards 0. Secondly, relying on the stability properties of Lemmas 3.1 and 3.2 and on Corollary 3.3, the limit is proven.

*Proof of Theorem 2.1.* Let us first point out that the numerical cost of solving $(S^\epsilon)$ is independent of $\epsilon$. The resolution is an explicit procedure thanks to (2.6) and the discretization parameters can be chosen independently of $\epsilon$ ensuring that $(S^\epsilon)$ is an efficient numerical scheme in practice.

**Step 1.** Let $x \in \mathbb{R}$, let us fix $\Delta t, \Delta v > 0$, $j \in \mathbb{Z}$ and $n \in [1, N_t]$. We recall that the quantities $M^{\epsilon,n}$ and $m^{\epsilon,n}$ have been introduced in (2.9) and (2.10) to ensure non-positivity in the exponential
terms of $(S^\epsilon_{\Delta t, \Delta v})$. Then, focusing on the first line of $(S^\epsilon_{\Delta t, \Delta v})$, we observe that, by definition (2.9), at least one exponential in the discrete integral is equal to one. More precisely, there exists $j_0 \in \mathbb{Z}$ such that $\varphi_{j_0}^{\epsilon,n}(x - \Delta t v_{j_0}) - m^{\epsilon,n}(x) = 0$. The following relation then holds:

$$ \left\langle e^{-[\varphi_j^{\epsilon,n}(x-\Delta t v_j) - m^{\epsilon,n}(x)]/\epsilon} \right\rangle_{\Delta v} = \Delta v + \sum_{j \in \mathbb{Z}, j \neq j_0} e^{-[\varphi_j^{\epsilon,n}(x-\Delta t v_j) - m^{\epsilon,n}(x)]/\epsilon} \Delta v, \quad (3.15) $$

where the remaining terms in the exponentials are non-positive. On the one hand, one obtains from (3.15)

$$ \epsilon \ln \left( \left\langle e^{-[\varphi_j^{\epsilon,n}(x-\Delta t v_j) - m^{\epsilon,n}(x)]/\epsilon} \right\rangle_{\Delta v} \right) \geq \epsilon \ln(\Delta v) \underset{\epsilon \to 0}{\to} 0. \quad (3.16) $$

On the other hand, we establish an upper bound for (3.15) by exploiting the quadratic growth of $\varphi_j^{\epsilon,n}$. From Corollary 3.3, we obtain in particular

$$ -M - m^{\epsilon,n}(x) + \frac{v_j^2}{2} - |v_j| \le \varphi_j^{\epsilon,n}(x - \Delta t v_j) - m^{\epsilon,n}(x) - |v_j| \le M - m^{\epsilon,n}(x) + \frac{v_j^2}{2} - |v_j|, \quad (3.17) $$

where the second inequality and (2.9) and with $j=0$ yields

$$ m^{\epsilon,n}(x) \le M. \quad (3.18) $$

Combining (3.17) and (3.18) gives

$$ -2M + \frac{v_j^2}{2} - |v_j| \le \varphi_j^{\epsilon,n}(x - \Delta t v_j) - m^{\epsilon,n}(x) - |v_j|. \quad (3.19) $$

Remark now that the set

$$ A := \{j \in \mathbb{Z} : -2M + \frac{v_j^2}{2} - |v_j| < 0\} \quad (3.20) $$

is finite, independently of $\epsilon$. Moreover, since $0 < \epsilon \le 1$, for any $j \in \bar{A}$, the complement set of $A$, one has

$$ -\frac{1}{\epsilon} (\varphi_j^{\epsilon,n}(x - \Delta t v_j) - m^{\epsilon,n}(x)) \le -|v_j|. \quad (3.21) $$

Thus, splitting the sum in (3.15) according to $A$ and $\bar{A}$ gives

$$ \left\langle e^{-[\varphi_j^{\epsilon,n}(x-\Delta t v_j) - m^{\epsilon,n}(x)]/\epsilon} \right\rangle_{\Delta v} = \Delta v \sum_{j \in A} e^{-[\varphi_j^{\epsilon,n}(x-\Delta t v_j) - m^{\epsilon,n}(x)]/\epsilon} + \Delta v \sum_{j \in \bar{A}} e^{-[\varphi_j^{\epsilon,n}(x-\Delta t v_j) - m^{\epsilon,n}(x)]/\epsilon}. $$

The first sum is finite and bounded independently of $\epsilon$ by $\Delta v \# A$. Thanks to (3.21), the second sum is bounded from above by a convergent series that is also independent of $\epsilon$. Consequently, the quantity in (3.15) is bounded by a constant $C > 0$ for all $\epsilon$, and

$$ \epsilon \ln \left( \left\langle e^{-[\varphi_j^{\epsilon,n}(x-\Delta t v_j) - m^{\epsilon,n}(x)]/\epsilon} \right\rangle_{\Delta v} \right) \le \epsilon \ln(C) \underset{\epsilon \to 0}{\to} 0. \quad (3.22) $$

Finally, combining (3.22) with (3.16), we conclude that

$$ \forall j \in \mathbb{Z}, \forall \Delta t, \Delta v > 0, \quad \epsilon \ln \left( \left\langle e^{-[\varphi_j^{\epsilon,n}(x-\Delta t v_j) - m^{\epsilon,n}(x)]/\epsilon} \right\rangle_{\Delta v} \right) \underset{\epsilon \to 0}{\to} 0. \quad (3.23) $$

Using a similar argument, we deal with the logarithmic term in the second line of $(S^\epsilon_{\Delta t, \Delta v})$. By definition (2.10), at least one of the two exponentials is necessarily equal to one. Distinguishing two cases depending on the value of (2.10), the argument of the logarithm reduces to

$$ \begin{cases} c_{\Delta v}^{\epsilon} + (1 - e^{-\Delta t/\epsilon}) e^{-\left[\frac{v_j^2}{2} + \mu^{\epsilon,n+1}(x) - M_j^{\epsilon,n}(x)\right]/\epsilon}, & \text{if } \varphi_j^{\epsilon,n}(x - \Delta t v_j) + \Delta t \le \frac{v_j^2}{2} + \mu^{\epsilon,n+1}(x), \\ c_{\Delta v}^{\epsilon} e^{-[\varphi_j^{\epsilon,n}(x - \Delta t v_j) + \Delta t - M_j^{\epsilon,n}(x)]/\epsilon} + (1 - e^{-\Delta t/\epsilon}), & \text{else.} \end{cases} \quad (3.24) $$

In both cases, the remaining exponential terms have negative arguments and in the limit $\epsilon \to 0$, the dominating term is at most $c_{\Delta v}^{\epsilon} + 1$. Then, recalling that $\Delta v$ is fixed, we remark that $\epsilon \ln(c_{\Delta v}^{\epsilon})$ converges towards 0 thanks to the relation $\epsilon \ln(c_{\Delta v}^{\epsilon}) = \epsilon \ln(\sqrt{2\pi\epsilon} + O(\Delta v))$, where we used the error estimate on first order quadratures (2.12). Combining (3.24) with this observation, we obtain, for
all $j \in \mathbb{Z}$ and $\Delta t, \Delta v > 0$, that the logarithmic terms in the second line of $(S^\epsilon_{\Delta t, \Delta v})$ also vanishes in the limit $\epsilon \to 0$:

$$ \epsilon \ln \left( c_{\Delta v}^{\epsilon} e^{-[\varphi_{j}^{\epsilon,n}(x-\Delta tv_{j})+\Delta t-M_{j}^{\epsilon,n}(x)]/\epsilon} \right. \\ \left. + \left(1-e^{-\Delta t/\epsilon}\right) e^{-\left[\frac{v_{j}^{2}}{2}+\mu^{\epsilon,n+1}(x)-M_{j}^{\epsilon,n}(x)\right]/\epsilon} \right) + \epsilon \ln(c_{\Delta v}^{\epsilon}) \underset{\epsilon \to 0}{\to} 0. \quad (3.25) $$

**Step 2.** Let $\Delta t, \Delta v$, and $j \in \mathbb{Z}$ be fixed. The next step in the proof is to pass to the limit in the terms $m^{\epsilon,n}(x)$ and $M_j^{\epsilon,n}(x)$. We begin with $m^{\epsilon,n}(x)$. Let $I \subset \mathbb{R}$ an interval and let $x \in I$.

First, by applying Corollary 3.3, which provides a Lipschitz bound propagation, one may apply Ascoli's theorem on the family of functions $(y \mapsto \varphi_j^{\epsilon,n}(y) - v_j^2/2)_{\epsilon>0}$ for any $j \in \mathbb{Z}$. Therefore, for all $j \in \mathbb{Z}$, this yields convergence, up to the extraction of a subsequence, uniformly on $I$ towards a function $y \mapsto \varphi_j^n(y) - v_j^2/2$.

As a consequence, for all $j \in \mathbb{Z}$, the limit function $\varphi_j^n$ satisfies Assumption 1. Using the same reasoning as in the proof of Lemma 3.4, the minimum $\min_{j \in \mathbb{Z}} \{\varphi_j^n(x - \Delta tv_j)\}$ is well-defined and attained within a bounded set $[-J, J]$ for some $J \in \mathbb{N}$.

Secondly, let $A > 0$ and define $I = [-A - \Delta tv_{-J}, A + \Delta tv_{J}]$ and fix $j \in [-J, J]$. The family $(y \mapsto \varphi_j^{\epsilon,n}(y) - v_j^2/2)_{\epsilon>0, n \in \{1, N_t\}}$ converges uniformly on $I$ thanks to the above discussion. In particular, for $x \in [-A, A]$, one has $\varphi_j^{\epsilon,n}(x - \Delta tv_j) \to \varphi_j^n(x - \Delta tv_j)$, up to a subsequence. In what follows, the notation $\epsilon \to 0$ will refer to a subsequence for which these convergences are true for all $j \in [-J, J]$.

Since the minimum of Lipschitz functions is Lipschitz and the index set $[-J, J]$ is finite, the minimum is taken over a set of uniformly bounded and Lipschitz functions. Consequently, there exists a function $x \mapsto m^n(x)$ such that

$$ m^{\epsilon,n}(x) \underset{\epsilon \to 0}{\to} m^n(x), $$

uniformly on $[-A, A]$. It remains to identify the limit $m^n(x)$. As, for all $j \in [-J, J]$, the convergence of $(\varphi_j^{\epsilon,n}(\cdot - \Delta tv_j))_\epsilon$ towards $\varphi^n(\cdot - \Delta tv_j)$ is uniform on $[-A, A]$, for any $\eta > 0$, there exists $\epsilon_0 > 0$ such that for all $\epsilon < \epsilon_0$, $j \in [-J, J]$, $n \in \{1, N_t\}$, and $x \in [-A, A]$,

$$ |\varphi_j^{\epsilon,n}(x - \Delta tv_j) - \varphi_j^n(x - \Delta tv_j)| \leq \eta. $$

Taking the minimum over $j \in [-J, J]$ then gives

$$ \left| m^{\epsilon,n}(x) - \min_{j \in [-J, J]} \{\varphi_j^n(x - \Delta tv_j)\} \right| \leq \eta. $$

Therefore, $m^{\epsilon,n}(x) = \min_{j \in \Delta v \mathbb{Z}} \{\varphi_j^n(x - \Delta tv_j)\}$. This convergence, together with (3.23), implies

$$ \mu^{\epsilon,n+1}(x) \underset{\epsilon \to 0}{\to} \min_{j \in \Delta v \mathbb{Z}} \{\varphi_j^n(x - \Delta tv_j)\}, \quad (3.26) $$

locally uniformly towards $\mu^{n+1}(x)$.

We now turn our attention to the limit of $M_j^{\epsilon,n}(x)$. Using the same argument as before, uniform continuity of $(x, y) \mapsto \min\{x, y\}$ preserves the uniform convergence on $x \in [-A, A]$, and we can pass to the limit in (2.10) giving for all $x \in [-A, A]$, $j \in [-J, J]$ and $n \in \{1, N_t\}$:

$$ M_j^{\epsilon,n}(x) \underset{\epsilon \to 0}{\to} \min \left\{ \varphi_j^n(x - \Delta tv_j) + \Delta t, \frac{v_j^2}{2} + \mu^{n+1}(x) \right\}, \quad (3.27) $$

which is the definition of $\varphi_j^{n+1}(x)$. Finally, combining the convergence of (3.25), (3.26) and (3.27) we can conclude on the convergence of $(S^\epsilon_{\Delta t, \Delta v})$, as $\epsilon \to 0$, towards the semi-discrete limit scheme:

$$ \left\{ \begin{array}{l} \mu^{n+1}(x) = \min_{j \in \mathbb{Z}} \{\varphi_j^n(x - \Delta tv_j)\}, \\ \\ \varphi_j^{n+1}(x) = \min \left\{ \varphi_j^n(x - \Delta tv_j) + \Delta t, \frac{v_j^2}{2} + \mu^{n+1}(x) \right\}. \end{array} \right. \quad (S^0_{\Delta t, \Delta v}) $$

and this convergence is locally uniform in $x \in \mathbb{R}$. It remains to remove the restriction “up to a subsequence” in the convergence as $\epsilon \to 0$. Since $(S^0_{\Delta t, \Delta v})$ is an explicit scheme, and since, for all $j \in \mathbb{Z}$ and $x \in \mathbb{R}$, the initial data $\varphi_{in,j}(x)$ is assumed to be independent of $\epsilon$, this restriction
can be lifted. If, instead, the initial data does depend on $\varepsilon$ — that is, if one further assume that $\varphi_{\text{in},j}^{\varepsilon}(x) \to \varphi_{\text{in},j}(x)$ locally uniformly in $x$ as $\varepsilon \to 0$ — then all subsequences satisfy the same explicit induction relation. Consequently, the sequences $\varphi_j^{\varepsilon, n+1}$ and $\mu_j^{\varepsilon, n+1}$ converge locally uniformly when $\varepsilon \to 0$ and their limits satisfy $(S_{\Delta t, \Delta v}^0)$. $\square$

The proof of Theorem 2.1 can be extended to the fully discrete setting following Remark 1 by using the properties of the linear interpolation (2.13).

### 4. PROPERTIES OF THE ASYMPTOTIC SCHEME

In this section, we again focus on the semi discrete scheme $(S_{\Delta t, \Delta v}^0)$. The fully discrete scheme $(S^0)$ mainly adds technicalities. We start by showing stability properties satisfied by $(S_{\Delta t, \Delta v}^0)$ in Section 4.1. Then, we present several equivalent reformulations of $(S_{\Delta t, \Delta v}^0)$. In particular, following the continuous setting, a discrete representation formula is introduced in Section 4.3.

**4.1. Stability properties of the limit scheme.** In this section, we prove stability properties of the limit scheme $(S_{\Delta t, \Delta v}^0)$, such as the preservation of the equilibrium and maximum principle that are classical in the study of numerical schemes for Hamilton-Jacobi equations [18]. Note that the same procedures can be adapted to the fully discrete case thanks to the monotonicity property of the first order interpolation (2.13), and we refer back to Remark 1 for more details. The section is concluded with two Lemmas that emphasize the link between $(S_{\Delta t, \Delta v}^0)$ and $(\mathcal{P}^0)$.

**Lemma 4.1.** The semi-discrete scheme $(S_{\Delta t, \Delta v}^0)$ preserves the equilibrium $(v_j^2/2)_{j \in \mathbb{Z}}$; if $\varphi_{\text{in}}(x, v) = v^2/2$ for all $x \in \mathbb{R}$ and $v \in \mathbb{R}$, then

$$ \forall n \in [1, N_t], \forall j \in \mathbb{Z}, \forall x \in \mathbb{R}, \quad \varphi_j^n(x) = \frac{v_j^2}{2}. $$

**Lemma 4.2.** Let $(n, j) \in [1, N_t] \times \mathbb{Z}$ and $x \in \mathbb{R}$. Let $(\mu^n(x), \varphi_j^n(x))_{n \in [1, N_t], j \in \mathbb{Z}}$ be a solution to $(S_{\Delta t, \Delta v}^0)$ with initial data (2.11) and denote $(\psi^n(x))_{n \in [1, N_t], j \in \mathbb{Z}}$ another solution to $(S_{\Delta t, \Delta v}^0)$. For any function $\psi_j^n(x)$, we use the notation

$$ \nu^n(x) = \min_{j \in \mathbb{Z}} \{ \psi_j^{n-1}(x - \Delta t v_j) \}. $$

Then, the following properties hold:

i) **Monotonicity.** Suppose that $\varphi_{\text{in}} \le \psi_{\text{in}}$. Then, for all $j \in \mathbb{Z}$, $x \in \mathbb{R}$ and for all $n \in [1, N_t]$,
$$ \mu^n(x) \le \nu^n(x), $$
$$ \varphi_j^n(x) \le \psi_j^n(x). $$

ii) **Invariance by translation.** Let $y \in \mathbb{R}$ and suppose that $\psi_{\text{in}} = \varphi_{\text{in}}(\cdot + y, \cdot)$. Then, for all $j \in \mathbb{Z}$, $x \in \mathbb{R}$ and for all $n \in [1, N_t]$,
$$ \nu^n(x) = \mu^n(x + y), $$
$$ \psi_j^n(x) = \varphi_j^n(x + y). $$

iii) **Commutation with constants.** Let $K \in \mathbb{R}$ and suppose that $\psi_{\text{in}} = \varphi_{\text{in}} + K$. Then, for all $j \in \mathbb{Z}$, $x \in \mathbb{R}$ and for all $n \in [1, N_t]$,
$$ \nu^n(x) = \mu^n(x) + K, $$
$$ \psi_j^n(x) = \varphi_j^n(x) + K. $$

From Lemma 4.2, one can then deduce the exact propagation of Lipschitz bounds by the semi-discrete numerical scheme $(S_{\Delta t, \Delta v}^0)$ (see [18]). In addition, one has:

**Corollary 4.3.** For all $n \in [1, N_t]$, the solution to $(S_{\Delta t, \Delta v}^0)$ satisfies Assumption 1.

The proofs of Lemmas 4.2 and 4.1 are done by induction, using the definition of $(S_{\Delta t, \Delta v}^0)$ and are therefore not detailed here.

We conclude this section with two lemmas linking $(S_{\Delta t, \Delta v}^0)$ and $(\mathcal{P}^0)$. First, we prove that $(S_{\Delta t, \Delta v}^0)$ also satisfies a discrete version of the condition $\partial_t \min_{v \in \mathbb{R}^{d_v}} \{\varphi\} \le 0$ from $(\mathcal{P}^0)$. To show this result, we start by introducing the following property of $(S_{\Delta t, \Delta v}^0)$.
**Lemma 4.4.** Let $\Delta t, \Delta x, \Delta v > 0$. Let $(\varphi_j^n(x))_{n \in [1, N_t], j \in \mathbb{Z}}$ and $(\mu^n(x))_{n \in [1, N_t]}$ be solution to $(S^0_{\Delta t, \Delta v})$. Then,

$$ \forall n \in [1, N_t], \quad \mu^n(x) = \min_{j \in \mathbb{Z}} \{\varphi_j^n(x)\}. \qquad (4.1) $$

*Proof.* Starting from the second line of $(S^0_{\Delta t, \Delta v})$, the minimum in $j$ is taken, and since $\Delta t \ge 0$, one obtains

$$ \begin{aligned} \min_{j \in \mathbb{Z}} \{\varphi_j^{n+1}(x)\} &= \min \left\{ \min_{j \in \mathbb{Z}} \{\varphi_j^n(x - \Delta t v_j)\} + \Delta t, \min_{j \in \mathbb{Z}} \left\{ \frac{v_j^2}{2} \right\} + \mu^{n+1}(x) \right\} \\ &= \mu^{n+1}(x), \end{aligned} $$

□

Using Lemma 4.4, the following property, can now be stated.

**Lemma 4.5.** Let $x \in \mathbb{R}$ and $(\mu^n(x), \varphi_j^n(x))_{j \in \mathbb{Z}, n \in [0, N_t]}$ be a solution to $(S^0_{\Delta t, \Delta v})$. Then, under condition (2.19),

$$ \forall x \in \mathbb{R}, n \in [0, N_t - 1], \quad \mu^{n+1}(x) - \mu^n(x) \le 0. \qquad (4.2) $$

*Proof.* Let $x \in \mathbb{R}$. Starting from (4.1), it gives in particular that,

$$ \forall j \in \mathbb{Z} \quad \mu^{n+1}(x) \le \varphi_j^{n+1}(x). $$

Then, using the second line of $(S^0_{\Delta t, \Delta v})$, we deduce that for all $j \in \mathbb{Z}$,

$$ \varphi_j^{n+1}(x) \le \mu^{n+1}(x) + \frac{v_j^2}{2}. \qquad (4.3) $$

In particular, combining the first line of $(S^0_{\Delta t, \Delta v})$ and (4.3) at time $n$, one obtains, for $j = 0$

$$ \varphi_0^{n+1}(x) \le \mu^{n+1}(x) \le \varphi_0^n(x) \le \mu^n(x), $$

which yields (4.2).

□

To conclude this section, we present an equivalent reformulation of $(S^0_{\Delta t, \Delta v})$ that more closely mirrors the structure of the continuous equation $(\mathcal{P}^0)$.

**Lemma 4.6.** Let $(\mu^n(x), \varphi_j^n(x))_{j \in \mathbb{Z}, x \in \mathbb{R}, n \in [0, N_t]}$ be solution to $(S^0_{\Delta t, \Delta v})$. Then, it is solution to:

$$ \left\{ \begin{array}{l} \mu^{n+1}(x) = \min_{j \in \mathbb{Z}} \{\varphi_j^n(x - \Delta t v_j)\}, \\ \\ \max \left\{ \frac{\varphi_j^{n+1}(x) - \varphi_j^n(x - \Delta t v_j)}{\Delta t} + 1, \varphi_j^{n+1}(x) - \frac{v_j^2}{2} - \mu^{n+1}(x) \right\} = 0. \end{array} \right. $$

4.2. **Towards a discrete representation formula.** To derive a semi-discrete representation formula for $(\mathcal{P}^0)$, we begin by introducing the necessary notations and definitions. These are based on the discussion in the introduction concerning the duality between Eulerian and Lagrangian perspectives. We first recall the notion of minimizing curves for $(\mathcal{P}^0)$ defined in [12], along with their associated continuous action. We then introduce corresponding concepts in the semi-discrete setting and recall the representation formula for $(\mathcal{P}^0)$.

**Definition 4.7** (Continuous minimizing curves). Let $\Sigma^{0,T}$ be the space of piecewise constant, càdlàg functions defined over the time interval $(0, T]$ taking values in $\mathbb{R}$. For any $\sigma \in \Sigma^{0,T}$, denote $s^0 = 0$ and let $\mathbf{s} = (s^k)_{0 \le k \le N_*}$ be the $N_* \ge 1$ times, of discontinuity of $\sigma$ in $(0, T]$, such that

$$ \sigma = w^0 \mathbf{1}_{(0, s_1)} + \sum_{k=1}^{N_*-1} w^k \mathbf{1}_{[s^k, s^{k+1})} + w^{N_*} \mathbf{1}_{[s^{N_*}, T]}. $$

Then, for any starting point $y \in \mathbb{R}$, the piecewise linear curve $\gamma \in \mathbb{R} \times \Sigma^{0,T}$ is defined for all $\tau \in [0, T]$ as

$$ \gamma(\tau) = y + \int_{0}^{\tau} \sigma(\tau') d\tau'. $$

In addition, we denote by $\mathbf{w} = (w^0, \dots, w^{N_*})$ the associated tuple of velocity jumps. It allows to characterize the path $\gamma$ by its starting point and the sequences $\mathbf{s}$ and $\mathbf{w} : \gamma = (y, \mathbf{s}, \mathbf{w}) \in \mathbb{R} \times \Sigma^{0,T}$.
From this definition, the action $\bar{A}_0^T[\gamma]$ of a path $\gamma$ is defined as follows.

**Definition 4.8** (Action of a piecewise linear curve (continuous case)). Let $\gamma = (y, \mathbf{s}, \mathbf{w}) \in \mathbb{R} \times \Sigma^{0,T}$ be defined as in Definition 4.7. The continuous action of $\gamma$ on $(0, T]$ is defined as follows:

$$ \bar{A}_0^T[\gamma] = \frac{1}{2} \sum_{k=1}^{N_*} |w^k|^2 + \sum_{k=0}^{N_*} (s^{k+1} - s^k) \mathbf{1}_{w^k \neq 0}, \quad (4.4) $$

with the convention $s^0 = 0$ and $s^{N_*+1} = T$.

The action $\bar{A}$ describes the cost of a moving individual. The first part corresponds to a change in velocity that is instantaneous and induces a cost $|w|^2/2$, where $w$ is the new velocity. In addition, when the individual is subject to free transport, there is a running cost that depends only on the time spent traveling and not on the velocity at which it occurs.

With these definitions, we recall the following result from [12, Theorem 1.8], which introduces a representation formula and establishes the uniqueness of the viscosity solution to $(\mathcal{P}^0)$.

**Theorem 4.9** ((1.8 in [12]), Kinetic Hopf-Lax formula). Let $\varphi_{\text{in}}$ be a continuous function such that $\varphi_{\text{in}} - v^2/2 \in L^{\infty}(\mathbb{R} \times \mathbb{R})$. Then, the following representation formula

$$ \varphi(t, x, v) = \inf_{\substack{\bar{\gamma}=(y, \mathbf{s}, \mathbf{w}) \in \mathbb{R} \times \Sigma^{0,T} \\ \bar{\gamma}(t)=x, \bar{\gamma}(t)=v}} \{\varphi_{\text{in}}(y, w^0) + \bar{A}_0^t[\bar{\gamma}]\} \quad (4.5) $$

is the viscosity solution of $(\mathcal{P}^0)$ with initial data $\min_{v \in \mathbb{R}} \{\varphi_{\text{in}}, \min_{v \in \mathbb{R}} \{\varphi_{\text{in}}\} + \frac{v^2}{2}\}$.

In the discrete setting, analogous definitions can be introduced. A semi-discrete analog of (4.5) will be derived in the next section. We begin with time-discrete curves, where the velocity jump times are restricted to the grid $(n\Delta t)_{n \in \llbracket 0, N_t \rrbracket}$. The following definition can be viewed as a discrete counterpart — or a special case — of Definition 4.7. It serves as a technical bridge between fully continuous minimizing curves and their time-velocity semi-discrete counterparts.

**Definition 4.10** (Time-discrete minimizing curves). Let $\Sigma_{\Delta t}^{0,T}$ be the subset of $\Sigma^{0,T}$ such that for any $\sigma \in \Sigma_{\Delta t}^{0,T}$, the $N_*$ times $\mathbf{s} = (s^k)_{1 \le k \le N_*}$ of discontinuity of $\sigma \in (0, T]$ are located on the time grid: for all $k \in \llbracket 1, N_* \rrbracket$, $s^k \in \Delta t\mathbb{N}$. Then, for any starting point $y \in \mathbb{R}$ and for any tuple $\mathbf{w} = (w^0, \dots, w^{N_*})$ of velocity jumps, the piecewise linear curve $\gamma = (y, \mathbf{s}, \mathbf{w}) \in \mathbb{R} \times \Sigma_{\Delta t}^{0,T}$ is defined on the time grid $(t^n)_{0 \le n \le N_t}$ as in Definition 4.7.

As the discontinuities of the elements of $\Sigma_{\Delta t}^{0,T}$ are located on the grid, the element $\sigma \in \Sigma_{\Delta t}^{0,T}$ can also be described with all the values it takes on time steps $[k\Delta t, (k+1)\Delta t)$, for all $k \in \llbracket 0, N_t - 1 \rrbracket$, and at time $T = t^{N_t}$. Introducing $\mathbf{w} = (w^0, \dots, w^{N_t}) \in \mathbb{R}^{[0, N_t]}$, $\sigma$ is defined as

$$ \sigma = w^0 \mathbf{1}_{(0, \Delta t)} + \sum_{k=1}^{N_t-1} w^k \mathbf{1}_{[k \Delta t, (k+1) \Delta t)} + w^{N_t} \mathbf{1}_{\{T\}}, $$

and for any initial point $y \in \mathbb{R}$, the associated piecewise linear curves $\gamma = (y, \mathbf{w}) \in \mathbb{R} \times \mathbb{R}^{[1, N_t]}$ is defined on the time grid for all $n \in \llbracket 0, N_t \rrbracket$ by

$$ \gamma(t^n) = y + \Delta t \sum_{k=1}^{n} w^{k-1}. $$

In the definition above, $\mathbf{w}$ may contain repetitions of velocities on adjacent intervals, while $\mathbf{w}$ does not. The first convention will be used when the focus is on the times of discontinuities of the paths, whereas the second will be used when the focus will be on the sequence of velocities arising in the path. In what follows, the statement $\gamma = (y, \mathbf{s}, \mathbf{w}) \in \mathbb{R} \times \Sigma_{\Delta t}^{0,T}$ refers to the first convention, while $\gamma = (y, \mathbf{w}) \in \mathbb{R} \times \mathbb{R}^{[0, N_t]}$ refers to the second.

Further discretizing the velocity variable on the grid $\Delta v\mathbb{Z}$ leads to a semi-discrete in time and velocity setting. In this framework, one can define a specific class of discrete paths that form another special case of Definition 4.7. These paths are particularly suited to build a discrete analog of the representation formula (4.5).

**Definition 4.11** (Time-velocity discrete minimizing curves). Let $\Sigma_{\Delta t, \Delta v}^{0,T}$ be the subset of $\Sigma_{\Delta t}^{0,T}$ such that for any $\sigma \in \Sigma_{\Delta t, \Delta v}^{0,T}$, the $N_*+1$ velocities $\mathbf{w} = (w^k)_{0 \le k \le N_*}$ are located on the velocity
grid : for all $k \in [0, N_*)$, $w^k \in \Delta v \mathbb{Z}$. Then for any starting point $y \in \mathbb{R}$, the piecewise linear curve $\gamma = (y, s, w) \in \mathbb{R} \times \Sigma_{\Delta t, \Delta v}^{0,T}$ is defined on the time grid as in Definition 4.7.

As in Definition 4.10, the path $\gamma$ can also be defined by its initial point $y \in \mathbb{R}$ and the sequence of velocities $w = (w^k)_{0 \le k \le N_t} \in \Delta v \mathbb{Z}^{[0, N_t]}$. We will denote $\gamma = (y, w) \in \mathbb{R} \times \Delta v \mathbb{Z}^{[0, N_t]}$ when we use this convention. With both conventions, the velocities $w^k$ arising in the time-velocitie discrete minimizing curves of Definition 4.11 lie on the velocity grid. Therefore, for all $k$, we introduce the index $j(k) \in \mathbb{Z}$, such that $w^k = v_{j(k)}$.

Since Definitions 4.10 and 4.11 are particular cases of Definition 4.7, the continuous action $\bar{A}$ can still be evaluated along discrete paths. However, the goal of the next section is to construct a fully discrete counterpart to the representation formula (4.5), which is equivalent to the semi-discrete system $(S_{\Delta t, \Delta v}^0)$. To this end, we introduce a discrete analog of Definition 4.8, tailored specifically to time-discrete paths and adapted to the structure of $(S_{\Delta t, \Delta v}^0)$.

**Definition 4.12** (Action of a piecewise linear curve (discrete case)). Let $\gamma = (y, w) \in \mathbb{R} \times \mathbb{R}^{[0, N_t]}$ be defined as in Definition 4.10. The discrete action of $\gamma$ on $(0, T]$ is defined as follows:

$$ \tilde{A}_0^T[\gamma] = \sum_{k=0}^{N_t-1} A_{w^k}^{w^{k+1}}, $$

where for all $(v, w) \in \mathbb{R}^2$, $A_w^v$ denotes the cost of a jump from velocity $w$ to $v$, and is defined as

$$ A_w^v = \begin{cases} \Delta t \mathbb{1}_{v \neq 0} & \text{if } v = w, \\ \frac{|v|^2}{2} & \text{if } v \neq w. \end{cases} \qquad (4.6) $$

This discrete action represents the cost of transitioning from velocity $w$ to velocity $v$ on a discrete time interval of size $\Delta t$. The first case indicates that if there is no change in velocity, the cost is simply the time spent in free transport, except when $v = 0$. The second case shows that the cost of changing velocity depends only on the new velocity and not on the initial one.

While the continuous action (4.4) and the discrete one (4.6) share similarities, such as the cost of a velocity jump and free transport, there is a key difference between the two descriptions. In the continuous setting, a jump occurs instantaneously, whereas in the discrete case, it takes an entire time interval $\Delta t$ to get to a new velocity. This discrepancy introduces an error of order $\Delta t$ per jump (if it is not at the end of the interval) when comparing the cost of the same path under the continuous and discrete actions. This is the content of the following lemma.

**Lemma 4.13.** Let $\gamma = (y, s, w) \in \mathbb{R} \times \Sigma_{\Delta t}^{0,T}$, and let $N_*$ the number of velocity jumps in $\gamma$. Therefore, with the convention $s^0 = 0$, $s = (s^0, \dots, s^{N_*})$ and $w = (w^0, \dots, w^{N_*})$ are such that

$$ \forall k \in [1, N_*], \quad s^k - s^{k-1} \ge \Delta t, \text{ and } w^k \ne w^{k-1}. $$

The action $\tilde{A}[\gamma]$ of Definition 4.12 can be rewritten

$$ \tilde{A}[\gamma] = \sum_{k=1}^{N_*} \frac{|w^k|^2}{2} + \sum_{k=0}^{N_*-2} (s^{k+1} - s^k - \Delta t) \mathbb{1}_{w^{k+1} \neq 0} + \max \{T - s^{N_*} - \Delta t, 0\} \mathbb{1}_{w^{N_*} \neq 0}. \quad (4.7) $$

Moreover,

$$ |\tilde{A}[\gamma] - \tilde{A}[\gamma]| \le N_* \Delta t. $$

*Proof.* Equation (4.7) is a consequence of the definition of $\tilde{A}[\gamma]$, obtained by rewriting the velocities $w$ of $\gamma$ as a sequence $w \in \mathbb{R}^{[0, N_t]}$ defined on the time grid. The last term is taken apart to distinguish the cases $s^{N_*} = T$ and $s^{N_*} \ne T$. The difference between $\tilde{A}$ and $\tilde{A}$ follows immediately. $\square$

4.3. **Derivation of the semi-discrete representation formula.** We now introduce a semi-discrete representation formula for $(S_{\Delta t, \Delta v}^0)$. In line with the continuous setting, this formulation involves the discrete action (4.6) along with the initial datum evaluated at the foot of a specific characteristic curve. We have the following proposition:

**Proposition 4.14** (Semi-discrete representation formula). Assume that the time step $\Delta t$ and the velocity step $\Delta v$ satisfy the relation (2.19), that is $\Delta t \le \Delta v^2/2$, and assume moreover that the
initial data satisfies Assumption 1. Then, the solution $(\varphi_j^{N_t}(x))_{j \in \mathbb{Z}}$ to $(S_{\Delta t, \Delta v}^0)$, whose existence is known, is also solution to the following discrete representation formula:

$$ \varphi_j^{N_t}(x) = \inf_{\substack{\gamma=(y,w) \in \mathbb{R} \times \Delta v \mathbb{Z}[0,N_t] \\ \gamma(T)=x, w^{N_t}=v_j}} \{ \varphi_{\text{in},j(0)}(y) + \tilde{A}_0^T[\gamma] \}, \quad (4.8) $$

where $\gamma$ is defined in Definition 4.11, the discrete action $\tilde{A}_0^T$ is given in Definition 4.8, and we use the notation $v_j(0) = w^0$.

The assumption on the time step can be understood as the requirement that, at the discrete level, free transport and a velocity jump cannot occur simultaneously. This condition is purely discrete and specific to both the model and this reformulation. It highlights a key distinction from the continuous setting, where jumps are instantaneous and independent of transport.

Moreover, Assumption 1 on the initial data, together with the definition (4.6) of the discrete action, ensure coercivity of the problem. In particular, the infimum in (4.8) is well defined.

The proof of Proposition 4.14 proceeds by induction, using the following lemma.

**Lemma 4.15.** Assume that the time step $\Delta t$ and the velocity step $\Delta v$ satisfy the relation (2.19). Assume moreover that the initial data satisfies Assumption 1. Then, the solution $(\varphi_j^{n+1}(x))_{j \in \mathbb{Z}}$ to $(S_{\Delta t, \Delta v}^0)$ is solution to the following discrete representation formula:

$$ \varphi_j^{n+1}(x) = \inf_{w^n \in \Delta v \mathbb{Z}} \{ \varphi_j^n(x - \Delta t w^n) + A_{w^n}^{v_j} \}, \quad (4.9) $$

where $A_{w^n}^{v_j}$ is defined in (4.6) and $j(n)$ in Definition 4.11.

*Proof.* The proof consists in showing that under condition (2.19), $\varphi_j^{n+1}(x)$ defined in (4.9) coincides with $\varphi_j^{n+1}(x)$ defined with $(S_{\Delta t, \Delta v}^0)$. Start by considering $(U_j^n(x))_{j \in \mathbb{Z}, x \in \mathbb{R}}$ satisfying

$$ U_j^{n+1}(x) = \inf_{w^n \in \Delta v \mathbb{Z}} \{ U_j^n(x - \Delta t w^n) + A_{w^n}^{v_j} \}. \quad (4.10) $$

Let us now assume that, for all $j \in \mathbb{Z}$ and $x \in \mathbb{R}$, we have $U_j^n(x) = \varphi_j^n(x)$, and that these functions satisfy Assumption 1. Our goal is to show that $U_j^{n+1}(x) = \varphi_j^{n+1}(x)$. Let $x \in \mathbb{R}$ be fixed. By separating the cases $v_j = w^n$ ($j(n) = j$) and $v_j \neq w^n$ ($j(n) \neq j$), equation (4.10) can be reformulated as follows:

$$ U_j^{n+1}(x) = \min \left\{ \varphi_j^n(x - \Delta t v_j) + \Delta t \mathbf{1}_{\{v_j \neq 0\}}, \inf_{w^n \neq v_j} \left\{ \varphi_{j(n)}^n(x - \Delta t w^n) + \frac{v_j^2}{2} \right\} \right\}. \quad (4.11) $$

We now consider separately the cases $v_j = 0$ and $v_j \neq 0$.

• Case $v_j \neq 0$. In that case, one has:

$$ U_j^{n+1}(x) = \min \left\{ \varphi_j^n(x - \Delta t v_j) + \Delta t, \inf_{w^n \neq v_j} \left\{ \varphi_{j(n)}^n(x - \Delta t w^n) + \frac{v_j^2}{2} \right\} \right\}. \quad (4.12) $$

Using condition (2.19) and since $v_j \neq 0$, the following inequality holds

$$ \Delta t \leq \frac{v_j^2}{2}. \quad (4.13) $$

Moreover, the case $w^n = v_j$ can be included in the infimum in (4.12), as it does not attain the minimum: the first term in the expression is the minimizer, as ensured by (4.13). This yields:

$$ U_j^{n+1}(x) = \min \left\{ \varphi_j^n(x - \Delta t v_j) + \Delta t, \inf_{w^n \in \Delta v \mathbb{Z}} \left\{ \varphi_{j(n)}^n(x - \Delta t w^n) + \frac{v_j^2}{2} \right\} + \frac{v_j^2}{2} \right\}. \quad (4.14) $$

Since $\varphi^n$ satisfies Assumption 1, the infimum in (4.14) is in fact a minimum. Consequently, as $v_j \neq 0$, we observe that both $\varphi_j^{n+1}(x)$ and $U_j^{n+1}(x)$ satisfy the same induction relation. It follows that $\varphi_j^{n+1}(x) = U_j^{n+1}(x)$ for all $j \in \mathbb{Z}^*$ and $x \in \mathbb{R}$.
* Case $v_j = 0$ ($j = 0$). On the one hand, using the notations from Definition 4.11 for the discrete velocities, the first line of $(S_{\Delta t, \Delta v}^0)$ rewrites

$$ \mu^{n+1}(x) = \min_{w^n \in \Delta v \mathbb{Z}} \{\varphi_{j(n)}^n(x - \Delta t w^n)\}. \quad (4.15) $$

On the other hand, evaluating (4.11) at $j = 0$ gives

$$ U_0^{n+1}(x) = \min \left\{ U_0^n(x), \inf_{w^n \neq 0} \left\{ U_{j(n)}^n(x - \Delta t w^n) \right\} \right\}. \quad (4.16) $$

As $U^n$ satisfies Assumption 1, the infimum is also attained meaning it is actually a minimum. Then, comparing (4.15) and (4.16), one concludes that $U_0^{n+1}(x) = \varphi_0^{n+1}(x)$ for all $x \in \mathbb{R}$.

□

*Proof of Proposition 4.14.* Thanks to Lemma 4.15, the following relation is obtained by induction:

$$ \varphi_j^{N_t}(x) = \inf_{w \in \Delta v \mathbb{Z}^{[0,N_t]}} \left\{ \varphi_{j(0)}^0 \left( x - \Delta t \sum_{k=0}^{N_t-1} w^k \right) + \sum_{k=0}^{N_t-1} A_{w^k}^{w^{k+1}} \right\}, $$

with the notation $w^{N_t} = v_j$. The proof is completed by introducing the piecewise linear curve $\gamma$ as in Definition 4.11 to recover the reformulation (4.8). □

## 5. CONVERGENCE TOWARDS THE VISCOSITY SOLUTION

We conclude our analysis by proving Theorem 2.2, namely the convergence of the solution to $(S^0)$ towards the viscosity solution of $(\mathcal{P}^0)$ and the quantification of the error. As mentioned at the beginning of this work, we first address the semi-discrete case, using the representation formula (4.8), that is equivalent to the semi-discrete scheme $(S_{\Delta t, \Delta v}^0)$ thanks to Proposition 4.14. We establish that the limit of (4.8), when $\Delta t, \Delta v \to 0$ is the unique viscosity solution to $(\mathcal{P}^0)$ by directly comparing the continuous and semi-discrete representation formulas, given by Proposition 4.14 and Theorem 4.9. Finally, we analyze the convergence of the solution to $(S^0)$ towards the one of the semi discrete scheme $(S_{\Delta t, \Delta v}^0)$. The lower part of Figure 2 summarizes these proofs.

### 5.1. Reduction of the discrete action.
We start by showing that the minimal path $\gamma \in \mathbb{R} \times \mathbb{R}^{[0,N_t]}$, defined in Proposition 4.14 admits at most two non-zero velocity jumps. This property is a discrete analog of the following proposition showed in [12].

**Proposition 5.1.** (5.1 in [12]) Consider the representation formula (4.5) with the action $\bar{A}_0^T[\gamma]$ of a minimal path $\gamma \in \mathbb{R} \times \Sigma^{0,T}$, respectively defined in (4.4) and Definition 4.7. In arbitrary dimension, there is at most three states with one intermediate non-trivial one ($s^1, w^1$) with non-zero velocity $w_1 \neq 0$ and at most one flat state with zero velocity. Hence, the action reduces to the following problem:

$$ \left\{ \begin{array}{ll} \text{either } w \equiv w^0 = \frac{x}{T}, \text{ then } \bar{A}_0^T[\gamma] = T \mathbf{1}_{w^0 \neq 0} & (\text{free transport}) \\ \text{or: } \bar{A}_0^T[\gamma] = \frac{|w^2|^2}{2} + \min_{\substack{s^0 w^0 + s^1 w^1 + s^2 w^2 = x \\ s^0, s^1, s^2 \ge 0 \\ s^0 + s^1 + s^2 \le T}} \left\{ \frac{|w^1|^2}{2} + s^0 + s^1 + s^2 \right\}. & \end{array} \right. $$

In the discrete setting, it reads

**Lemma 5.2 (Reduction of the semi-discrete in time action).** Suppose that $\varphi_{\text{in}}$ satisfies Assumption 1. Assume moreover that $\varphi_{\text{in}} - v^2/2$ is L-Lipschitz in $v$. There exists a minimal path $\gamma \in \mathbb{R} \times \mathbb{R}^{[0,N_t]}$, as defined in Definition 4.10, to the minimization problem

$$ \varphi_j^{N_t}(x) = \inf_{\substack{\gamma=(y,w) \in \mathbb{R} \times \mathbb{R}^{[0,N_t]} \\ \gamma(T)=x, w^{N_t}=v_j}} \left\{ \varphi_{\text{in}}(y, w^0) + \bar{A}_0^T[\gamma] \right\}. $$
Moreover, $\gamma$ admits at most one non-trivial velocity state. Denoting $\tilde{s}^k = \tilde{N}^k \Delta t$ the free transport duration at velocity $\tilde{w}^k$, $k = 0, 1, 2$, $\tilde{A}_0^T[\gamma]$ reduces to

$$ \left\{ \begin{array}{ll} \text{either } v_j \equiv w^0 = \frac{x}{T}, \text{ then } \tilde{A}_0^T[\gamma] = T \mathbb{1}_{w^0 \neq 0} \quad (\text{free transport}) \\ \\ \text{or: } \tilde{A}_0^T[\gamma] = \frac{|w^2|^2}{2} + \min_{\substack{(\tilde{s}^0, \tilde{s}^1, \tilde{s}^2) \ge 0 \\ \tilde{s}^0 w^0 + \tilde{s}^1 w^1 + \tilde{s}^2 w^2 = x \\ \tilde{s}^0 + \tilde{s}^1 + \tilde{s}^2 \le T}} \left\{ \frac{|w_1|^2}{2} + \tilde{s}^0 + \tilde{s}^1 + \tilde{s}^2 \right\}. & (5.1) \end{array} \right. $$

The initial velocity $w^0$ and the intermediate one $w^1$ are bounded by $v_j$.

*Proof.* The proof follows closely the proof of the continuous case recalled in Proposition 5.1. The only difference is that one has to deal with time intervals of fixed size $\Delta t$ instead of variable ones.

Let us now consider a path $\gamma = (y, w = (w^0, \dots, w^{N_t})) \in \mathbb{R} \times \mathbb{R}^{[0, N_t]}$. Without changing the cost of $\gamma$ under the discrete action $\tilde{A}$, one can merge all 0 velocity sections ($w^k = 0$) and put them at the end of the time interval yielding a path $\gamma_0$. Then, instead of making redundant jumps, segments with the same velocity are merged, as free transport is generally more cost-effective. More precisely, this holds when $w^k$ satisfies $(w^k)^2 \ge 2N^k \Delta t$ where $N^k$ is the number of time interval spent in free transport. The resulting path is then denoted $\gamma_1$. The case $(w^k)^2 \le 2N^k \Delta t$, which corresponds to relatively small velocities, is considered at the end of the proof, but it does not affect the final results.

Finally, consider a sub-path $\gamma_2 \subset \gamma_1$ consisting of two consecutive segments, $(N^\alpha \Delta t, w^\alpha)$ and $(N^\beta \Delta t, w^\beta)$, with $w^\alpha \ne w^\beta$. We denote by $s^\alpha$ and $s^\beta$ the starting times of these two segments. The sub-path $\gamma_2$ is simplified to a single segment $\gamma_3$ by averaging the velocities over the total duration $(N_\alpha + N_\beta) \Delta t$ into

$$ \bar{w}^\alpha = \frac{N_\alpha \Delta t w^\alpha + (N_\beta - 1) \Delta t w^\beta}{(N_\alpha + N_\beta - 1) \Delta t}. $$

The action of $\gamma_2$ is given by

$$ \tilde{A}_{s^{\alpha}-\Delta t}^{s^{\beta}+N^{\beta}}[\gamma_2] = (N^{\alpha} + N^{\beta} - 1)\Delta t + \frac{(w^{\alpha})^2}{2} + \frac{(w^{\beta})^2}{2} - \Delta t. \quad (5.2) $$

Note that the starting time is $t = s^\alpha - \Delta t$ to account for the fact that $\alpha > 0$: while there is no cost associated with the initial velocity, each subsequent jump does incur a cost. By convexity, one has

$$ |\bar{w}^\alpha|^2 \le \frac{N_\alpha \Delta t}{(N_\alpha + N_\beta - 1) \Delta t} |w^\alpha|^2 + \frac{N_\beta \Delta t}{(N_\alpha + N_\beta - 1) \Delta t} |w^\beta|^2 \le |w^\alpha|^2 + |w^\beta|^2. \quad (5.3) $$

Combining (5.2) and (5.3) one obtains

$$ \tilde{A}_{s^{\alpha}-\Delta t}^{s^{\beta}+N^{\beta}}[\gamma_2] \ge \tilde{A}_{s^{\alpha}-\Delta t}^{s^{\beta}+N^{\beta}}[\gamma_3] - \Delta t. $$

An illustration of this procedure is provided in Figure 3 for clarity. The solid blue path represents an initial trajectory whose cost we aim to reduce. The dotted red line depicts a rearranged path where all segments with zero velocity are shifted to the end, without altering the cost, while the remaining velocity segments are ordered in descending magnitude. This reordering decreases the number of jumps by merging the segments with the same velocities. Finally, the dashed green line shows a further reduction in the number of jumps by introducing a single one at an average velocity, up to introducing an additional $-\Delta t$ term in the cost of the action.

Note that this approach ensures the cost is reduced to at most two jumps: an intermediate one and a final one. However, depending on the values of $x$ and $v$, it may be more advantageous to perform a single large jump, relatively to $x$ and $v$, and remain stationary rather than accumulating the running cost of free transport. This could for example occur for small target velocities $v$.

The existence of a minimal path follows from the reduced form (5.1), which admits a solution. By Assumption 1 and the continuity of $\varphi_{\text{in}}$, the minimization problem involves a continuous coercive function over a closed domain. In particular, this ensures coercivity with respect to $w^0$, while $x$ and $w^2$ act merely as parameters in (5.1).

Furthermore, the inclusion of the initial datum $\varphi_{\text{in}}$ allows one to determine the initial velocity $w^0$. The boundedness of both $w^0$ and $w^1$ by $v_j$ then follows from the definition of the cost of a velocity jump in (4.6): overshooting the target velocity and subsequently dropping down to it would incur a strictly higher cost. $\square$
<figure><img src="image_7.png" /><figcaption>FIGURE 3. Diagram of the reduction of a time discrete path with 4 non-trivial intermediate velocities. Original path (solid blue line), reorganized path (dotted red line), reduced path (dashed green line).</figcaption></figure>

**Remark 2.** (1) Using the same argument as in the discrete case, the initial and intermediate velocities in Proposition 5.1 are also bounded by the final one.

(2) In the proof of Lemma 5.2, the velocities $w$ are defined in $\mathbb{R}^{[0,N_t]}$. This is crucial, as it ensures that the average of two velocities remains within the admissible velocities. We emphasize that this property does not hold when velocities are restricted to a discrete grid. As a result, the reduction of the cost $\tilde{A}[\gamma]$ for a path $\gamma \in \mathbb{R} \times \Delta v \mathbb{Z}^{[0,N_t]}$ cannot be achieved using this method.

(3) Remark that the cost $v^2/2$ of the final jump is always paid except if $x = y + w^0 T$ and $w^0 = v_j$.

(4) The minimal path is not necessarily unique. Indeed, one can shift the location of the zero-velocity segment and any possible non-trivial portion without affecting the cost.

**5.2. Semi-discrete convergence.** As mentioned at the beginning of the section, we first address the semi-discrete case and show that, as $\Delta t, \Delta v \to 0$, the solution given by the representation formula (4.8) converges towards the viscosity solution of $(\mathcal{P}^0)$.

**Lemma 5.3.** Let $x \in \mathbb{R}$ and $(\varphi_j^n(x))_{j \in \mathbb{Z}, n \in [1, N_t]}$ be given by (4.8). Suppose that $\varphi_{\text{in}}$ satisfies Assumption 1 and assume moreover that $\varphi_{\text{in}} - v^2/2$ is L-Lipschitz in $v$. Let $j \in \mathbb{Z}$ be fixed. Then, in the limit $\Delta t, \Delta v \to 0$, $(\varphi_j^n(x))_{n \in [1, N_t]}$ converges towards the viscosity solution of $(\mathcal{P}^0)$. Moreover, there exists $C > 0$, depending on $j$, such that

$$ \forall n \in [1, N_t], \quad |\varphi(t^n, x, v_j) - \varphi_j^n(x)| \le C(\Delta t + \Delta v). \qquad (5.4) $$

*Proof.* The proof of this result is a consequence of the identification of the minimal paths introduced in Section 5.1. We start by recalling the following inclusions:

$$ \mathbb{R} \times \Sigma_{\Delta t, \Delta v}^{0,T} \subset \mathbb{R} \times \Sigma_{\Delta t}^{0,T} \subset \mathbb{R} \times \Sigma^{0,T}. \qquad (5.5) $$

Inequality (5.4) is then proven in two steps.

**Step 1.** We start from (4.8) and use the first inclusion in (5.5) to obtain

$$ \varphi_j^{N_t}(x) \ge \inf_{\substack{\gamma=(y,w) \in \mathbb{R} \times \Sigma_{\Delta t}^{0,T} \\ \gamma(T)=x, w^{N_t}=v_j}} \left\{ \varphi_{\text{in}}(y, w^0) + \tilde{A}_0^T[\gamma] \right\}, $$

which is a consequence of minimizing over a larger set. Using Lemma 5.2, this infimum is in particular reached at a minimal path $\tilde{\gamma} \in \mathbb{R} \times \Sigma_{\Delta t}^{0,T}$ with at most two jumps. We denote by $\tilde{w}^0$ the initial velocity and $\tilde{w}$ the possible non-trivial velocity. It yields:

$$ \varphi_j^{N_t}(x) \ge \varphi_{\text{in}}(y, \tilde{w}^0) + \tilde{A}_0^T[\tilde{\gamma}], $$
<figure><img src="image_8.png" /><figcaption>FIGURE 4. Illustration a minimal path $\bar{\gamma}$ under the action $\bar{A}$ defined in (4.4) (solid red line) and its approximation by a path $\tilde{\gamma}$ defined on the time and velocity grids (dashed blue line).</figcaption></figure>

and we introduce $\bar{A}$ defined in Definition 4.8 to obtain

$$ \varphi_j^{N_t}(x) \geq \inf_{\substack{\gamma=(y,\mathbf{w})\in\mathbb{R}\times\Sigma^{0,T} \\ \gamma(T)=x, \dot{\gamma}(T)=v_j}} \{\varphi_{\text{in}}(y, w^0) + \bar{A}_0^T[\gamma]\} + (\tilde{A}_0^T[\tilde{\gamma}] - \bar{A}_0^T[\tilde{\gamma}]) \\ = \varphi(T, x, v_j) + \tilde{A}_0^T[\tilde{\gamma}] - \bar{A}_0^T[\tilde{\gamma}], $$

One can then apply Lemma 4.13 to the difference $\tilde{A}_0^T[\tilde{\gamma}] - \bar{A}_0^T[\tilde{\gamma}]$. Since $\tilde{\gamma}$ admits at most two velocity jumps, one obtains

$$ \varphi_j^{N_t}(x) \geq \varphi(T, x, v_j) - 2\Delta t. \qquad (5.6) $$

**Step 2.** We start from (4.5) in which we choose a minimal path $\bar{\gamma} \in \mathbb{R} \times \Sigma^{0,T}$ with $\bar{\gamma}(0) = \bar{y}$. Thanks to Proposition 5.1, it has at most two jumps. In the following, we assume that $\bar{\gamma}$ admits two jumps but the same analysis can be done for one or no jumps. We denote by $\bar{w}^0$ and $\bar{w}^1$ the associated velocities. One has:

$$ \varphi(T, x, v_j) = \varphi_{\text{in}}(\bar{y}, \bar{w}^0) + \bar{A}_0^T[\bar{\gamma}]. \qquad (5.7) $$

In what follows, we define $\tilde{\gamma} \in \mathbb{R} \times \Sigma^{0,T}_{\Delta t, \Delta v}$, such that it has a similar cost compared to $\bar{\gamma} \in \mathbb{R} \times \Sigma^{0,T}$ using the action $\bar{A}$. It is chosen with the same number of velocity jumps as $\bar{\gamma}$, that is at most 2 thanks to Proposition 5.1. The velocities in $\tilde{\gamma}$ are the velocities of $\bar{\gamma}$ shifted towards the nearest neighbor on the grid $\Delta v\mathbb{Z}$. Similarly, the jump times in $\tilde{\gamma}$ are the jump times in $\bar{\gamma}$ shifted to their nearest neighbor on the grid $\{n\Delta t, 0 \le n \le N_t\}$. We denote $\tilde{w}^0 = \tilde{\gamma}(0) \in \Delta v\mathbb{Z}$ and $\tilde{w}^1 \in \Delta v\mathbb{Z}$ the initial and intermediate velocity of this path. Its starting point is denoted $\tilde{y} = \tilde{\gamma}(0) \in \mathbb{R}$. Without loss of generality, let us assume that the trivial state with zero velocity is located at the end of the time interval. Then, by construction of $\tilde{\gamma}$, one has the following estimates:

$$ |\tilde{s}^1 - \bar{s}^1| \le \Delta t, \quad |\tilde{s}^0 - \bar{s}^0| \le \Delta t, \quad |\tilde{w}^1 - \bar{w}^1| \le \Delta v, \text{ and } |\tilde{w}^0 - \bar{w}^0| \le \Delta v. \qquad (5.8) $$

The following observation can then be made: the definition of $\tilde{\gamma}$ introduces an error on $y^0$, $\bar{w}^0$, $\bar{w}^1$ and on the jump times $\bar{s}^0$ and $\bar{s}^1$. Figure 4 illustrates this remark, and the notations. To estimate these errors, one follows the characteristic curves backward at velocities $\tilde{w}^1$ and $\bar{w}^1$ for the durations $\tilde{s}_1 = \tilde{s}^1 - \tilde{s}^0$ and $\bar{s}_1 = \bar{s}^1 - \bar{s}^0$ respectively. We denote $\tilde{x}_1 = x - \tilde{w}^1\tilde{s}_1$ and $\bar{x}_1 = x - \bar{w}^1\bar{s}_1$ the feet of the intermediate states of $\tilde{\gamma}$ and $\bar{\gamma}$, so that

$$ |\bar{x}_1 - \tilde{x}_1| \le 2v_j\Delta t + T\Delta v, $$

thanks to (5.8) and $\bar{w}^1 \le v_j$ (Remark 2). As $\tilde{\gamma}$ and $\bar{\gamma}$ have at most two jumps, the same reasoning can be applied to the sections with velocities $\tilde{w}^0$ and $\bar{w}^0$, so that

$$ |\tilde{y} - \bar{y}| = |\tilde{\gamma}(0) - \bar{\gamma}(0)| \le C(\Delta t + \Delta v), $$

where the constant $C$ depends on $T$ and $j \in \mathbb{Z}$.
Introducing $\tilde{\gamma}$ in (5.7) yields

$$ \varphi(T, x, v_j) = \varphi_{\text{in}}(\tilde{y}, \tilde{w}^0) + \tilde{A}_0^T[\tilde{\gamma}] + (\varphi_{\text{in}}(\tilde{y}, \bar{w}^0) - \varphi_{\text{in}}(\tilde{y}, \tilde{w}^0)) + (\bar{A}_0^T[\tilde{\gamma}] - \bar{A}_0^T[\tilde{\gamma}]) + (\bar{A}_0^T[\tilde{\gamma}] - \tilde{A}_0^T[\tilde{\gamma}]). $$

Using the inclusion (5.5), one obtains

$$ \begin{aligned} \varphi(T, x, v_j) &\ge \inf_{\substack{\gamma=(y, \mathbf{w}) \in \mathbb{R} \times \Sigma_{\Delta_t, \Delta v}^{0, T} \\ \gamma(T)=x, w^{N_t}=v}} \left\{ \varphi_{\text{in}}(y, w^0) + \tilde{A}_0^T[\gamma] \right\} \\ &\quad + (\varphi_{\text{in}}(y, \bar{w}^0) - \varphi_{\text{in}}(\tilde{y}, \tilde{w}^0)) + (\bar{A}_0^T[\tilde{\gamma}] - \bar{A}_0^T[\tilde{\gamma}]) + (\bar{A}_0^T[\tilde{\gamma}] - \tilde{A}_0^T[\tilde{\gamma}]) \\ &= \varphi_j^{N_t}(x) + (\varphi_{\text{in}}(\tilde{y}, \bar{w}^0) - \varphi_{\text{in}}(\tilde{y}, \tilde{w}^0)) + (\bar{A}_0^T[\tilde{\gamma}] - \bar{A}_0^T[\tilde{\gamma}]) + (\bar{A}_0^T[\tilde{\gamma}] - \tilde{A}_0^T[\tilde{\gamma}]). \end{aligned} \quad (5.9) $$

We now estimate the three differences in (5.9). Firstly, the L-Lipschitz continuity of $\varphi_{\text{in}} - v^2/2$ and inequalities (5.8) yield

$$ |\varphi_{\text{in}}(\tilde{y}, \bar{w}^0) - \varphi_{\text{in}}(\tilde{y}, \tilde{w}^0)| \le C L(\Delta t + \Delta v) + \frac{1}{2} \Delta v |\bar{w}^0 + \tilde{w}^0|. $$

As both $\bar{w}^0$ and $\tilde{w}^0$ are less than the target velocity $v_j$ (Lemma 5.2 and Remark 2), there exists a constant $C$, that depends on $j$, such that

$$ |\varphi_{\text{in}}(\tilde{y}, \bar{w}^0) - \varphi_{\text{in}}(\tilde{y}, \tilde{w}^0)| \le C(\Delta t + \Delta v). \quad (5.10) $$

Secondly, we deal with the difference of actions $\bar{A}_0^T[\tilde{\gamma}] - \bar{A}_0^T[\tilde{\gamma}]$. By definition, the paths $\tilde{\gamma}$ and $\tilde{\tilde{\gamma}}$ admit two jumps with the final one being shared. The difference reads,

$$ |\bar{A}_0^T[\tilde{\gamma}] - \bar{A}_0^T[\tilde{\tilde{\gamma}}]| \le |\bar{s}^0 + \bar{s}^1 - \bar{s}^0 - \bar{s}^1| + \left| \frac{(\bar{w}^1)^2}{2} - \frac{(\tilde{w}^1)^2}{2} \right| + |\mathbf{1}_{\{v_j \ne \frac{x}{T}, \bar{w}^0 \ne v_j\}} - \mathbf{1}_{\{v_j \ne \frac{x}{T}, \tilde{w}^0 \ne v_j\}}| \frac{v_j^2}{2}, $$

where the indicator functions take the value 1 if the condition is satisfied and 0 otherwise. Using (5.8), the differences in time and velocity are of order $\Delta t$ and $\Delta v$, respectively. Regarding the difference of indicator functions encoding the last jump, we separate two cases. First, if $v_j \ne \frac{x}{T}$, the difference is zero because a final jump must be paid to reach the target velocity both in $\tilde{\gamma}$ and $\tilde{\tilde{\gamma}}$, and they have the same cost. Secondly, if $v_j = \frac{x}{T}$, the cases $\bar{w}^0 \ne v_j$ and $\tilde{w}^0 \ne v_j$ again both require final a jump. Lastly, assume that $\bar{w}^0 = v_j \in \Delta v\mathbb{Z}$, then $\tilde{w}^0 = v_j$. Consequently, the indicator functions both equal zero at the same time and the last jump never introduces an error. It yields

$$ |\bar{A}_0^T[\tilde{\gamma}] - \bar{A}_0^T[\tilde{\tilde{\gamma}}]| \le 2\Delta t + (\bar{w}^1 + \tilde{w}^1) \frac{\Delta v}{2} \le 2\Delta t + 2v_j \Delta v. \quad (5.11) $$

Thirdly, the difference $\bar{A}_0^T[\tilde{\gamma}] - \bar{A}_0^T[\tilde{\tilde{\gamma}}]$ is bounded by $2\Delta t$ thanks to Lemma 4.13. Combining this bound with (5.11) and (5.10), there exists a constant $C$, depending on the final time $T$ and on the target velocity $v_j$, such that

$$ \varphi(T, x, v_j) \ge \varphi_j^{N_t}(x) - C(\Delta t + \Delta v). \quad (5.12) $$

Finally, by combining (5.6) and (5.12), and by equivalence of (4.8) and $(S^0_{\Delta t, \Delta v})$, we conclude the proof of Lemma 5.3. $\square$

**5.3. Fully discrete convergence.** In this section, we prove the convergence of the solution to $(S^0)$ towards the one of $(S^0_{\Delta t, \Delta v})$ when $\Delta x \to 0$ with fixed $\Delta t$ and $\Delta v$.

**Lemma 5.4.** Let $(\mu_i^n, \varphi_{ij}^n)_{i \in \mathbb{Z}, j \in \mathbb{Z}, n \in [1, N_t]}$ be given by $(S^0)$. Suppose that $\varphi_{\text{in}}$ satisfies Assumption 1 and that $\varphi_{\text{in}} - v^2/2$ is L-Lipschitz in $v$. Then, in the limit $\Delta x \to 0$, the solution to $(S^0)$ converges towards the solution to the semi-discrete scheme $(S^0_{\Delta t, \Delta v})$: for all $(i, j) \in \mathbb{Z}^2$, and all $n \in [1, N_t]$,

$$ |\mu^n(x_i) - \mu_i^n| \le 2LT \frac{\Delta x}{\Delta t}, \quad (5.13) $$

$$ |\varphi_j^n(x_i) - \varphi_{ij}^n| \le 2LT \frac{\Delta x}{\Delta t}. \quad (5.14) $$
Before beginning the proof, it is worth mentioning that the transport in $x$ is exactly resolved, meaning that no interpolations are needed, if $\Delta t\Delta v = \Delta x$. In this case, the results of Lemma 5.4 are true for $(S^0)$. However, condition (2.19) then implies a huge computational cost, as $\Delta x \le \Delta v^3/2$.

*Proof.* Proceeding by induction, we start by estimating the propagation of interpolation errors over a single time step. Assume that for some $n \ge 0$, and for all $(i, j) \in \mathbb{Z}^2$,

$$e_{\mu,i}^n := \mu^n(x_i) - \mu_i^n, \quad |e_{\mu,i}^n| \le 2Ln\Delta x, \qquad (5.15)$$

$$e_{\varphi,ij}^n := \varphi_j^n(x_i) - \varphi_{ij}^n, \quad |e_{\varphi,ij}^n| \le 2Ln\Delta x, \qquad (5.16)$$

where $(\mu_i^n, \varphi_{ij}^n)$ is obtained from $(S^0)$ and $(\mu^n(x_i), \varphi_j^n(x_i))$ from $(S^0_{\Delta t, \Delta v})$. Remark that the interpolation (2.13) applied to a function $h = h(x)$ that is $L$-Lipschitz yields an error of size $2L\Delta x$. Applying this result to $\bar{\varphi}_{ij,j}^n$, given by (2.13), and Corollary 4.3 on the propagation of Lipschitz bounds, one obtains:

$$|\varphi_j^n(x_i - \Delta t v_j) - \bar{\varphi}_{ij,j}^n| \le 2L\Delta x, \qquad (5.17)$$

where $L$ is the Lipschitz constant of $\varphi_{in}$. In particular, introducing $\varphi_j^n(x_{ij})$ in the interpolation (2.13), one has

$$\bar{\varphi}_{ij,j}^n = \alpha_j \varphi_j^n(x_{ij}) + (1 - \alpha_j) \varphi_j^n(x_{ij-1}) + \alpha_j e_{\varphi,ij,j}^n + (1 - \alpha_j) e_{\varphi,ij-1,j}^n \qquad (5.18)$$

Let us now show the propagation of (5.15). Introducing $\varphi_j^n(x_i - \Delta t v_j)$ and (5.18) in the update relation on $\mu_i^{n+1}$ in $(S^0)$ one obtains

$$\mu_i^{n+1} = \min_{j \in \mathbb{Z}} \left\{ (\alpha_j \varphi_j^n(x_{ij}) + (1 - \alpha_j) \varphi_j^n(x_{ij-1}) - \varphi_j^n(x_i - \Delta t v_j)) + \alpha_j e_{\varphi,ij,j}^n + (1 - \alpha_j) e_{\varphi,ij-1,j}^n + \varphi_j^n(x_i - \Delta t v_j) \right\}.$$

Using the interpolation error estimate (5.17), it yields

$$\mu^{n+1}(x_i) - 2L\Delta x - \sup_{(i,j) \in \mathbb{Z}^2} \{|e_{\varphi,ij}^n|\} \le \mu_i^{n+1} \le \mu^{n+1}(x_i) + 2L\Delta x + \sup_{(i,j) \in \mathbb{Z}^2} \{|e_{\varphi,ij}^n|\},$$

from which one obtains (5.15) at step $n+1$ using the induction hypothesis (5.16) on $e_{\varphi,ij}^n$. The same ideas are then applied to the equation on $\varphi_{ij}^{n+1}$ in $(S^0)$. Let

$$A := \varphi_j^n(x_i - \Delta t v_j) + (\bar{\varphi}_{ij,j}^n - \varphi_j^n(x_i - \Delta t v_j)) + \Delta t,$$

$$B := \frac{v_j^2}{2} + \mu^{n+1}(x_i) + (\mu_i^{n+1} - \mu^{n+1}(x_i)),$$

such that $\varphi_{ij}^n = \min\{A, B\}$. Applying estimate (5.17) to $A$ and (5.15) at step $n+1$ to $B$, one has in particular

$$A \le \varphi_j^n(x_i - \Delta t v_j) + \Delta t + 2L\Delta x, \qquad (5.19)$$

$$A \ge \varphi_j^n(x_i - \Delta t v_j) + \Delta t - 2L\Delta x, \qquad (5.20)$$

$$B \le \frac{v_j^2}{2} + \mu^{n+1}(x_i) + 2L(n+1)\Delta x, \qquad (5.21)$$

$$B \ge \frac{v_j^2}{2} + \mu^{n+1}(x_i) - 2L(n+1)\Delta x. \qquad (5.22)$$

Therefore, combining (5.19) and (5.21), one obtains the upper bound

$$\varphi_{ij}^{n+1} \le \varphi_j^{n+1}(x_i) + 2L(n+1)\Delta x. \qquad (5.23)$$

Similarly, estimates (5.20) and (5.22) yield the lower bound

$$\varphi_{ij}^{n+1} \ge \varphi_j^{n+1}(x_i) - 2L(n+1)\Delta x. \qquad (5.24)$$

Finally, (5.16) is obtained at step $n+1$ combining (5.23) and (5.24). $\square$

Using Lemmas 5.3 and 5.4, one can conclude the proof of Theorem 2.2.

**Remark 3.** On the generalization to higher dimensions. The generalization of Theorems 2.1 and 2.2 to higher dimensions in velocity is straightforward, as all velocity-specific tools naturally generalize to this setting. However, when considering the position variable, the monotonicity of the interpolation (2.13) played a crucial role in the analysis of both $(S^\epsilon)$ and $(S^0)$. Therefore, the multidimensional interpolation must be carefully selected to ensure this property is preserved. On
suitably structured meshes, a linear interpolation in each coordinate direction should be sufficient. Additionally, since the unknown is expected to exhibit only Lipschitz regularity, this imposes an extra constraint on the interpolation method. Finally, the commutation of the interpolation (2.13) with constants and translations allowed the exact preservation of Lipschitz bounds. Without these last two properties, this preservation may deteriorate over time.

## 6. NUMERICAL RESULTS

In this section, we present numerical experiments validating the properties of schemes ($S^\epsilon$) and ($S^0$), as established in Theorems 2.1 and 2.2. While these results were proven in an unbounded phase space, practical considerations require to restrict the computational domain.

For the spatial discretization, for the sake of simplicity, we impose periodic boundary conditions on the domain $[-x_*, x_*]$, with a uniform grid spacing defined as

$$ \Delta x = \frac{2x_*}{N_x}, \quad x_i = -x_* + \left(i - \frac{1}{2}\right) \Delta x, \quad i \in [1, N_x]. $$

A similar discretization is applied in velocity but no boundary condition needs to be defined as there is no transport in this variable:

$$ \Delta v = \frac{2v_*}{N_v}, \quad v_j = -v_* + \left(j - \frac{1}{2}\right) \Delta v, \quad j \in [1, N_v]. $$

Unless stated otherwise, we set the discretization parameters as follows:

$$ N_x = 64, \quad N_v = 61, \quad x_* = 10, \quad v_* = 10. $$

Note that to ensure that $v = 0$ belongs to the grid, we choose an odd number of points in velocity. The time step is chosen to satisfy the stability condition (2.19): $\Delta t = 0.9\Delta v/2$. In the following, we consider the initial data

$$ \varphi_{\text{in}}(t, x, v) = \min \{3(v-3)^2 + 5, 5(v+7)^2 + 2\} + 0.9 \cos \left( \frac{4\pi xv}{x_* v_*} \right), \quad (6.1) $$

which enjoys Assumption 1, and is such that $\varphi_{\text{in}} - v^2/2$ is L-Lipschitz in the $v$ variable. To illustrate the asymptotic-preserving property of ($S^\epsilon$), we compare it with a standard explicit finite difference scheme on ($\mathcal{P}^\epsilon$). Denoting $f_{ij}^n$ an approximation of $f(t^n, x_i, v_j)$, the scheme is given by

$$ \frac{f_{ij}^{n+1} - f_{ij}^n}{\Delta t} + (v \partial_x f^n)_{ij} = \frac{1}{\epsilon} (\rho_i^n \mathcal{M}_j^\epsilon - f_{ij}^n), \quad (6.2) $$

where

$$ (v \partial_x f^n)_{ij} = v_j^+ \frac{f_{ij}^n - f_{i-1,j}^n}{\Delta x} + v_j^- \frac{f_{i+1,j}^n - f_{ij}^n}{\Delta x} $$

is a classical upwind scheme for the transport term, $\rho_i^n = \langle f_i^n \rangle_{\Delta v}$, and $\mathcal{M}_j^\epsilon = c\mathcal{M}^\epsilon(v_j)$ with $c = \langle \mathcal{M}^\epsilon \rangle_{\Delta v}^{-1}$ being a normalization constant. To compare the results of scheme (6.2) to the ones of ($S^\epsilon$), we introduce the Hopf-Cole transform of the solution to (6.2): for all $(n, i, j) \in [0, N_t] \times [1, N_x] \times [1, N_v]$,

$$ \psi_{ij}^n = -\epsilon \log(f_{ij}^n). $$

Since the scheme (6.2) is explicit, its stability is guaranteed only under a CFL condition of the form $\Delta t \le \epsilon C \Delta x$, which becomes excessively restrictive as $\epsilon \to 0$. Moreover, due to the concentration phenomena in velocity as $\epsilon \to 0$, this scheme is expected to poorly resolve this variable in that regime. Another challenge arises from the use of the logarithm in the Hopf-Cole transform: for small values of $\epsilon$, the distribution tails approach the machine epsilon, leading to significant numerical errors. Consequently, the scheme (6.2) is expected to approach accurately the solution to ($\mathcal{P}^\epsilon$) only in the regime $\epsilon \sim 1$.

**Remark 4.** With the discretization parameters presented above, the exact preservation of the equilibrium $v^2/2$ is observed. More precisely, setting $\varphi_{\text{in}} = v^2/2$, the quantities $\max_{1 \le n \le N_t} \{\varphi^\epsilon - v^2/2\}$ and $\max_{1 \le n \le N_t} \{\varphi - v^2/2\}$ are exactly equal to 0.
<figure><img src="image_9.png" /><figcaption>FIGURE 5. Evolution of the solution to (6.2) (Kinetic, top) and ($S^ε$) (bottom) compared to the solution to ($S^0$) for $v = 0$ at $t = 0.07$ (left), 0.15 (center), and 1.5 (right) for different values of $ε$.</figcaption></figure>

6.1. **AP property.** We start by evaluating the AP property of ($S^ε$). This is done through a comparison with both the limiting scheme ($S^0$) and the explicit (non-AP) scheme (6.2).

Figures 5 and 6 illustrate the evolution of $\varphi^ε$, $\psi$ and $\varphi$ under various values of the scaling parameter $ε$. In Figure 5, $\varphi^ε$, $\psi$ and $\varphi$ are shown as a function of $x$ with fixed velocity $v = 0$. We observe that the explicit scheme (top row) fails to capture the correct limiting behavior as $ε \to 0$, whereas the AP scheme ($S^ε$) (bottom row) is able to align with the asymptotic solution. In addition, it is worth noticing that, as $ε$ decreases, sharp spatial peaks emerge, and a Lipschitz-type regularity in $x$ becomes apparent in the limiting profile. This behavior is typical of Hamilton-Jacobi equation.

Figure 6 presents the same simulation but now with fixed spatial coordinate $x = 0$, displaying $\varphi^ε$, $\psi$ and $\varphi$ as a function of $v$. Once again the AP scheme approximates the limit well, while the explicit scheme exhibits notable discrepancies. In particular, the explicit method performs poorly for large velocities. This can be attributed to the procedure it follows—first computing the distribution function and then taking the logarithm to recover $\psi$. At high velocities, this approach suffers from numerical instabilities due to values near machine precision, leading to inaccuracies in the logarithm.

Finally, the convergence of the AP scheme as $ε \to 0$ is quantified in Figure 7. The pointwise $L^∞$-error between $\varphi^ε$ and $\varphi$ converges as $ε \to 0$, as predicted in Theorem 2.1. While Theorem 2.1 does not provide an explicit rate of convergence with respect to $ε$, a rate of 1 can be observed for this test case.

6.2. **Convergence.** We now illustrate the convergence behavior of schemes ($S^ε$) and ($S^0$) with respect to the discretization parameters $Δt$, $Δv$, and $Δx$. The error associated to a given triplet of parameters $Δ = (Δt, Δv, Δx)$ is defined by

$$E(\varepsilon, \Delta) = \frac{|\varphi_{\text{Ref}}^{\varepsilon} - u_{\Delta}^{\varepsilon}|_{\infty}}{|\varphi_{\text{Ref}}^{\varepsilon}|_{\infty}},$$

where $E(0, Δ)$ corresponds to the error of the limiting scheme and $u_{\Delta}^{\varepsilon} = \varphi^{\varepsilon}$, $\varphi$ or $\psi$ computed with ($S^ε$), ($S^0$) and (6.2) respectively. The reference solution $\varphi_{\text{Ref}}^{\varepsilon}$ is computed using a sufficiently fine discretization, while $u_{\Delta}^{\varepsilon}$ is obtained using varying values of $Δ$. All convergence tests use (6.1) as
<figure><img src="image_10.png" /><figcaption>FIGURE 6. Comparison of the solutions to ($S^\epsilon$), ($S^0$) and (6.2), for $x = 0$, at $t = 0.075$ (top left), 0.15 (top right), 1.5 (bottom left) and 3.0 (bottom right) for different values of $\epsilon$.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>FIGURE 7. Convergence in $\epsilon$. $L^\infty$-error between ($S^\epsilon$) and ($S^0$), as a function of $\epsilon$ for different times</figcaption></figure>

initial data. The convergence curves associated to the schemes ($S^\epsilon$), ($S^0$), and the explicit scheme (6.2) are displayed in Figure 8.

To study the convergence with respect to the spatial step $\Delta x$, we fix $N_v = 201$, $\Delta t = 2.5 \times 10^{-5}$, and $T = 0.01$. The reference solution is computed using $N_x = 2^{15}$. The AP scheme ($S^\epsilon$) for $\epsilon = 1$ and limit scheme ($S^0$) both exhibit an observed convergence order of approximately 0.7. The deterioration of the expected order 1 of these two schemes can be explained by the Lipschitz (and not $C^1$) regularity of the approximated solution. Nevertheless, the schemes do converge. For scheme ($S^0$), the results should be interpreted in light of the error estimate (2.20) from Theorem 2.2, where the ratio $\Delta x/\Delta t$ appears. As a result, the influence of $\Delta x$ also cannot be fully separated from that of $\Delta t$.
<figure><img src="image_12.png" /><figcaption>FIGURE 8. Convergence study of the schemes ($S^\epsilon$), ($S^0$) and (6.2). Values of $E(1, \Delta)$ and $E(0, \Delta)$ as a function of $\Delta t$ (top left), $E(1, \Delta)$ and $E(0, \Delta)$ as a function of $\Delta v$ (top right) and $E(1, \Delta)$ and $E(0, \Delta)$ as a function of $\Delta x$ (bottom).</figcaption></figure>

The convergence with respect to the velocity step $\Delta v$ is investigated by setting $N_x = 256$, $\Delta t = 2.5 \times 10^{-5}$, and $T = 0.01$. A reference solution is computed using $N_v = 3^{10}$. In this case, the error curves for all schemes are superimposed, and an order of 1 is observed, in agreement with the theoretical prediction from (2.20).

To analyze the convergence with respect to the time step $\Delta t$, we fix $N_v = 101$ and $T = 0.5$, and compute the reference solution using $N_t = 2^{11}$. As noted in the spatial convergence study, $\Delta t$ also appears in the denominator in (2.20), so letting $\Delta t \to 0$ requires adapting $\Delta x$ accordingly to avoid divergence of the ratio $\Delta t/\Delta x$. As discussed in Section 5.3, if the relation $\Delta t\Delta v = \Delta x$ is satisfied, the transport step defined in (2.13) becomes exact, and the term $\Delta t/\Delta x$ in (2.20) vanishes. We adopt this setting and determine the number of grid points from $N_x = \lfloor 2x_*/\Delta x \rfloor$. With this choice, the convergence in $\Delta t$ can be properly assessed, and the predicted first-order behavior in time is observed. Note that the errors of ($S^\epsilon$) and ($S^0$) actually coincide. In addition, the naive scheme (6.2) fails to achieve first-order accuracy in time for this test case, highlighting the need for an appropriate discretization of the evolution equation for $\phi^\epsilon$.

### 6.3. Decreasing of the minimum.
We now investigate the decay of the minimum of $\phi$. Recall that the limit system ($P^0$) satisfies the condition $\partial_t \min_v \{\phi\} \le 0$. As shown in Lemma 4.5, the discrete system ($S^0$) also preserves a discrete analogue of this property. Figure 9 displays the evolution of $\min_j \{\phi^\epsilon\}$ for the limit scheme. For comparison, the evolution of $\min_j \{\phi^\epsilon\}$ computed with ($S^\epsilon$) for $\epsilon = 1$ is also displayed in Figure 9. For the limit scheme, we observe a monotonic decrease of the minimum, and the solution appears to reach a non-trivial stationary state. In contrast, the simulation with $\epsilon = 1$ shows that the minimum can increase, confirming that this monotonicity is indeed a behavior specific to the asymptotic regime.

### 6.4. Persistence of the initial data.
In this final numerical experiment, we examine the long-time behavior of the solution. It is well known (see, e.g., [22]) that for $\epsilon = 1$, the solution $f^\epsilon$ to ($P^\epsilon$) converges, as $t \to \infty$, towards a stationary state given by

$$f_\infty = \rho_\infty \mathcal{M}(v),$$
<figure><img src="image_13.png" /><figcaption>FIGURE 9. Evolution of $\min_v \{\varphi\} (t, x)$ (left) and $\min_v \{\varphi^\varepsilon\} (t, x)$ for $\varepsilon = 1$ (right), as functions of $x$ at different times.</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>FIGURE 10. Long-time (left) and short-time (right) spatial amplitude evolution of the solution, evaluated at $v = 0$, corresponding to $(S^\varepsilon)$, $(S^0)$ and (6.2), for various values of $\varepsilon$.</figcaption></figure>

where $\rho_\infty$ is a positive constant that depends on the initial mass and the length of the periodic spatial domain $\mathbb{T}$. This convergence is driven by the mixing effect resulting from the interplay between the transport operator $v \cdot \nabla_x$ and the relaxation term $\rho M - f$. Notably, this mixing mechanism remains active even when $\varepsilon > 0$, i.e., after the rescaling in $(P^\varepsilon)$. However, as already hinted at in Figure 1, a striking observation is that in the limit $\varepsilon \to 0$, this behavior no longer holds: the stationary profile in position is no longer uniform but retains features of the initial data. Mathematically, this reflects the non-commutation of the limits $\varepsilon \to 0$ and $t \to \infty$.

To illustrate this phenomenon, we consider the initial data (6.1). Figure 11 compares, at fixed velocity $v = 0$, the long-time behavior of the solution to $(S^0)$ with that of $(S^\varepsilon)$ for $\varepsilon = 1$ and $\varepsilon = 0.1$. As $\varepsilon$ decreases, the mixing effect becomes progressively weaker. Another way to observe and quantify this weakening as $\varepsilon \to 0$ is through the spatial amplitude of the solution. The left panel of Figure 10 shows the evolution of

$$ \max_{x \in \mathbb{T}} u(t, x, 0) - \min_{x \in \mathbb{T}} u(t, x, 0), \qquad (6.3) $$

computed from the solutions to $(S^\varepsilon)$, $(S^0)$, and (6.2). The weakening of the mixing effect is evident: as $\varepsilon$ decreases, it takes increasingly longer for the amplitude to decay to a given value. Furthermore, the asymptotic-preserving property of $(S^\varepsilon)$ is clearly visible, as it correctly captures the convergence toward the limiting behavior. In contrast, the classical scheme (6.2) does converge to a stationary profile with non-zero amplitude, but this profile does not correspond to that of the limiting scheme. On the right panel of Figure 10, only $(S^\varepsilon)$ and $(S^0)$ are compared. In particular, one can observe that at short times there appears to be a characteristic point of inflexion in the curve (between $t = 1$ and $t = 1.5$), where the dynamics intersect. This highlights the change in the system's behavior induced by the scaling (1.2).
<figure><img src="image_15.png" /><figcaption>FIGURE 11. Evolution of $\varphi^\epsilon(t, x, 0)$ for $\epsilon = 1$ (top left), $\varphi^\epsilon(t, x, 0)$ for $\epsilon = 0.1$ (top right), and $\varphi(t, x, 0)$ (bottom), as functions of $x$ at different times.</figcaption></figure>

An even more striking illustration of this property is obtained by looking at the fundamental solution of (P⁰), starting—at the distributional level—from a Dirac mass in position and a Gaussian in velocity. Numerically, since computing the logarithm of a Dirac mass is not feasible for initializing $\varphi$, the initial data $\varphi_{\text{in}}$ is instead set to zero on a few grid cells centered around $x = 0$, and to 100 elsewhere.

The asymptotic behavior as $t \to \infty$ was characterized in the one-dimensional case in [12]. The action between a point $x$ and the origin is given by

$$ \bar{A}_0^t[\bar{\gamma}] = \begin{cases} \frac{3}{2}|x|^{2/3}, & \text{if } |x| \le t^{3/2}, \\ \frac{|x|^2}{2t^2} + t, & \text{if } |x| \ge t^{3/2}. \end{cases} \qquad (6.4) $$

The evolution of $\varphi$ is shown in Figure 12. As $t \to \infty$, the spatial profile converges to the cusp-shaped function

$$ x \mapsto \frac{3}{2}|x|^{2/3}, $$

as predicted by the first case of (6.4).

Figure 12 also presents the evolution from two Dirac masses. In this case, the solution converges to the superposition of two cusps, and a clear imprint of the initial data remains visible in the long-time regime. A similar phenomenon is illustrated in Figure 11: for $\epsilon = 0$, the oscillatory structure of the initial condition (6.1) persists as $t \to \infty$, whereas for $\epsilon > 0$, the solution tends to relax and lose memory of the initial state due to mixing. It is also worth noting that, as $\epsilon$ decreases, relaxation towards a constant occurs more slowly—larger times are required to observe this behavior, as seen in the top left and top right panels of Figure 11. This persistent structure in the limiting regime is further illustrated in Figure 13, which shows the long-time behavior of the solution to (S⁰) in the $(x, v)$-phase space: the profile converges to the parabolic shape $v^2/2$ in velocity, while remaining non-trivial and structured in position.
<figure><img src="image_16.png" /><figcaption>FIGURE 12. Evolution of the logarithm of one (left) and two (right) Dirac masses, as functions of $x$ at different times.</figcaption></figure>

## ACKNOWLEDGMENT

The authors would like to warmly thank Vincent Calvez for the fruitful discussions and insights on this project.

The authors have received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation program (grant agreement No 865711).

## REFERENCES

[1] ALMEIDA, L., PERTHAME, B., AND RUAN, X. An asymptotic preserving scheme for capturing concentrations in age-structured models arising in adaptive dynamics. *J. Comput. Phys.* **464** (2022), 26. Id/No 111335.

[2] BARDOS, C. W., GOLSE, F., AND LEVERMORE, D. M. Fluid dynamic limits of kinetic equations. i. formal derivations. *J. Statist. Phys.* **63** (1991), 323–344.

[3] BARLES, G. An introduction to the theory of viscosity solutions for first-order Hamilton-Jacobi equations and applications. In *Hamilton-Jacobi equations: approximations, numerical analysis and applications. Based on the lectures of the CIME summer school, Cetraro, Italy, August 29–September 3, 2011*. Berlin: Springer; Florence: Fondazione CIME, 2013, pp. 49–109.

[4] BARLES, G., EVANS, L. C., AND SOUGANDIS, P. E. Wavefront propagation for reaction-diffusion systems of PDE. *Duke Math. J.* **61**, 3 (1990), 835–858.

[5] BARLES, G., MIRRAHIMI, S., AND PERTHAME, B. Concentration in Lotka-Volterra parabolic or integral equations: a general convergence result. *Methods Appl. Anal.* **16**, 3 (2009), 321–340.

[6] BARLES, G., AND PERTHAME, B. Concentrations and constrained Hamilton-Jacobi equations arising in adaptive dynamics. In *Recent developments in nonlinear partial differential equations. Proceedings of the second symposium on analysis and PDEs, West Lafayette, IN, USA, June 7–10, 2004*. Providence, RI: American Mathematical Society (AMS), 2007, pp. 57–68.

[7] BARLES, G., AND SOUGANIDIS, P. E. Convergence of approximation schemes for fully nonlinear second order equations. *Asymptotic Anal.* **4**, 3 (1991), 271–283.

[8] BENNOUNE, M., LEMOU, M., AND MIEUSSENS, L. Uniformly stable numerical schemes for the Boltzmann equation preserving the compressible Navier-Stokes asymptotics. *J. Comput. Phys.* **227** (2008).

[9] BENSOUSSAN, A., LIONS, J.-L., AND PAPANICOLAU, G. C. Boundary Layers and Homogenization of Transport Processes. *Publ. Res. I. Math. Sci.* (1979).

[10] BOUIN, E. A Hamilton-Jacobi approach for front propagation in kinetic equations. *Kinet. Relat. Models* **8**, 2 (2015), 255–280.

[11] BOUIN, E., AND CALVEZ, V. A kinetic eikonal equation. *C. R., Math., Acad. Sci. Paris 350*, 5-6 (2012), 243-248.

[12] BOUIN, E., CALVEZ, V., GRENIER, E., AND NADIN, G. Large-scale asymptotics of velocity-jump processes and nonlocal Hamilton-Jacobi equations. *J. Lond. Math. Soc., II. Ser.* **108**, 1 (2023), 141–189.

[13] BOUIN, E., CALVEZ, V., AND NADIN, G. Propagation in a kinetic reaction-transport equation: travelling waves and accelerating fronts. *Arch. Ration. Mech. Anal.* **217**, 2 (2015), 571–617.

[14] CAILLERIE, N. Large deviations of a velocity jump process with a Hamilton-Jacobi approach. *C. R., Math., Acad. Sci. Paris 355*, 2 (2017), 170-175.

[15] CALVEZ, V., HENDERSON, C., MIRRAHIMI, S., TURANOVA, O., AND DUMONT, T. Non-local competition slows down front acceleration during dispersal evolution. *Ann. Henri Lebesgue* **5** (2022), 1–71.

[16] CALVEZ, V., HIVERT, H., AND YOLDAŞ, H. Concentration in Lotka-Volterra parabolic equations: an asymptotic-preserving scheme. *Numer. Math.* **154**, 1-2 (2023), 103–153.

[17] CRANDALL, M. G., ISHII, H., AND LIONS, P.-L. User's guide to viscosity solutions of second order partial differential equations. *Bull. Am. Math. Soc., New Ser.* **27**, 1 (1992), 1–67.

[18] CRANDALL, M. G., AND LIONS, P.-L. Two approximations of solutions of Hamilton-Jacobi equations. *Math. Comput.* **43** (1984), 1–19.
<figure><img src="image_17.png" /><figcaption>FIGURE 13. Evolution of the solution to (S⁰) in the (x, v)-phase-space.</figcaption></figure>
[19] DEGOND, P., GOUDON, T., AND POUPAUD, F. Diffusion limit for non homogeneous and non-micro-reversible processes. *Indiana Univ. Math. J.* (2000).

[20] DIEKMANN, O., JABIN, P.-E., MISCHLER, S., AND PERTHAME, B. The dynamics of adaptation: An illuminating example and a Hamilton–Jacobi approach. *Theor. Popul. Biol.* 67, 4 (2005), 257–271.

[21] DIMARCO, G., AND PARESCHI, L. Numerical methods for kinetic equations. *Acta Numer.* 23 (2014), 369–520.

[22] DOLBEAULT, J., MOUHOT, C., AND SCHMEISER, C. Hypocoercivity for linear kinetic equations conserving mass. *Trans. Am. Math. Soc.* 367, 6 (2015), 3807–3828.

[23] DUDERSTADT, J. J., AND MARTIN, W. R. Transport theory. A Wiley-Interscience Publication. New York etc.: John Wiley & Sons. X, 613 p. £ 24.00 (1979), 1979.

[24] EVANS, L. C., AND SOUGANIDIS, P. E. Differential games and representation formulas for solutions of Hamilton-Jacobi-Isaacs equations. *Indiana Univ. Math. J.* 33 (1984), 773–797.

[25] EVANS, L. C., AND SOUGANIDIS, P. E. A PDE approach to geometric optics for certain semilinear parabolic equations. *Indiana Univ. Math. J.* 38, 1 (1989), 141–172.

[26] FATHI, A. Weak KAM from a PDE point of view: viscosity solutions of the Hamilton-Jacobi equation and Aubry set. *Proc. R. Soc. Edinb., Sect. A, Math.* 142, 6 (2012), 1193–1236.

[27] FREIDLIN, M. I. Geometric optics approach to reaction-diffusion equations. *SIAM J. Appl. Math.* 46 (1986), 222–232.

[28] GAUDEUL, B., AND HIVERT, H. Numerical approximation of a class of constrained Hamilton-Jacobi equations. Preprint, arXiv:2403.12557 (2024), 2024.

[29] HIVERT, H. A first-order asymptotic preserving scheme for front propagation in a one-dimensional kinetic reaction-transport equation. *J. Comput. Phys.* 367 (2018), 253–278.

[30] JIN, S. Efficient asymptotic-preserving (AP) schemes for some multiscale kinetic equations. *SIAM J. Sci. Comput.* 21, 2 (1999), 441–454.

[31] JIN, S. Asymptotic-preserving schemes for multiscale physical problems. *Acta Numer.* 31 (2022), 415–489.

[32] KLAR, A. A numerical method for kinetic semiconductor equations in the drift-diffusion limit. *SIAM J. Sci. Comput.* 20, 5 (1999), 1696–1712.

[33] LEMOU, M., AND MIEUSSENS, L. A new asymptotic preserving scheme based on micro-macro formulation for linear kinetic equations in the diffusion limit. *SIAM J. Sci. Comput.* 31, 1 (2008), 334–368.

[34] LIU, J.-G., AND MIEUSSENS, L. Analysis of an asymptotic preserving scheme for linear kinetic equations in the diffusion limit. *SIAM J. Numer. Anal.* 48, 4 (2010), 1474–1491.

[35] LORZ, A., MIRRAHIMI, S., AND PERTHAME, B. Dirac mass dynamics in multidimensional nonlocal parabolic equations. *Commun. Partial Differ. Equations* 36, 4-6 (2011), 1071–1098.

[36] SAINT-RAYMOND, L. *Hydrodynamic limits of the Boltzmann equation*, vol. 1971 of *Lect. Notes Math*. Berlin: Springer, 2009.

[37] SOUGANIDIS, P. E. Approximation schemes for viscosity solutions of Hamilton-Jacobi equations. *J. Differ. Equations* 59 (1985), 1–43.

(Hélène Hivert) **UNIV. RENNES, INRIA, GÉOSCIENCES RENNES - UMR 6118, IRMAR - UMR 6625, F-35000 RENNES, FRANCE**

(Tino Laidin) **UNIV BREST, CNRS UMR 6205, LABORATOIRE DE MATHÉMATIQUES DE BRETAGNE ATLANTIQUE, F-29200 BREST, FRANCE**

Email address: tino.laidin@univ-brest.fr

Email address: helene.hivert@inria.fr