# Convergence of the incremental projection method using conforming approximations

Robert Eymard, David Maltese

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

## To cite this version:

Robert Eymard, David Maltese. Convergence of the incremental projection method using conforming approximations. Computational Methods in Applied Mathematics, 2023, 10.1515/cmam-2023-0038. hal-03982033v2

HAL Id: hal-03982033

https://cnrs.hal.science/hal-03982033v2

Submitted on 28 Jun 2023

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.

Public Domain
# Convergence of the incremental projection method using conforming approximations

R. Eymard and D. Maltese*

**Keywords:** Incompressible Navier-Stokes equations, incremental projection scheme, conforming scheme, convergence analysis

## Abstract

We prove the convergence of an incremental projection numerical scheme for the time-dependent incompressible Navier-Stokes equations, without any regularity assumption on the weak solution. The velocity and the pressure are discretized in conforming spaces, whose compatibility is ensured by the existence of an interpolator for regular functions which preserves approximate divergence free properties. Owing to a priori estimates, we get the existence and uniqueness of the discrete approximation. Compactness properties are then proved, relying on a Lions-like lemma for time translate estimates. It is then possible to show the convergence of the approximate solution to a weak solution of the problem. The construction of the interpolator is detailed in the case of the lowest degree Taylor-Hood finite element.

## 1 Introduction

The Navier-Stokes equations for a homogeneous incompressible fluid can be written in a strong form as:

$$ \partial_t u + (u \cdot \nabla)u - \Delta u + \nabla p = f \quad \text{in } (0, T) \times \Omega, \quad (1a) $$

$$ \operatorname{div} u = 0 \quad \text{in } (0, T) \times \Omega, \quad (1b) $$

where the density and the viscosity are set to one for the sake of simplicity, and where

$$ T > 0, \text{ and } \Omega \text{ is a connected, open and bounded subset of } \mathbb{R}^d, d \in \{2, 3\}, \\ \text{with a Lipschitz boundary } \partial \Omega. \quad (2) $$

The variables $u$ and $p$ are respectively the velocity and the pressure of the fluid, and (1a) and (1b) respectively model the momentum conservation and the mass conservation of an incompressible fluid. This system is supplemented with the homogeneous Dirichlet boundary condition

$$ u = 0 \quad \text{on } (0, T) \times \partial \Omega, \quad (3) $$

and the initial condition

$$ u(0) = u_0 \quad \text{in } \Omega. \quad (4) $$

The function $u_0$ is the initial datum for the velocity and the function $f$ is the source term.

While this system of equations is coupling the velocity and the pressure, projection numerical schemes, introduced in [9] and [37], enable the successive resolution of decoupled elliptic equations for the velocity

<footnote>*Université Gustave Eiffel, LAMA, (UMR 8050), UPEM, UPEC, CNRS, F-77454, Marne-la-Vallée (France), robert.eymard, david.maltese@univ-eiffel.fr</footnote>
and the pressure. This leads to cheaper and smaller computations than those issued from a coupled approximation. Nevertheless, the projection method as originally presented has some drawbacks: the precision is limited for the $L^2$ norm of the pressure and for the $H^1$ norm of the velocity. This limited precision is a consequence of the artificial boundary condition for the pressure. The use of a higher-order time discretization does not improve the convergence rate. We refer to [31] for error estimates for the semi-discrete case (see also [33] and [34]). Different choices of the pressure boundary condition have been discussed to improve the efficiency of this method. The first major improvement of the original projection scheme was proposed in [19] (see also [38]) and now referred to (following [22]) as incremental projection schemes. Such schemes can be written as follows in the case of a continuous-in-space discrete-in-time scheme:

* **Prediction step**: Find a function $\tilde{u}^{n+1} \in H_0^1(\Omega)^d$, which is therefore regular in space and respects the boundary conditions but is not divergence-free, such that the following linearized momentum equation holds in the homogeneous Dirichlet weak sense:

$$ \frac{1}{\delta t}(\tilde{u}^{n+1} - u^n) + $$

$$ (\tilde{u}^n \cdot \nabla)\tilde{u}^{n+1} + \frac{1}{2}\mathrm{div}\tilde{u}^n \tilde{u}^{n+1} - \Delta\tilde{u}^{n+1} + \nabla p^n = f^{n+1} \quad \text{in } \Omega. \quad (5) $$

This step involves $d$ decoupled resolution of elliptic problems for each of the components of the velocity.

* **Correction step**: Find the new pressure field $p^{n+1} \in H^1(\Omega) \cap L_0^2(\Omega)$ (denoting by $L_0^2(\Omega)$ the space of functions of $L^2(\Omega)$ with null average on $\Omega$), and a corrected velocity $u^{n+1} \in V(\Omega)$ (denoting by $V(\Omega)$ the space of $L^2$-divergence-free functions and null normal trace, precisely defined by (11)) therefore with a lower regularity in space and a weaker boundary condition, such that, in the homogeneous Neumann weak sense for the unknown $p^{n+1} - p^n \in H^1(\Omega) \cap L_0^2(\Omega)$,

$$ \frac{1}{\delta t}(u^{n+1} - \tilde{u}^{n+1}) + \nabla(p^{n+1} - p^n) = 0 \quad \text{in } \Omega, \quad (6a) $$

$$ \mathrm{div}u^{n+1} = 0 \quad \text{in } \Omega. \quad (6b) $$

This step involves the resolution of an elliptic problem for finding the new pressure.

Scheme (5)-(6) is called an “incremental projection scheme”, since it is obtained by introducing the previous pressure gradient in the prediction step and by solving the increment of the pressure in the correction step, which is a projection step of the predicted velocity on the divergence-free functions. Such a scheme seems to be much more efficient from a computational point of view than non-incremental projection schemes and has been the object of several error analyses, under some regularity assumptions on the solution of the continuous problem in the semi-discrete setting, see [33] and [34]. The incremental schemes have been the object of some works in the fully discrete setting (see [25], [24], [20], [23]). For more details concerning projection methods we refer to [22] and references therein.

Some recent papers [27, 18] propose convergence proofs for fully discrete projection schemes. In [18], the incremental projection scheme is considered without any assumption of regularity assumptions on the exact solution. The proof of its convergence is done both for the semi-discrete scheme and for a fully discrete scheme using the Marker-And-Cell scheme (introduced in the seminal paper [26]). In this paper, we extend such a proof to the case where the spatial discretization is done using conforming methods (in most of the cases, finite element methods or spectral methods).

For any integer $N \ge 1$, we define the time step by $\delta t_N = T/N$ and we consider approximation spaces $X_N \subset H_0^1(\Omega)^d$ for the predicted velocity, $M_N \subset H^1(\Omega) \cap L_0^2(\Omega)$ for the approximate pressure. Then the space $V_N$, defined as the space of functions which are $L^2$-orthogonal to the gradient of the elements of
$M_N$ (it is therefore not a subset of $V(\Omega)$) is used for the approximation of the corrected velocity. The fully discrete incremental projection scheme is given in Section 2, and is shown to have a unique solution $(\tilde{u}_N^n, p_N^n, u_N^n)_{n=1,...,N} \in (X_N \times M_N \times V_N)^N$.

In order to prove the convergence of the incremental projection scheme, we need to prove the compactness of the sequence of the approximate velocities in $L^2$. The main difficulty relies in the time translate estimates; this difficulty is solved in the coupled case by the famous Lions' Lemma [29, Lemme 5.1 p.58], bounding the norm in $L^2$ by a combination of the norm in $H^1$ and of a semi-norm similar to a $H^{-1}$ norm. In the case of the incremental projection method, the difficulty relies in the fact that there are two different approximations for the velocities. We handle this difficulty in Section 3.3, by closely following the method first introduced in [18]. We can then conclude the convergence proof to a solution, considered in the weak sense given by Definition 1.1, in Section 3.4.

It is noticeable that, for this existence and uniqueness result, no compatibility condition between $M_N$ and $X_N$ is required. This is no longer the case in Section 3, where the convergence of the method is proved. The compatibility condition is expressed through the existence of an interpolator $\Pi_N : W \to X_N \cap V_N \cap L^\infty(\Omega)^d$, whose purpose is to build approximate values in the discrete velocity spaces of some regular functions. This interpolator replaces in this paper a Fortin interpolator (see [32, 12, 13, 17] for more properties on such operator). In the present paper, the requested properties are slightly different (see Section 3.1):

* The set $W$ is a dense subset of the set $H_0^1(\Omega)^d \cap V(\Omega)$ for the $H_0^1(\Omega)^d$ norm (see (21)).

* For any $\varphi \in W$, the convergence of $\Pi_N \varphi$ to $\varphi$ is assumed to hold in $H_0^1(\Omega)^d$, while the $L^\infty$ norm of $\Pi_N \varphi$ remains bounded (see (24)).

So the hypotheses which are given on the space discretization in this paper are not shown to be equivalent to an inf-sup property. In Section 4, we consider the example of the lowest degree Taylor-Hood finite element. In this section we provide all the calculations that prove the convergence property which is requested on $\Pi_N$, since the properties which are studied in the literature are generally different. We show that constructing $\Pi_N$ and checking its properties is much simplified by the fact that we only need to apply it on regular divergence free functions with compact support. It is in some way simpler than the verification of the inf-sup condition for these spaces (we provide in this section simple examples where there is no inf-sup stability with the lowest degree Taylor-Hood finite element, but where $\Pi_N$ is nevertheless defined). Note that the study done in this paper also applies to the semi-discrete scheme (see Remark 3.2). We refer to [18] for the convergence of the semi-discrete scheme using another discretization of the convective term.

Throughout the paper, we shall assume that the data $f$ and $u_0$ satisfy

$$f \in L^2(0, T; H^{-1}(\Omega)^d) \text{ and } u_0 \in L^2(\Omega)^d. \qquad (7)$$

We denote in the whole paper

$$\forall v, w \in L^2(\Omega)^d, (v, w) := \int_{\Omega} v(x) \cdot w(x) \, dx \qquad (8)$$

and

$$\forall g \in H^{-1}(\Omega)^d, \forall v \in H_0^1(\Omega)^d, \langle g, v \rangle := g(v). \qquad (9)$$

We recall that we denote by $L_0^2(\Omega)$ the space of $L^2$ functions with null average on $\Omega$, which can be defined by

$$L_0^2(\Omega) = \{q \in L^2(\Omega) \text{ such that } \int_{\Omega} q(x) \, dx = 0\}, \qquad (10)$$

and by $V(\Omega)$ the space of $L^2$-divergence-free functions, which can be defined by

$$V(\Omega) = \{u \in L^2(\Omega)^d \text{ such that } (u, \nabla \xi) = 0 \text{ for any } \xi \in H^1(\Omega)\}. \qquad (11)$$
Recall that $V(\Omega)$ (sometimes also denoted by $L^2_\sigma(\Omega)$) is the closure in $L^2(\Omega)$ of $C_{c,\sigma}^\infty(\Omega) := \{v \in C_c^\infty(\Omega)^d; \text{div}v = 0\}$ and that we can write as well [35, Eq. 3.5.13], [4, Definition IV.3.2 p. 248]

$$V(\Omega) = \{u \in H_{\text{div}}(\Omega) \text{ such that } \text{div}u = 0 \text{ and } \gamma_{\nu}u = 0\},$$

denoting by $\gamma_\nu$ the normal trace of $u$ on $\partial\Omega$. Then $H_0^1(\Omega)^d \cap V(\Omega)$ (sometimes also denoted by $W_{0,\sigma}^{1,2}(\Omega)$) is the subset of $H_0^1(\Omega)^d$ of divergence-free functions, and is the closure in $H^1(\Omega)$ of $C_{c,\sigma}^\infty(\Omega)$ [35, p.157], [4, Lemma IV.3.4 p.249]. Let us define the weak solutions of Problem (1)-(4) in the sense of Leray [28].

**Definition 1.1 (Weak solution):** Under the assumptions (2) and (7), a function $u \in L^2(0, T; H_0^1(\Omega)^d \cap V(\Omega)) \cap L^\infty(0, T; L^2(\Omega)^d)$ is a weak solution of the problem (1)-(4) if

$$- \int_0^T (u, \partial_t v) \, dt + \int_0^T ((u \cdot \nabla)u, v) \, dt + \int_0^T \int_\Omega \nabla u : \nabla v \, dx \, dt \\ = (u_0, v(0, \cdot)) + \int_0^T \langle f, v \rangle \, dt \quad (12)$$

for any $v \in C_c^\infty(\Omega \times [0, T))^d$, such that $\text{div}v = 0$ a.e. in $\Omega \times (0, T)$.

# 2 The projection scheme using a conforming method

## 2.1 Space and time discretizations

We consider a partition of the time interval $[0, T]$, which we suppose uniform to alleviate the notations, so that the assumptions read (omitting in the whole paper to recall that $N$ is assumed to be an integer):

$$N \ge 1, \quad \partial_t N = \frac{T}{N}, \quad t_N^n = n \partial_t N \text{ for } n \in [0, N]. \quad (13)$$

We consider a sequence of velocity-pressure approximations in $H_0^1(\Omega)^d$ and $H^1(\Omega) \cap L_0^2(\Omega)$ respectively. For the approximation of the predicted velocity, let

$$(X_N)_{N \ge 1} \text{ be a sequence of closed subspaces of } H_0^1(\Omega)^d. \quad (14)$$

For the approximation of the pressure, let

$$(M_N)_{N \ge 1} \text{ be a sequence of closed subspaces of } H^1(\Omega) \cap L_0^2(\Omega). \quad (15)$$

For the approximation of the corrected velocity we define the space of weakly divergence free functions by

$$V_N = \{v \in L^2(\Omega)^d \text{ such that } (v, \nabla q) = 0 \text{ for any } q \in M_N\}. \quad (16)$$

We denote by $\mathcal{P}_{V_N} : L^2(\Omega)^d \to V_N$ the orthogonal projection in $L^2(\Omega)^d$ onto the space $V_N$. We denote by $\mathcal{P}_{V(\Omega)} : L^2(\Omega)^d \to V(\Omega)$ the orthogonal projection in $L^2(\Omega)^d$ onto the space $V(\Omega)$.

**Remark 2.1:** The choice $X_N = H_0^1(\Omega)^d$, $M_N = H^1(\Omega) \cap L_0^2(\Omega)$ (which implies $V_N = V(\Omega)$) yields the continuous-in-space discrete-in-time scheme.
## 2.2 The projection scheme

Let $a: H_0^1(\Omega)^d \times H_0^1(\Omega)^d \to \mathbb{R}$ be the coercive bilinear form defined by

$$a(\mathbf{u}, \mathbf{v}) = \int_{\Omega} \nabla \mathbf{u} : \nabla \mathbf{v} \, d\mathbf{x},$$

and let $b: H_0^1(\Omega)^d \times H_0^1(\Omega)^d \times H_0^1(\Omega)^d \to \mathbb{R}$ be the “skew-symmetric” trilinear form defined by

$$b(\mathbf{u}, \mathbf{v}, \mathbf{w}) = \frac{1}{2} \int_{\Omega} \left( (\mathbf{u} \cdot \nabla)\mathbf{v} \cdot \mathbf{w} - (\mathbf{u} \cdot \nabla)\mathbf{w} \cdot \mathbf{v} \right) d\mathbf{x}.$$

Note that we have for any $(\mathbf{u}, \mathbf{v}, \mathbf{w}) \in H_0^1(\Omega)^d \times H_0^1(\Omega)^d \times H_0^1(\Omega)^d$

$$b(\mathbf{u}, \mathbf{v}, \mathbf{w}) = \int_{\Omega} (\mathbf{u} \cdot \nabla)\mathbf{v} \cdot \mathbf{w} \, d\mathbf{x} + \frac{1}{2} \int_{\Omega} \mathrm{div}\mathbf{u}(\mathbf{v} \cdot \mathbf{w}) \, d\mathbf{x}. \quad (17)$$

Note that the fact that $b$ is continuous over $H_0^1(\Omega)^d \times H_0^1(\Omega)^d \times H_0^1(\Omega)^d$ is a consequence of the Sobolev continuous embedding of $H_0^1(\Omega)$ in $L^4(\Omega)$ for $d \in \{2, 3\}$. Using the above definitions, the first order time fully discrete incremental projection scheme reads:

Initialization:

$$\text{Let } \mathbf{u}_N^0 = \mathcal{P}_{V_N} \mathbf{u}_0, \quad \tilde{\mathbf{u}}_N^0 = 0 \text{ and } p_N^0 = 0. \qquad (18a)$$

Solve for $0 \le n \le N - 1$:

Prediction step: Find $\mathbf{u}_N^{n+1} \in \mathbf{X}_N$ such that

$$\frac{1}{\partial_N} \left( (\tilde{\mathbf{u}}_N^{n+1}, \mathbf{v}) - (\mathbf{u}_N^n, \mathbf{v}) \right) + b(\tilde{\mathbf{u}}_N^n, \tilde{\mathbf{u}}_N^{n+1}, \mathbf{v}) \qquad (18b)$$

$$+ a(\tilde{\mathbf{u}}_N^{n+1}, \mathbf{v}) + (\nabla p_N^n, \mathbf{v}) = \langle \mathbf{f}_N^{n+1}, \mathbf{v} \rangle, \text{ for any } \mathbf{v} \in \mathbf{X}_N,$$

Correction step 1: Find $p_N^{n+1} \in M_N$ such that

$$(\nabla(p_N^{n+1} - p_N^n), \nabla q) = \frac{1}{\partial_N} (\tilde{\mathbf{u}}_N^{n+1}, \nabla q), \text{ for any } q \in M_N(\Omega). \qquad (18c)$$

Correction step 2: Define $\mathbf{u}_N^{n+1} \in L^2(\Omega)^d$ as follows

$$\mathbf{u}_N^{n+1} = \tilde{\mathbf{u}}_N^{n+1} - \partial_N \nabla (p_N^{n+1} - p_N^n). \qquad (18d)$$

The term $\mathbf{f}_N^{n+1} \in H^{-1}(\Omega)^d$ is given by $\mathbf{f}_N^{n+1} = \frac{1}{\partial_N} \int_{t^n}^{t^{n+1}} \mathbf{f}(t) \, dt$.

**Remark 2.2:** Note that $\mathbf{u}_N^{n+1} \in \mathbf{X}_N + \nabla M_N \subset L^2(\Omega)^d$. We refer to [20] for other discretizations of the projection step. Multiplying (18d) by the gradient of an element of $M_N$ and using (18d) we obtain for any $n \in [0, N-1]$

$$(\mathbf{u}_N^{n+1}, \nabla q) = 0 \text{ for any } q \in M_N, \qquad (19)$$

which means that $\mathbf{u}_N^{n+1} \in V_N$. Hence, as stated in the introduction, $\mathbf{u}_N$ satisfies the discrete divergence free constraint, but not the regularity one. We therefore obtain a discrete Helmholtz-Leray decomposition of the function $\tilde{\mathbf{u}}_N^{n+1}$ by writing $\tilde{\mathbf{u}}_N^{n+1} = \mathbf{u}_N^{n+1} + \partial_N \nabla (p_N^{n+1} - p_N^n)$. Although $\tilde{\mathbf{u}}_N$ satisfies the regularity constraint and the boundary conditions, but not the discrete divergence free constraint, Theorem 3.8 proves that both approximate fields converge to the same weak solution of the continuous problem, simultaneously satisfying all these constraints.
**Remark 2.3:** Note that the numerical method can be expressed without using the sequence of corrected velocity $(\tilde{u}_N^n)_{n \in [0, N]}$. Letting $p_N^0 = 0$, the predicted velocity is such that

$$ \frac{1}{\partial_N} \left( (\tilde{u}_N^1, v) - (P_{V_N} u_0, v) \right) + a(\tilde{u}_N^1, v) = \langle \tilde{f}_N^1, v \rangle, \text{ for any } v \in \mathbf{X}_N, $$

and, for $n \in [1, N-1]$,

$$ \frac{1}{\partial_N} \left( (\tilde{u}_N^{n+1}, v) - (\tilde{u}_N^n, v) \right) + b(\tilde{u}_N^n, \tilde{u}_N^{n+1}, v) \\ + a(\tilde{u}_N^{n+1}, v) + (\nabla(2p_N^n - p_N^{n-1}), v) = \langle \tilde{f}_N^{n+1}, v \rangle, \text{ for any } v \in \mathbf{X}_N. $$

The pressure $p_N^{n+1} \in M_N$ is such that

$$ (\nabla(p_N^{n+1} - p_N^n), \nabla q) = \frac{1}{\partial_N} (\tilde{u}_N^{n+1}, \nabla q) \text{ for any } q \in M_N. $$

**Remark 2.4:** The choice for the initial conditions $\tilde{u}_N^0$ and $p_N^0$ plays no significant role for the convergence study. It suffices to consider any sequences $(\tilde{u}_N^0)_{N \ge 1}$ bounded in $H^1(\Omega)^d$ and $(\tilde{p}_N^0)_{N \ge 1}$ bounded in $H^1(\Omega)$.

## 2.3 Existence of a solution to the projection scheme

The following existence result allows to define the approximate solutions obtained by the projection scheme (18).

**Lemma 2.1 (Approximate solutions):** Under Assumptions (2), (7), (13), (14), (15), for any $N \ge 1$, there exists a unique sequence $(\tilde{u}_N^n, u_N^n, p_N^n)_{n \in [1, N]} \subset \mathbf{X}_N \times \mathbf{V}_N \times M_N$ satisfying (18). We then define the functions $u_N: (0, T) \to \mathbf{V}_N$ and $\tilde{u}_N: (0, T) \to \mathbf{X}_N$ by

$$ u_N(t) = \sum_{n=0}^{N-1} 1_{(t^n_N, t^{n+1}_N]}(t)u_N^n, \quad \tilde{u}_N(t) = \sum_{n=0}^{N-1} 1_{(t^n_N, t^{n+1}_N]}(t)\tilde{u}_N^{n+1}. \quad (20) $$

**PROOF.** For given functions $p_N^n \in M_N$, $\tilde{u}_N^n \in \mathbf{X}_N$ and $u_N^n \in \mathbf{V}_N$, Problem (18b) is under the form:

$$ \tilde{u}_N^{n+1} \in \mathbf{X}_N, \quad \forall v \in \mathbf{X}_N, \quad a_N^n(\tilde{u}_N^{n+1}, v) = \langle g_N^n, v \rangle, $$

where we define the bilinear form $a_N^n: H_0^1(\Omega)^d \times H_0^1(\Omega)^d \to \mathbb{R}$, $(u, v) \mapsto \frac{1}{\partial_N}(u, v) + b(\tilde{u}_N^n, u, v) + a(u, v)$, and where $g_N^n \in H^{-1}(\Omega)^d$. The continuity of this bilinear form is a consequence of that of b. Owing to the property $b(\cdot, v, v) = 0$, this continuous bilinear form is such that

$$ \forall v \in \mathbf{X}_N, \quad a_N^n(v, v) \ge \|v\|_{H_0^1(\Omega)^d}^2. $$

Hence the existence and uniqueness of the solution $\tilde{u}_N^{n+1} \in \mathbf{X}_N$ (which is closed by assumption) is given by the Lax-Milgram theorem.

For given functions $p_N^n \in M_N$, $\tilde{u}_N^{n+1} \in \mathbf{X}_N$, Problem (18c) is under the form:

$$ p_N^{n+1} \in M_N, \quad \forall q \in M_N, \quad (\nabla p_N^{n+1}, \nabla q) = c_N(q). $$

Hence the existence and uniqueness of the solution $p_N^{n+1} \in M_N$ (again closed by assumption) is given by the fact that the bilinear form $(p, q) \mapsto (\nabla p, \nabla q)$ is coercive on $H^1(\Omega) \cap L_0^2(\Omega)$. $\square$
# 3 Convergence study

When the Stokes problem is formulated as a coupled mixed problem whose unknowns are the velocity and the pressure, the approximation spaces used in the case of a coupled scheme must satisfy a compatibility condition (called the inf-sup condition) to obtain existence and uniqueness of a solution (see [6]). Although such condition is not needed for proving Lemma 2.1, the numerical simulations presented in [24] suggest that some compatibility conditions are necessary to obtain satisfactory convergence properties. We present in Section 3.1 some sufficient compatibility conditions, which are not equivalent to the standard inf-sup condition, as detailed in Section 4.

## 3.1 Convergence assumptions

The following assumptions are done on the sequences $(X_N)_{N \ge 1}$ and $(M_N)_{N \ge 1}$ in order to prove the convergence of the projection scheme as $N \to \infty$. These assumptions, provided in terms of the existence of a family $(\Pi_N)_{N \ge 1}$ of interpolators on some spaces, replace the standard uniform inf-sup hypothesis. We assume that

$$W \text{ is a dense subset of } V(\Omega) \cap H_0^1(\Omega)^d \text{ for the norm of } H_0^1(\Omega)^d. \qquad (21)$$

We define

$$\forall N \ge 1, \quad E_N = X_N \cap V_N \cap L^\infty(\Omega)^d, \qquad (22)$$

$$\forall N \ge 1, \forall \mathbf{u} \in E_N, \quad \|\mathbf{u}\|_{E} = \|\nabla \mathbf{u}\|_{L^2(\Omega)^{d \times d}} + \|\mathbf{u}\|_{L^\infty(\Omega)^d}. \qquad (23)$$

We assume that

$$\text{For any } N \ge 1 \text{ there exists a mapping } \Pi_N : W \to E_N \text{ such that,} \\ \text{for any } \varphi \in W, (\Pi_N \varphi)_{N \ge 1} \text{ converges to } \varphi \text{ in } H_0^1(\Omega)^d \\ \text{and the sequence } (\|\Pi_N \varphi\|_E)_{N \ge 1} \text{ is bounded.} \qquad (24)$$

The sequence $(M_N)_{N \ge 1}$ is assumed to be such that

$$(\mathcal{P}_{M_N}(p))_{N \ge 1} \text{ converges to } p \text{ in } H^1(\Omega), \\ \text{for any } p \in H^1(\Omega) \cap L_0^2(\Omega), \qquad (25)$$

where we denote by $\mathcal{P}_{M_N} : H^1(\Omega) \cap L_0^2(\Omega) \to M_N$ the orthogonal projection onto the space $M_N$ in the Hilbert space $H^1(\Omega) \cap L_0^2(\Omega)$ equipped with the scalar product $(p, q) \mapsto (\nabla p, \nabla q)$.

**Remark 3.1:** The density of $C_c^\infty(\Omega)^d \cap V(\Omega)$ in $V(\Omega)$ for the $L^2(\Omega)^d$ norm is proved in particular in [4, Lemma IV.3.5 p. 249]. Since the density of $W$ in $H_0^1(\Omega)^d \cap V(\Omega)$ enables to approximate any element of $C_c^\infty(\Omega)^d \cap V(\Omega)$ as closely as desired for the $L^2(\Omega)^d$ norm, we get that $W$ is dense as well in $V(\Omega)$ for the $L^2(\Omega)^d$ norm.

**Remark 3.2:** If, for all $N \ge 1$, we consider as in Remark 2.1 the semi-discrete case $X_N = H_0^1(\Omega)^d$, $M_N = H^1(\Omega) \cap L_0^2(\Omega)$ (recall that this yields $V_N = V(\Omega)$), we let $W = V(\Omega) \cap C_c^\infty(\Omega)^d$ and $\Pi_N = \text{Id}$. Then all the assumptions of this section are satisfied.

**Remark 3.3:** We provide in Section 4 the construction of $\Pi_N$ in the case of the lowest degree Taylor-Hood finite element. We show in Theorem 4.14 that it satisfies the assumptions given in this section, for a regular family of meshes in the standard sense.
## 3.2 Space estimate

**Lemma 3.1:** Under Assumptions (2), (7), (13), (14), (15), for any $N \ge 1$, the following relation holds for $n \in [0, N-1]$:

$$ \frac{1}{2\partial_N} \left( \|\mathbf{u}_N^{n+1}\|_{L^2(\Omega)^d}^2 - \|\mathbf{u}_N^n\|_{L^2(\Omega)^d}^2 \right) + \frac{\partial_N}{2} \left( \|\nabla p_N^{n+1}\|_{L^2(\Omega)^d}^2 - \|\nabla p_N^n\|_{L^2(\Omega)^d}^2 \right) \\ + \frac{1}{2\partial_N} \|\tilde{\mathbf{u}}_N^{n+1} - \mathbf{u}_N^n\|_{L^2(\Omega)^d}^2 + \|\tilde{\mathbf{u}}_N^{n+1}\|_{H_0^1(\Omega)^d}^2 = \langle \mathbf{f}_N^{n+1}, \tilde{\mathbf{u}}_N^{n+1} \rangle. \quad (26) $$

Consequently, there exists $C_1$ depending only on $|\Omega|$, on the $L^2(\Omega)^d$-norm of $\mathbf{u}_0$ and on the $L^2(0, T; H^{-1}(\Omega)^d)$-norm of $\mathbf{f}$ such that the functions $\tilde{\mathbf{u}}_N$ and $\mathbf{u}_N$ defined by (20) satisfy:

$$ \|\tilde{\mathbf{u}}_N\|_{L^2(0,T;H_0^1(\Omega)^d)} \le C_1 \quad \text{and} \quad \|\mathbf{u}_N\|_{L^\infty(0,T;L^2(\Omega)^d)} \le C_1. \quad (27) $$

$$ \|\tilde{\mathbf{u}}_N\|_{L^\infty(0,T;L^2(\Omega)^d)} \le C_1 \quad \text{and} \quad \|\mathbf{u}_N - \tilde{\mathbf{u}}_N\|_{L^2(0,T;L^2(\Omega)^d)} \le C_1 \sqrt{\partial_N}. \quad (28) $$

**PROOF.** We take $\tilde{\mathbf{u}}_N^{n+1}$ in (18b) and we obtain for $n \in [0, N-1]$

$$ \frac{1}{2\partial_N} \|\tilde{\mathbf{u}}_N^{n+1}\|_{L^2(\Omega)^d}^2 - \frac{1}{2\partial_N} \|\mathbf{u}_N^n\|_{L^2(\Omega)^d}^2 + \frac{1}{2\partial_N} \|\tilde{\mathbf{u}}_N^{n+1} - \mathbf{u}_N^n\|_{L^2(\Omega)^d}^2 \\ + (\nabla p_N^n, \tilde{\mathbf{u}}_N^{n+1}) + \|\tilde{\mathbf{u}}_N^{n+1}\|_{H_0^1(\Omega)^d}^2 = \langle \mathbf{f}_N^{n+1}, \tilde{\mathbf{u}}_N^{n+1} \rangle. \quad (29) $$

Squaring the relation $\mathbf{u}_N^{n+1} + \partial_N \nabla p_N^{n+1} = \tilde{\mathbf{u}}_N^{n+1} + \partial_N \nabla p_N^n$, integrating over $\Omega$, and owing to $\mathbf{u}_N^{n+1} \in V_N$, we get that for $n \in [0, N-1]$,

$$ \frac{1}{2\partial_N} \|\mathbf{u}_N^{n+1}\|_{L^2(\Omega)^d}^2 + \frac{\partial_N}{2} \|\nabla p_N^{n+1}\|_{L^2(\Omega)^d}^2 = \frac{1}{2\partial_N} \|\tilde{\mathbf{u}}_N^{n+1}\|_{L^2(\Omega)^d}^2 \\ + \frac{\partial_N}{2} \|\nabla p_N^n\|_{L^2(\Omega)^d}^2 + (\tilde{\mathbf{u}}_N^{n+1}, \nabla p_N^n). $$

Summing the latter relation with (29) yields (26) for $n \in [0, N-1]$. We then get Relations (27) by summing over the time steps, using the Cauchy-Schwarz and Poincaré inequalities. $\square$

**Lemma 3.2:** Under Assumptions (2), (7), let $(\partial_N, \mathbf{X}_N, M_N)_{N \ge 1}$ be a sequence such that (13), (14), (15) and (25) hold. Then there exists $\mathbf{u} \in L^\infty(0, T; L^2(\Omega)^d) \cap L^2(0, T; H_0^1(\Omega)^d \cap V(\Omega))$ and a subsequence again denoted by $(\partial_N, \mathbf{X}_N, M_N)_{N \ge 1}$ such that, letting $(\tilde{\mathbf{u}}_N, \mathbf{u}_N)_{N \ge 1}$ be defined by Lemma 2.1, then $(\tilde{\mathbf{u}}_N)_{N \ge 1}$ weakly converges to $\mathbf{u}$ in $L^2(0, T; H_0^1(\Omega)^d)$ and $(\mathbf{u}_N)_{N \ge 1}$ weakly-\* converges to $\mathbf{u}$ in $L^\infty(0, T; L^2(\Omega)^d)$.

**PROOF.** Owing to (27), we get the existence of $\tilde{\mathbf{u}} \in L^2(0, T; H_0^1(\Omega)^d)$ such that, up to the extraction of a subsequence, the sequence $(\tilde{\mathbf{u}}_N)_{N \ge 1}$ weakly converges to $\tilde{\mathbf{u}}$ in $L^2(0, T; H_0^1(\Omega)^d)$, and of $\mathbf{u} \in L^\infty(0, T; L^2(\Omega)^d)$ such that $(\mathbf{u}_N)_{N \ge 1}$ weakly converges to $\mathbf{u}$ in $L^\infty(0, T; L^2(\Omega)^d)$ for the weak star topology. Let $\xi \in H^1(\Omega) \cap L_0^2(\Omega)$ and $\varphi \in C_c^\infty([0, T])$ be given.

We notice that (19) yields

$$ \int_0^T (\mathbf{u}_N, \nabla P_{M_N} \xi) \varphi(t) \, dt = 0, \text{ for any } N \ge 1. $$

Using the convergence of the sequence $(P_{M_N} \xi)_{N \ge 1}$ in $H^1(\Omega)$, the weak convergence of the sequence $(\mathbf{u}_N)_{N \ge 1}$ in $L^2((0, T) \times \Omega)^d$ and passing to the limit in the previous identity gives

$$ \int_0^T (\mathbf{u}, \nabla \xi) \varphi(t) \, dt = 0, \text{ for any } \xi \in H^1(\Omega). $$

We then obtain that $\mathbf{u} \in L^2(0, T; H_0^1(\Omega)^d \cap V(\Omega))$. Using (28), we get that $\tilde{\mathbf{u}} = \mathbf{u}$, which concludes the proof. $\square$
## 3.3 Time estimates

The weak convergence property given by Lemma 3.2 is not sufficient for passing to the limit in the scheme, owing to the presence of the nonlinear convection term. Hence we need some stronger compactness property on one of the subsequences $(\tilde{u}_N)_{N \ge 1}$ or $(\tilde{u}_N)_{N \ge 1}$, obtained through the application of Theorem 3.7 below. We therefore need to prove an estimate on the time translates which fulfills (41). This is done by first proving an estimate for the time translates with respect to a semi-norm which plays the role of a discrete $H^{-1}$ seminorm. This seminorm, denoted by $|\cdot|_{*,N}$, makes the pressure gradient in (18d) vanish by considering test functions $v \in E_N \subset V_N$. It is defined for any $w \in L^2(\Omega)^d$ by

$$|w|_{*,N} = \sup\{(w, v), v \in E_N, \|v\|_E \le 1\}. \qquad (30)$$

Recall that $E_N$ and $\|v\|_E$ are defined by (22)-(23) and that the above definition remains meaningful even if $E_N = \{0\}$.

**Lemma 3.3 (A first estimate on the time translates):** Under Assumptions (2), (7), (13), (14), (15), there exists $C_2$ only depending on $|\Omega|$, $\|u_0\|_{L^2(\Omega)^d}$ and $\|f\|_{L^2(0,T;H^{-1}(\Omega)^d)}$ such that for any $N \ge 1$ and for any $\tau \in (0,T)$,

$$\int_{0}^{T-\tau} |\tilde{u}_N(t + \tau) - \tilde{u}_N(t)|^2_{*,N} dt \le C_2 \tau (\tau + \delta_N).$$

**Remark 3.4:** The proof of Lemma 3.3 relies on the $L^2(0, T; H_0^1(\Omega)^d)$ estimates (27)-(28) proved in Lemma 3.1.

**PROOF.** Let $N \ge 2$ and $\tau \in (0, T)$ (for $N = 1$ the quantity we have to estimate is zero). Let $(\chi_{N,\tau}^n)_{n \in \llbracket 1, N-1 \rrbracket}$ be the family of measurable functions defined for $n \in \llbracket 1, N-1 \rrbracket$ and $t \in \mathbb{R}$ by $\chi_{N,\tau}^n(t) = \mathbb{1}_{(t_N^n - \tau, t_N^n]}(t)$, then

$$\tilde{u}_N(t + \tau) - \tilde{u}_N(t) = \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t)(\tilde{u}_N^{n+1} - \tilde{u}_N^n), \quad \text{for any } t \in (0, T - \tau). \qquad (31)$$

Hence, owing to (18b), we have for any $v \in X_N$ and for any $t \in (0, T - \tau)$ the following identity

$$(\tilde{u}_N(t + \tau) - \tilde{u}_N(t), v) = -\delta_N \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) a(\tilde{u}_N^{n+1}, v) \\ - \delta_N \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) b(\tilde{u}_N^n, \tilde{u}_N^{n+1}, v) - \delta_N \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) (\nabla(2p_N^n - p_N^{n-1}), v) \\ + \delta_N \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \langle f_N^{n+1}, v \rangle.$$

Let $v \in E_N$, and define $A(t) = (\tilde{u}_N(t + \tau) - \tilde{u}_N(t), v)$. Using the previous identity we obtain

$$A(t) = A_d(t) + A_c(t) + A_p(t) + A_f(t),$$

with

$$A_d(t) = - \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \delta_N \int_{\Omega} \nabla \tilde{u}_N^{n+1} : \nabla v \, dx,$$

$$A_c(t) = - \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \delta_N b(\tilde{u}_N^n, \tilde{u}_N^{n+1}, v),$$
$$
\begin{align*}
A_p(t) &= \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \partial_t N \int_{\Omega} (2p_N^n - p_N^{n-1}) \mathrm{div} \, \mathbf{v} \, \mathrm{d}\mathbf{x}, \\
A_f(t) &= \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \partial_t N \langle \mathbf{f}_N^{n+1}, \mathbf{v} \rangle.
\end{align*}
$$

Using (23) we have

$$
A_d(t) \le \|\mathbf{v}\|_{\mathcal{E}} \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \partial_t N \|\tilde{\mathbf{u}}_N^{n+1}\|_{H_0^1(\Omega)^d}. \quad (32)
$$

Using the identity (17) and the estimates (27)-(28) we have

$$
\begin{align*}
A_c(t) = & -\sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \partial_t N ((\tilde{\mathbf{u}}_N^n \cdot \nabla) \tilde{\mathbf{u}}_N^{n+1}, \mathbf{v}) \\
& -\frac{1}{2} \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \partial_t N (\operatorname{div} \tilde{\mathbf{u}}_N^n \tilde{\mathbf{u}}_N^{n+1}, \mathbf{v}) \\
\le & \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \partial_t N \|\tilde{\mathbf{u}}_N^n\|_{L^2(\Omega)^d} \|\tilde{\mathbf{u}}_N^{n+1}\|_{H_0^1(\Omega)^d} \|\mathbf{v}\|_{L^\infty(\Omega)^d} \\
& + \frac{1}{2} \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \partial_t N \|\operatorname{div} \tilde{\mathbf{u}}_N^n\|_{L^2(\Omega)^d} \|\tilde{\mathbf{u}}_N^{n+1}\|_{L^2(\Omega)^d} \|\mathbf{v}\|_{L^\infty(\Omega)^d} \\
\le & C_1 \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \partial_t N \|\tilde{\mathbf{u}}_N^{n+1}\|_{H_0^1(\Omega)^d} \|\mathbf{v}\|_{L^\infty(\Omega)^d} \\
& + \frac{1}{2} C_1 \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \partial_t N \|\tilde{\mathbf{u}}_N^n\|_{H_0^1(\Omega)^d} \|\mathbf{v}\|_{L^\infty(\Omega)^d} \\
\le & C_1 \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \partial_t N (\|\tilde{\mathbf{u}}_N^{n+1}\|_{H_0^1(\Omega)^d} + \|\tilde{\mathbf{u}}_N^n\|_{H_0^1(\Omega)^d}) \|\mathbf{v}\|_{L^\infty(\Omega)^d} \\
\le & C_1 \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \partial_t N (\|\tilde{\mathbf{u}}_N^{n+1}\|_{H_0^1(\Omega)^d} + \|\tilde{\mathbf{u}}_N^n\|_{H_0^1(\Omega)^d}) \|\mathbf{v}\|_{\mathcal{E}}. \quad (33)
\end{align*}
$$

We get from $E_N \subset V_N$ that $A_p(t) = 0$. Next, we note that

$$
A_f(t) \le C_{\text{sob}} \|\mathbf{v}\|_{\mathcal{E}} \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \partial_t N \|\mathbf{f}_N^{n+1}\|_{H^{-1}(\Omega)^d}. \quad (34)
$$

where $C_{\text{sob}}$ is such that $\|\mathbf{v}\|_{L^2(\Omega)^d} \le C_{\text{sob}} \|\mathbf{v}\|_{H_0^1(\Omega)^d}$ for any $\mathbf{v} \in H_0^1(\Omega)^d$. Summing Equations (32), (33), (34), we obtain

$$
A(t) \le C \|\mathbf{v}\|_{\mathcal{E}} \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \partial_t N (\|\tilde{\mathbf{u}}_N^{n+1}\|_{H_0^1(\Omega)^d} + \|\tilde{\mathbf{u}}_N^n\|_{H_0^1(\Omega)^d} + \|\mathbf{f}_N^{n+1}\|_{H^{-1}(\Omega)^d})
$$

where $C = \frac{3}{2}C_1 + C_{\text{sob}} + 1$. This implies

$$
|\tilde{\mathbf{u}}_N(t + \tau) - \tilde{\mathbf{u}}_N(t)|_{*,N} \le C \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \partial_t N (\|\tilde{\mathbf{u}}_N^{n+1}\|_{H_0^1(\Omega)^d} + \|\tilde{\mathbf{u}}_N^n\|_{H_0^1(\Omega)^d} + \|\mathbf{f}_N^{n+1}\|_{H^{-1}(\Omega)^d}).
$$
Since $\sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \partial_t N \le \tau + \partial_t N$ for any $t \in (0, T - \tau)$ we then obtain

$$|\tilde{u}_N(t + \tau) - \tilde{u}_N(t)|_{*,N}^2 \le 3C^2(\tau + \partial_t N) \sum_{n=1}^{N-1} \chi_{N,\tau}^n(t) \partial_t N (\|\tilde{u}_N^{n+1}\|_{H_0^1(\Omega)^d}^2 + \|\tilde{u}_N^n\|_{H_0^1(\Omega)^d}^2 + \|f_N^{n+1}\|_{H^{-1}(\Omega)^d}^2).$$

Noting that $\int_0^{T-\tau} \chi_{N,\tau}^n(t) dt \le \tau$ for any $n \in [1, N-1]$ yields

$$\begin{aligned} & \int_0^{T-\tau} |\tilde{u}_N(t + \tau) - \tilde{u}_N(t)|_{*,N}^2 dt \\ & \le 3C^2(\tau + \partial_t N) \sum_{n=1}^{N-1} \partial_t N (\|\tilde{u}_N^{n+1}\|_{H_0^1(\Omega)^d}^2 + \|\tilde{u}_N^n\|_{H_0^1(\Omega)^d}^2 + \|f_N^{n+1}\|_{H^{-1}(\Omega)^d}^2) \int_0^{T-\tau} \chi_{N,\tau}^n(t) dt \\ & \le 3C^2(\tau + \partial_t N) \tau (2\|\tilde{u}_N\|_{L^2(0,T;H_0^1(\Omega)^d)}^2 + \|f\|_{L^2(0,T;H^{-1}(\Omega)^d)}^2) \le C_2\tau(\tau + \partial_t N) \end{aligned}$$

which gives the expected result.

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

**Remark 3.5:** Note that the property $\|v\|_{L^\infty(\Omega)^d} \le \|v\|_E$ for any $v \in E_N$ is only used in (33) to estimate the discrete convective term.

Note that Lemma 3.3 is only the first step for proving an estimate on the time translates of the predicted velocity in the $L^2(L^2)$ norm. The next steps are based on the Lions-like result given below, for which we need the following result.

**Lemma 3.4:** Under Assumptions (2), (14), (15), (21), (24), (25), let $(v_N)_{N \ge 1}$ be a sequence of functions of $L^2(\Omega)^d$ such that $(v_N)_{N \ge 1}$ converges to $v$ in $L^2(\Omega)^d$. Then the sequence $(P_{V_N(\Omega)}v_N)_{N \ge 1}$ converges to $P_{V_N(\Omega)}v$ in $L^2(\Omega)^d$.

**PROOF.** Using the fact that $(v_N)_{N \ge 1}$ is bounded in $L^2(\Omega)^d$ we obtain that the sequence $(P_{V_N(\Omega)}v_N)_{N \ge 1}$ is bounded in $L^2(\Omega)^d$. Hence there exists a subsequence still denoted by $(P_{V_N(\Omega)}v_N)_{N \ge 1}$ that converges to a function $\tilde{v}$ weakly in $L^2(\Omega)^d$. Let $\xi \in H^1(\Omega) \cap L_0^2(\Omega)$. Using the fact that for any $N \ge 1$ we have $P_{V_N(\Omega)}v_N \in V_N$ and $P_{M_N}\xi \in M_N$ we obtain

$$(P_{V_N(\Omega)}v_N, \nabla P_{M_N}\xi) = 0, \text{ for any } N \ge 1.$$

Using the weak convergence of the sequence $(P_{V_N(\Omega)}v_N)_{N \ge 1}$ in $L^2(\Omega)^d$ and the strong convergence of the sequence $(P_{M_N}\xi)_{N \ge 1}$ in $H^1(\Omega)$ and passing to the limit in the previous identity gives

$$(\tilde{v}, \nabla \xi) = 0, \text{ for any } \xi \in H^1(\Omega).$$

We then obtain that $\tilde{v} \in V(\Omega)$. Let $\varphi \in W$ be given. Using the fact that $\Pi_N\varphi \in V_N$ for any $N \ge 1$, we obtain

$$(v_N, \Pi_N \varphi) = (P_{V_N(\Omega)}v_N, \Pi_N \varphi), \text{ for any } N \ge 1.$$

Using the convergence to $\varphi$ of the sequence $(\Pi_N\varphi)_{N \ge 1}$ in $L^2(\Omega)^d$, the weak convergence of the sequence $(P_{V_N(\Omega)}v_N)_{N \ge 1}$ in $L^2(\Omega)^d$, passing to the limit in the previous identity and using the fact that $(v, \varphi) = (P_v v, \varphi)$ for any $\varphi \in W$ give

$$(P_v v, \varphi) = (\tilde{v}, \varphi) \text{ for any } \varphi \in W.$$

Using the fact that $W$ is dense in $V(\Omega)$ for the $L^2(\Omega)^d$ norm we obtain that $\tilde{v} = P_v v$ and that the sequence $(P_{V_N(\Omega)}v_N)_{N \ge 1}$ converges to $P_{V_N(\Omega)}v$ weakly in $L^2(\Omega)^d$. We can write

$$\|P_{V_N(\Omega)}v_N\|_{L^2(\Omega)^d}^2 = (v_N, P_{V_N(\Omega)}v_N), \text{ for any } N \ge 1.$$
Using the convergence of the sequence $(v_N)_{N \ge 1}$ to $v$ in $L^2(\Omega)^d$ and the weak convergence of the sequence $(P_{V_N(\Omega)}v_N)_{N \ge 1}$ to $P_{V(\Omega)}v$ in $L^2(\Omega)^d$ we obtain

$$ \lim_{N \to +\infty} \| P_{V_N(\Omega)} v_N \|_{L^2(\Omega)^d}^2 = (v, P_{V(\Omega)} v) = \| P_{V(\Omega)} v \|_{L^2(\Omega)^d}^2. $$

The weak convergence of the sequence $(P_{V_N(\Omega)}v_N)_{N \ge 1}$ to $P_{V(\Omega)}v$ in $L^2(\Omega)^d$ and the convergence of the sequence $(\|P_{V_N(\Omega)}v_N\|_{L^2(\Omega)^d})_{N \ge 1}$ to $\|P_{V(\Omega)}v\|_{L^2(\Omega)^d}$ give the expected result. $\square$

**Remark 3.6:** Note that, in the preceding proof, we only use the assumptions weaker than (24) and (25) that $W$ is a dense subset of $V(\Omega)$ for the norm of $L^2(\Omega)^d$ and that, for any $N \ge 1$ there exists a mapping $\Pi_N: W \to V_N$ such that for any $\varphi \in W$ such that the sequence $(\Pi_N \varphi)_{N \ge 1}$ converges to $\varphi$ in $L^2(\Omega)^d$.

**Lemma 3.5 (Lions-like):** Under Assumptions (2), (14), (15), (21) (24), (25), we have

$$ \forall \varepsilon > 0, \exists C_{\varepsilon} > 0, \exists N_{\varepsilon} \ge 1, \forall N \ge N_{\varepsilon}, \forall w \in \mathbf{X}_N $$

$$ \|P_{V_N} w\|_{L^2(\Omega)^d} \le \varepsilon \|w\|_{H_0^1(\Omega)^d} + C_{\varepsilon} |w|_{*,N}. \quad (35) $$

**PROOF.** Let us assume that (35) does not hold. This means that there exists $\varepsilon > 0$ such that

$$ \forall C \ge 1, \forall M \ge 1, \exists N \ge M, \exists w \in \mathbf{X}_N $$

$$ \|P_{V_N} w\|_{L^2(\Omega)^d} > \varepsilon \|w\|_{H_0^1(\Omega)^d} + C |w|_{*,N}. \quad (36) $$

Let us set $\nu(0) = 0$, and let us build the infinite set $\mathcal{I} = \{\nu(n), n \in \mathbb{N}^*\} \subset \mathbb{N}^*$ and the sequence $(w_N)_{N \in \mathcal{I}}$ by induction. For $n \in \mathbb{N}^*$, we select $C = n$ and $M = \nu(n-1)+1$ in (36). We get the existence of a given $N \ge M$ and of a given $w \in \mathbf{X}_N(\Omega)$ such that

$$ \|P_{V_N} w\|_{L^2(\Omega)^d} > \varepsilon \|w\|_{H_0^1(\Omega)^d} + n |w|_{*,N}. $$

We then define $\nu(n) = N$ and $w_N = w/\|P_{V_N}w\|_{L^2(\Omega)^d}$. Let us denote $\mathcal{I} = \{\nu(n) \text{ such that } n \in \mathbb{N}^*\} \subset \mathbb{N}^*$. We denote by $\nu^{-1}$ the reciprocal function of $\nu: \mathbb{N}^* \to \mathcal{I}$. Since the mapping $\nu$ is strictly increasing, we then indeed get that $\mathcal{I}$ is infinite. We can then write, for any $N \in \mathcal{I}$,

$$ 1 > \varepsilon \|w_N\|_{H_0^1(\Omega)^d} + \nu^{-1}(N) |w_N|_{*,N}, \text{ for any } N \in \mathcal{I}. $$

It then follows from the latter inequality that the sequence $(w_N)_{N \in \mathcal{I}}$ is bounded in $H_0^1(\Omega)^d$ and that $|w_N|_{*,N} \to 0$ as $N \to \infty$ with $N \in \mathcal{I}$. Hence there exists an infinite set $\mathcal{J} \subset \mathcal{I}$ such that $(w_N)_{N \in \mathcal{J}}$ converges in $L^2(\Omega)^d$ to a function $w \in H_0^1(\Omega)^d$ while $|w_N|_{*,N} \to 0$ as $N \to \infty$ with $N \in \mathcal{J}$. We notice that $\|P_{V_N}w_N\|_{L^2(\Omega)^d} = 1$ for all $N \in \mathcal{J}$.

Using Lemma 3.4 we have $P_{V_N}w_N \to P_{V(\Omega)}w$ in $L^2(\Omega)^d$ as $N \to \infty$ with $N \in \mathcal{J}$, which therefore implies $\|P_{V}w\|_{L^2(\Omega)^d} = 1$.

Let $\varphi \in W$ be given. For any $N \in \mathcal{J}$, by definition of $|w_N|_{*,N}$, we have

$$ (w_N, \Pi_N \varphi) \le |w_N|_{*,N} \|\Pi_N \varphi\|_E. $$

Using the fact that, for any $N \in \mathcal{J}$, $\Pi_N \varphi \in V_N$, we then obtain

$$ (w_N, \Pi_N \varphi) = (P_{V_N} w_N, \Pi_N \varphi) \le |w_N|_{*,N} \|\Pi_N \varphi\|_E. $$

Since $\|\Pi_N \varphi\|_E$ remains bounded by assumption (24), letting $N \to \infty$ with $N \in \mathcal{J}$ in this inequality yields that

$$ (P_{V(\Omega)} w, \varphi) = 0, \text{ for any } \varphi \in W. $$

Using the fact that $W$ is dense in $V(\Omega)$ for the $L^2(\Omega)^d$ norm we can therefore let $\varphi \to P_{V(\Omega)}w$ in $L^2(\Omega)^d$. This yields

$$ \|P_{V(\Omega)} w\|_{L^2(\Omega)^d}^2 = 0, $$

which contradicts $\|P_{V(\Omega)}w\|_{L^2(\Omega)^d} = 1$. $\square$
**Remark 3.7:** Note that, in the preceding proof, we only used the assumptions weaker than (24) and (25) that $\mathbf{W}$ is a dense subset of $\mathbf{V}(\Omega)$ for the norm of $L^2(\Omega)^d$ and that, for any $N \ge 1$ there exists a mapping $\Pi_N: \mathbf{W} \to \mathbf{E}_N$ such that for any $\varphi \in \mathbf{W}$ such that the sequence $(\Pi_N \varphi)_{N \ge 1}$ converges to $\varphi$ in $L^2(\Omega)^d$ and the sequence $(\|\Pi_N \varphi\|_{\mathbf{E}})_{N \ge 1}$ is bounded.

Our aim is now to use Lemma 3.3 on the time translates of $(\tilde{u}_N)_{N \ge 1}$ for the $L^2(|\cdot|_{*,N})$ semi-norm and (35) in the above lemma, in order to obtain an estimate on the time translates for the $L^2(0, T; L^2(\Omega)^d)$ norm, as stated by the next lemma.

**Lemma 3.6 (L² estimate on the time translates):** Under Assumptions (2), (7), (13), (14), (15), (21), (24), (25), the sequence $(\tilde{u}_N)_{N \ge 1}$ satisfies

$$ \int_0^{T-\tau} \|\tilde{u}_N(t+\tau) - \tilde{u}_N(t)\|_{L^2(\Omega)^d}^2 dt \to 0 \text{ as } \tau \to 0, \text{ uniformly with respect to } N, \quad (37) $$

and is therefore relatively compact in $L^2(0, T; L^2(\Omega)^d)$.

**PROOF.** Let us show that $A_N(\tau) \to 0$ as $\tau \to 0$ uniformly with respect to $N$, where we define for any $\tau \in (0, T)$

$$ A_N(\tau) := \int_0^{T-\tau} \|\tilde{u}_N(t+\tau) - \tilde{u}_N(t)\|_2^2 dt. $$

For any $\mathbf{w} \in L^2(\Omega)^d$ and $\mathbf{v} \in \mathbf{V}_N$, we have

$$ \mathbf{w} = \mathbf{w} - \mathcal{P}_{\mathbf{V}_N}(\mathbf{w}) + \mathcal{P}_{\mathbf{V}_N}(\mathbf{w}) = \mathbf{w} - \mathbf{v} - \mathcal{P}_{\mathbf{V}_N}(\mathbf{w} - \mathbf{v}) + \mathcal{P}_{\mathbf{V}_N}(\mathbf{w}). $$

Using the fact that $\mathcal{P}_{\mathbf{V}_N}$ is an orthogonal projection, we have

$$ \|\mathbf{w} - \mathbf{v} - \mathcal{P}_{\mathbf{V}_N}(\mathbf{w} - \mathbf{v})\|_{L^2(\Omega)^d} \le \|\mathbf{w} - \mathbf{v}\|_{L^2(\Omega)^d}, $$

which leads to

$$ \|\mathbf{w}\|_{L^2(\Omega)^d} \le \|\mathbf{w} - \mathbf{v}\|_{L^2(\Omega)^d} + \|\mathcal{P}_{\mathbf{V}_N}(\mathbf{w})\|_{L^2(\Omega)^d}, $$

giving

$$ \forall \mathbf{w} \in L^2(\Omega)^d, \forall \mathbf{v} \in \mathbf{V}_N, \quad \|\mathbf{w}\|_{L^2(\Omega)^d}^2 \le 2\|\mathbf{w} - \mathbf{v}\|_{L^2(\Omega)^d}^2 + 2\|\mathcal{P}_{\mathbf{V}_N}(\mathbf{w})\|_{L^2(\Omega)^d}^2. $$

We apply the preceding inequality, letting $\mathbf{w} = \tilde{u}_N(t+\tau) - \tilde{u}_N(t)$ and $\mathbf{v} = \mathbf{u}_N(t+\tau) - \mathbf{u}_N(t)$ and we integrate on $t \in (0, T - \tau)$. Setting

$$ B_N(\tau) = \int_0^{T-\tau} \|(\tilde{u}_N - \mathbf{u}_N)(t+\tau) - (\tilde{u}_N - \mathbf{u}_N)(t)\|_2^2 dt $$

we obtain

$$ A_N(\tau) \le 2B_N(\tau) + 2 \int_0^{T-\tau} \|\mathcal{P}_{\mathbf{V}_N}(\tilde{u}_N(t+\tau) - \tilde{u}_N(t))\|_{L^2(\Omega)^d}^2 dt. \quad (38) $$

Let $\zeta > 0$ be given.

From (28), we have $B_N(\tau) \le 4C_1^2 \delta_{N}$. Let $N_B \ge 1$ be such that $4C_1^2 \delta_{N} \le \zeta$ for all $N \ge N_B$. Since $B_N(\tau) \to 0$ as $\tau \to 0$ for $N = 1, \dots, N_B$, we can choose $\tau_B$ such that, for any $0 \le \tau \le \tau_B$ and $N = 1, \dots, N_B$, we have $B_N(\tau) \le \zeta$.

This yields $B_N(\tau) \le \zeta$ for all $N \ge 1$ and any $0 \le \tau < \tau_B$.

Our aim is now to use (35) in Lemma 3.5, which implies

$$ \forall \varepsilon > 0, \exists C_{\varepsilon} > 0, \exists N_{\varepsilon} \ge 1, \forall N \ge N_{\varepsilon}, \forall \mathbf{w} \in \mathbf{X}_N $$

$$ \|\mathcal{P}_{\mathbf{V}_N}(\mathbf{w})\|_{L^2(\Omega)^d}^2 \le 2\varepsilon^2 \|\mathbf{w}\|_{H_0^1(\Omega)^d}^2 + 2C_{\varepsilon}^2 |\mathbf{w}|_{*,N}^2, $$
letting $w = \tilde{u}_N(t + \tau) - \tilde{u}_N(t)$.

From (27), we have

$$ \int_{0}^{T-\tau} \|\tilde{u}_N(t + \tau) - \tilde{u}_N(t)\|_{H_0^1(\Omega)^d}^2 \, dt \le 4C_1^2. \quad (39) $$

We then select $\varepsilon$ such that

$$ 2\varepsilon^2 4C_1^2 = \zeta. \quad (40) $$

Then there exists $C_\zeta > 0$ and $N_\zeta$ such that for any $N \ge N_\zeta$ and for any $\tau \in (0, T)$ and for any $t \in (0, T - \tau)$,

$$ \|\mathcal{P}_{V_N}(\tilde{u}_N(t + \tau) - \tilde{u}_N(t))\|_{L^2(\Omega)^d}^2 \le 2\varepsilon^2 \|\tilde{u}_N(t + \tau) - \tilde{u}_N(t)\|_{H_0^1(\Omega)^d}^2 + C_\zeta^2 |\tilde{u}_N(t + \tau) - \tilde{u}_N(t)|_{*,N}^2. $$

Integrating the previous relation provides, using (39) and (40),

$$ \int_{0}^{T-\tau} \|\mathcal{P}_{V_N}(\tilde{u}_N(t + \tau) - \tilde{u}_N(t))\|_{L^2(\Omega)^d}^2 \, dt \le \zeta + 2C_\zeta^2 \int_{0}^{T-\tau} |\tilde{u}_N(t + \tau) - \tilde{u}_N(t)|_{*,N}^2 \, dt. $$

Thus, owing to Lemma 3.3, we have for any $N \ge N_\zeta$ and for any $\tau \in (0, T)$

$$ \int_{0}^{T-\tau} \|\mathcal{P}_{V_N}(\tilde{u}_N(t + \tau) - \tilde{u}_N(t))\|_{L^2(\Omega)^d}^2 \, dt \le \zeta + 2C_\zeta^2 C_2 \tau (\tau + \delta_N), $$

and therefore, for any $N \ge N_\zeta$ and for any $\tau \in (0, T)$, using (38)

$$ A_N(\tau) \le 2B_N(\tau) + 2\zeta + 4C_\zeta^2 C_2 \tau (\tau + \delta_N). $$

We now choose $\tau_C > 0$, such that, for any $0 \le \tau \le \tau_C$, $4C_\zeta^2 C_2 \tau (\tau + T) \le \zeta$.

We then get, for any $N \ge N_\zeta$ and $0 \le \tau \le \min(\tau_B, \tau_C)$,

$$ A_N(\tau) \le 2\zeta + 2\zeta + \zeta = 5\zeta. $$

Since $A_N(\tau) \to 0$ as $\tau \to 0$ for $N = 1, \dots, N_\zeta$, we can choose $\tau_D$ such that, for any $0 \le \tau \le \tau_D$ and $N = 1, \dots, N_\zeta$, we have $A_N(\tau) \le 5\zeta$.

We then obtain that $A_N(\tau) \le 5\zeta$ for any $\tau \in [0, \min(\tau_B, \tau_C, \tau_D)]$ and $N \ge 1$. The proof of (37) is thus complete. $\square$

## 3.4 Convergence of the projection scheme to a weak solution

The convergence proof given below relies on the following result which is a consequence of the Kolmogorov theorem, as noticed in [16, Corollaire 4.41].

**Theorem 3.7:** Let $1 \le p < +\infty$. Let $B$ be a Banach space and let $X$ be a Banach space compactly embedded in $B$. Let $(u_N)_{N \ge 1}$ be a sequence of $L^p(0, T; B)$ satisfying the following conditions:

1. The sequence $(u_N)_{N \ge 1}$ is bounded in $L^p(0, T; B)$.

2. The sequence $(u_N)_{N \ge 1}$ is bounded in $L^1(0, T; X)$.

3. The sequence $(u_N)_{N \ge 1}$ satisfies

$$ \int_{0}^{T-\tau} \|u_N(t + \tau) - u_N(t)\|_{B}^{p} \, dt \to 0 \text{ as } \tau \to 0, \quad (41) $$

uniformly with respect to $N \ge 1$.
Then the sequence $(u_N)_{N \ge 1}$ is relatively compact in $L^p(0, T; B)$.

Let us now detail the way this result is used. By Lemma 3.2, up to a subsequence, the sequence $(\tilde{u}_N)_{N \ge 1}$ weakly converges to some limit $\bar{u}$ in $L^2(0, T; L^2(\Omega)^d)$, and owing to (28), so does the sequence $(u_N)_{N \ge 1}$. There remains to check that $\bar{u}$ is a weak solution to (1) in the sense of Definition 1.1, by passing to the limit in the scheme with suitable test functions. This is possible using the compactness theorem 3.7 and the time translate estimate lemma 3.6.

**Theorem 3.8:** Under Assumptions (2), (7), let $(\partial_t N, X_N, M_N, \Pi_N)_{N \ge 1}$ be a sequence such that (13), (14), (15), (21), (24), (25) hold. Then there exists a subsequence samely denoted, and a function $\bar{u} \in L^2(0, T; H_0^1(\Omega)^d \cap V(\Omega)) \cap L^\infty(0, T; L^2(\Omega)^d)$, such that the corresponding sequences $(\tilde{u}_N)_{N \ge 1}$ and $(u_N)_{N \ge 1}$ defined by the incremental projection scheme (18) converge to $\bar{u}$ strongly in $L^2(0, T; L^2(\Omega)^d)$. Moreover the function $\bar{u}$ is a weak solution to (1) in the sense of Definition 1.1.

**PROOF.** We proceed in several steps.

* **Step 1: compactness and convergence in $L^2$.** Using Lemma 3.6, we can apply Theorem 3.7. We get the existence of $\bar{u}$ in $L^2(0, T; L^2(\Omega)^d)$ and of a subsequence of $(\partial_t N, X_N, M_N, \Pi_N)_{N \ge 1}$, still denoted $(\partial_t N, X_N, M_N, \Pi_N)_{N \ge 1}$, such that the corresponding sequences $(u_N)_{N \ge 1}$ and $(\tilde{u}_N)_{N \ge 1}$ converge to $\bar{u}$ in $L^2(0, T; L^2(\Omega)^d)$. Owing to Lemma 3.2, we get that $\bar{u} \in L^\infty(0, T; L^2(\Omega)^d) \cap L^2(0, T; H_0^1(\Omega)^d \cap V(\Omega))$.

* **Step 2: convergence towards a weak solution.** There remains to show that $\bar{u}$ is a weak solution in the sense of Definition 1.1 and in particular that $\bar{u}$ satisfies (12). Let $\psi \in C_c^\infty([0, T])$ and $\varphi \in W$. Let $v$ be the function defined by $v(t, x) = \psi(t)\varphi(x)$ for any $(t, x) \in [0, T) \times \Omega$. Let $N \ge 1$. Let $(v_N^n)_{n \in \mathbb{N}}$ be the sequence of functions of $X_N \cap V_N$ defined by $v_N^n = \psi(t_N^n)\Pi_N \varphi$ for any $n \in [0, N]$. Let $v_N: (0, T) \to H_0^1(\Omega)^d$ be defined by

$$v_N(t) = \sum_{n=0}^{N-1} \mathbb{1}_{(t_N^n, t_{N+1}^n]}(t) v_N^{n+1}.$$

We remark that we have

$$\|v_N - v\|_{L^\infty(0,T; H_0^1(\Omega)^d)} \le \partial_t N \|\psi'\|_{L^\infty(0,T)} \|\Pi_N \varphi\|_E \\ + \|\psi\|_{L^\infty(0,T)} \|\Pi_N \varphi - \varphi\|_{H_0^1(\Omega)^d}, \text{ for any } N \ge 1. \quad (42)$$

Consequently we obtain using (24) that the sequence $(v_N)_{N \ge 1}$ converges to $v$ in $L^\infty(0, T; H_0^1(\Omega)^d)$. Multiplying (18b) by $\partial_t N v_N^{n+1}$, using the fact that $(\tilde{u}_N^{n+1}, v_N^{n+1}) = (u_N^{n+1}, v_N^{n+1})$ for any $n \in [0, N-1]$ and summing over $n \in [0, N-1]$ lead to

$$\sum_{n=0}^{N-1} (u_N^{n+1} - u_N^n, v_N^{n+1}) + \partial_t N \sum_{n=0}^{N-1} b(\tilde{u}_N^n, \tilde{u}_N^{n+1}, v_N^{n+1}) \\ + \int_0^T \int_\Omega \nabla \tilde{u}_N : \nabla v_N \,dx \,dt = \int_0^T \langle f, v_N \rangle \,dt. \quad (43)$$
Using the fact that $v_N^N = 0$ in $\Omega$ the first term of the left hand side reads

$$ 
\begin{aligned}
\sum_{n=0}^{N-1} (\boldsymbol{u}_N^{n+1} - \boldsymbol{u}_N^n, \boldsymbol{v}_N^{n+1}) &= -\sum_{n=0}^{N-1} (\boldsymbol{u}_N^n, \boldsymbol{v}_N^{n+1} - \boldsymbol{v}_N^n) - (\boldsymbol{u}_N^0, \boldsymbol{v}_N^0) \\
&= -\int_0^T \psi'(t)(\boldsymbol{u}_N(t), \Pi_N \varphi) \, dt - \psi(0)(\mathcal{P}_{\boldsymbol{v}_N} \boldsymbol{u}_0, \Pi_N \varphi) \\
&= -\int_0^T \psi'(t)(\boldsymbol{u}_N(t), \Pi_N \varphi) \, dt - \psi(0)(\boldsymbol{u}_0, \Pi_N \varphi).
\end{aligned}
$$

Since the sequence $(\boldsymbol{u}_N)_{N \ge 1}$ converges to $\tilde{\boldsymbol{u}}$ in $L^2(0, T; L^2(\Omega)^d)$ and the sequence $(\Pi_N \varphi)_{N \ge 1}$ converges to $\varphi$ in $L^2(\Omega)^d$, we can then write

$$ \lim_{N \to +\infty} \sum_{n=0}^{N-1} (\boldsymbol{u}_N^{n+1} - \boldsymbol{u}_N^n, \boldsymbol{v}_N^{n+1}) = - \int_0^T (\tilde{\boldsymbol{u}}, \partial_t \boldsymbol{v}) \, dt - (\boldsymbol{u}_0, \boldsymbol{v}(0, \cdot)). \quad (44) $$

Using the fact that the initial predicted velocity is zero we have

$$ 
\begin{aligned}
\partial_t \sum_{n=0}^{N-1} b(\tilde{\boldsymbol{u}}_N^n, \tilde{\boldsymbol{u}}_N^{n+1}, \boldsymbol{v}_N^{n+1}) &= \partial_t \sum_{n=1}^{N-1} b(\tilde{\boldsymbol{u}}_N^n, \tilde{\boldsymbol{u}}_N^{n+1}, \boldsymbol{v}_N^{n+1}) \\
&= \frac{1}{2} \int_{\partial t_N}^T ((\tilde{\boldsymbol{u}}_N(t - \partial_t N) \cdot \boldsymbol{\nabla}) \tilde{\boldsymbol{u}}_N(t), \boldsymbol{v}_N(t)) \, dt \\
&\quad - \frac{1}{2} \int_{\partial t_N}^T ((\tilde{\boldsymbol{u}}_N(t - \partial_t N) \cdot \boldsymbol{\nabla}) \boldsymbol{v}_N(t), \tilde{\boldsymbol{u}}_N(t)) \, dt.
\end{aligned}
$$

The first term of the discrete convective term can be expressed by

$$ 
\begin{aligned}
\int_{\partial t_N}^T ((\tilde{\boldsymbol{u}}_N(t - \partial_t N) \cdot \boldsymbol{\nabla}) \tilde{\boldsymbol{u}}_N(t), \boldsymbol{v}_N(t)) \, dt &= \int_0^T ((\tilde{\boldsymbol{u}}_N(t) \cdot \boldsymbol{\nabla}) \tilde{\boldsymbol{u}}_N(t), \boldsymbol{v}_N(t)) \, dt \\
&\quad - \int_0^{\partial t_N} ((\tilde{\boldsymbol{u}}_N(t) \cdot \boldsymbol{\nabla}) \tilde{\boldsymbol{u}}_N(t), \boldsymbol{v}_N(t)) \, dt \\
&\quad + \int_{\partial t_N}^T (((\tilde{\boldsymbol{u}}_N(t - \partial_t N) - \tilde{\boldsymbol{u}}_N(t)) \cdot \boldsymbol{\nabla}) \tilde{\boldsymbol{u}}_N(t), \boldsymbol{v}_N(t)) \, dt.
\end{aligned}
$$

Using Lemma 3.1, we get the following bound for the second term of the right hand-side of the previous identity:

$$ \int_0^{\partial t_N} ((\tilde{\boldsymbol{u}}_N(t) \cdot \boldsymbol{\nabla}) \tilde{\boldsymbol{u}}_N(t), \boldsymbol{v}_N(t)) \, dt \le C_1^2 \|\psi\|_{L^\infty(0,T)} \|\Pi_N \varphi\|_{L^\infty(\Omega)^d} \sqrt{\partial t_N}. $$

Using the fact that the sequence $(\|\Pi_N \varphi\|_{L^\infty(\Omega)^d})_{N \ge 1}$ is bounded we obtain that this term tends to zero as $N$ tends to infinity. Using Lemma 3.1, a bound for the third term of the right hand-side of the previous identity is given by

$$ 
\begin{aligned}
& \int_{\partial t_N}^T (((\tilde{\boldsymbol{u}}_N(t - \partial_t N) - \tilde{\boldsymbol{u}}_N(t)) \cdot \boldsymbol{\nabla}) \tilde{\boldsymbol{u}}_N(t), \boldsymbol{v}_N(t)) \, dt \\
&\le C_1 \|\psi\|_{L^\infty(0,T)} \|\Pi_N \varphi\|_{L^\infty(\Omega)^d} \left( \int_0^{T - \partial_t N} \|\tilde{\boldsymbol{u}}_N(t + \partial_t N) - \tilde{\boldsymbol{u}}_N(t)\|_{L^2(\Omega)^d}^2 \right)^{1/2}.
\end{aligned}
$$
Using Lemma 3.6 and the fact that $(\|\Pi_N\varphi\|_{L^\infty(\Omega)^d})_{N\ge 1}$ is bounded we obtain that the preceding term tends to zero as $N$ tends to infinity. Using the convergence of the sequence $(\tilde{u}_N)_{N\ge 1}$ in $L^2(0,T; L^4(\Omega)^d)$, the weak convergence of the sequence $(\tilde{u}_N)_{N\ge 1}$ in $L^2(0,T; H_0^1(\Omega)^d)$ and the convergence of the sequence $(v_N)_{N\ge 1}$ in $L^\infty(0,T; L^4(\Omega)^d)$ we obtain

$$ \lim_{N \to +\infty} \int_0^T ((\tilde{u}_N(t) \cdot \nabla) \tilde{u}_N(t), v_N(t)) \, dt = \int_0^T ((\bar{u} \cdot \nabla) \bar{u}, v) \, dt. $$

The second term of the discrete convective term is such that

$$ \int_{\partial N}^T ((\tilde{u}_N(t - \partial_t N) \cdot \nabla) v_N(t), \tilde{u}_N(t)) \, dt = \int_0^T ((\tilde{u}_N(t) \cdot \nabla) \tilde{u}_N(t), v_N(t)) \, dt - \int_0^{\partial t_N} ((\tilde{u}_N(t) \cdot \nabla) v_N(t), \tilde{u}_N(t)) \, dt + \int_{\partial N}^T (((\tilde{u}_N(t - \partial_t N) - \tilde{u}_N(t)) \cdot \nabla) v_N(t), \tilde{u}_N(t)) \, dt. $$

A bound for the second term of the previous identity is obtained by

$$ \int_0^{\partial t_N} ((\tilde{u}_N(t) \cdot \nabla) v_N(t), \tilde{u}_N(t)) \, dt \le \|\psi\|_{L^\infty(0,T)} \|\Pi_N\varphi\|_{H_0^1(\Omega)^d} \int_0^{\partial t_N} \|\tilde{u}_N(t)\|_{L^4(\Omega)^d}^2 \, dt. $$

Using the fact that the sequence $(\tilde{u}_N)_{N\ge 1}$ is bounded in $L^2(0,T; L^4(\Omega)^d)$ and the fact that the sequence $(\Pi_N\varphi\|_{H_0^1(\Omega)^d})_{N\ge 1}$ is bounded we obtain that this term tends to zero as $N$ tends to infinity. The third term of the previous equation can be bounded as follows

$$ \begin{aligned} & \int_{\partial N}^T (((\tilde{u}_N(t - \partial_t N) - \tilde{u}_N(t)) \cdot \nabla) v_N(t), \tilde{u}_N(t)) \, dt \\ & \le C_1 \|\psi\|_{L^\infty(0,T)} \|\Pi_N\varphi\|_{H_0^1(\Omega)^d} \left( \int_0^{T-\partial_N} \|\tilde{u}_N(t + \partial_N) - \tilde{u}_N(t)\|_{L^4(\Omega)^d}^2 \right)^{1/2} \\ & \qquad \times \left( \int_0^{T-\partial_N} \|\tilde{u}_N(t)\|_{L^4(\Omega)^d}^2 \right)^{1/2} \\ & \le C_1 \|\psi\|_{L^\infty(0,T)} \|\Pi_N\varphi\|_{H_0^1(\Omega)^d} \left( \int_0^{T-\partial_N} \|\tilde{u}_N(t + \partial_N) - \tilde{u}_N(t)\|_{L^2(\Omega)^d}^2 \right)^{1/8} \\ & \qquad \times \left( \int_0^{T-\partial_N} \|\tilde{u}_N(t + \partial_N) - \tilde{u}_N(t)\|_{L^6(\Omega)^d}^2 \right)^{3/8} \left( \int_0^{T-\partial_N} \|\tilde{u}_N(t)\|_{L^4(\Omega)^d}^2 \right)^{1/2} \end{aligned} $$

Using the estimates on the translation of the sequence $(\tilde{u}_N)_{N\ge 1}$ given in Lemma 3.6, the fact that the sequence $(\tilde{u}_N)_{N\ge 1}$ is bounded in $L^2(0,T; L^6(\Omega)^d)$ and the fact that $(\|\Pi_N\varphi\|_{H_0^1(\Omega)^d})_{N\ge 1}$ is bounded we obtain that this term tends to zero as $N$ tends to infinity. Using the convergence of the sequence $(\tilde{u}_N)_{N\ge 1}$ in $L^2(0,T; L^4(\Omega)^d)$ and the convergence of the sequence $(v_N)_{N\ge 1}$ in $L^\infty(0,T; H_0^1(\Omega)^d)$ we obtain

$$ \lim_{N \to +\infty} \int_0^T ((\tilde{u}_N(t) \cdot \nabla) v_N(t), \tilde{u}_N(t)) \, dt = \int_0^T ((\bar{u} \cdot \nabla) v, \bar{u}) \, dt. $$

We then obtain

$$ \lim_{N \to +\infty} \partial_t N \sum_{n=0}^{N-1} b(\tilde{u}_N^n, \tilde{u}_N^{n+1}, v_N^{n+1}) = \int_0^T ((\bar{u} \cdot \nabla) \bar{u}, v) \, dt. \quad (45) $$
Using the weak convergence of the sequence $(\nabla \tilde{u}_N)_{N \ge 1}$ in $L^2(0, T; L^2(\Omega)^d)$ and the convergence of the sequence $(\nabla v_N)_{N \ge 1}$ in $L^2(0, T; L^2(\Omega)^d)$ we obtain

$$ \lim_{N \to +\infty} \int_0^T \int_{\Omega} \nabla \tilde{u}_N : \nabla v_N \, dx \, dt = \int_0^T \int_{\Omega} \nabla \tilde{u} : \nabla v \, dx \, dt. \quad (46) $$

Using the convergence of the sequence $(v_N)_{N \ge 1}$ in $L^2(0, T; L^2(\Omega)^d)$ we obtain

$$ \lim_{N \to +\infty} \int_0^T \langle \mathbf{f}, v_N \rangle \, dt = \int_0^T \langle \mathbf{f}, v \rangle \, dt. \quad (47) $$

Using (44)-(47) and passing to the limit in (43) gives the weak sense (12) with $v$ as test function. We get that (12) holds for all test functions by density of the set of linear combinations of such functions.

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

## 4 Lowest-order Taylor-Hood finite element approximation

In this section, we study the case of the lowest-order Taylor-Hood element. This pair is of particular interest, since it is the lowest-order conforming finite element method which is available on simplicial meshes (triangles in 2D, tetrahedra in 3D). It is such that the discrete velocity space $X_N$ contains continuous piecewise quadratic functions vanishing at the boundary, and the discrete pressure space $M_N$ contains continuous piecewise affine functions.

Our aim is to show that it is possible to construct an operator $\Pi_N$ satisfying the convergence assumptions (21)-(25) done in Section 3.1, under the only standard regularity hypothesis on the meshes that the ratio between exterior and interior diameters of the triangles or tetrahedra is uniformly bounded. We emphasize that the hypotheses which must be verified in Section 3.1 are not equivalent to the inf-sup condition necessary for schemes coupling the velocities and the pressure (as in [21, 15]). Notice that the verification of this inf-sup stability is not easy in the case of the Taylor-Hood approximation (see [1], [3], [36], [5], [2]). Indeed, this stability condition does not hold on general 2D and 3D simplicial meshes, only assuming the standard regularity hypothesis as recalled above. A few pathological cases can be listed in 2D for meshes containing 1, 2 or 3 triangles and in 3D for meshes containing a series of tetrahedra whose edges are all on the boundary, but these meshes are not relevant in the case of a convergence study.

But there also exists standard situations which make the inf-sup condition fail, while considering a convergence study based on a sequence of 3D meshes whose size tends to zero, the regularity property being checked. Indeed, Figure 1 shows the case where the mesh includes at least one tetrahedron (denoted $K_0$) whose all the edges are located on the boundary, with one vertex (denoted by $y_0$) only belonging to $K_0$. Then the inf-sup condition cannot hold. Indeed, for any $v \in X_N$, $v = 0$ on $K_0$, and for $q \in M_N$ such that $q(y_0) = 1 + c$ and $q(y_i) = c$ for any vertex $y_i \ne y_0$, $c$ being adjusted for yielding null mean value for $q$, then $(\nabla q, v) = 0$.

Sufficient conditions on the mesh (which then prevent from the cases listed above) are then requested for the proof of the existence of Fortin operators, that are defined as bounded interpolation operators preserving the discrete divergence of a function. The existence of a Fortin operator is then equivalent to the inf-sup condition (indeed, in the example of Figure 1, no operator can keep the discrete divergence in the case of a function whose support is included in the corner tetrahedron). We refer to [8], [30], [14] for the construction of a Fortin operator in the two dimensional case. The extension of the 2D face-bubble Fortin operator to the 3D case has been done in [11]. It relies on edge-bubble functions, which are modified at the boundary in order that the interpolation vanishes on $\partial\Omega$, under the constraint that each
<figure><img src="image_4.png" /><figcaption>Fig. 1: An example leading to inf-sup condition failure for the lowest-order Taylor-Hood finite element in 3D.</figcaption></figure>

simplex as at least one interior node (this eliminates the above cases of inf-sup condition failure).

In this paper, instead of constructing a complete Fortin operator, we focus on the only properties needed for the operator $\Pi_N$ involved in the convergence of the projection scheme. We recall below the definition of the Lagrange interpolator $\Pi_L$ and its approximation properties (which are elementary compared to that of the Scott-Zhang interpolator, and lead to very short proofs). We then adapt and simplify the divergence correcting linear operator $\Pi_D$ provided by [11] to the case of compact support functions. This enables us to define the operator $\Pi_N$ from $\Pi_L$ and $\Pi_D$ in the case of any simplicial mesh (including all the cases where the inf-sup condition fails).

## 4.1 Regular simplicial mesh

Let us assume that the computational domain $\Omega$ is polygonal ($d = 2$) or polyhedral ($d = 3$). We define a simplex in dimension $d \ge 1$ as the interior of the convex hull of a given set of $d + 1$ points (called its vertices) which are not all contained in the same hyperplane (a simplex is a triangle if $d = 2$ and a tetrahedron if $d = 3$).

We consider a regular simplicial mesh $\mathcal{T}$ of $\Omega$ in the usual sense of the finite element literature [10], as depicted in Figure 2, which means that:

* The family containing all the vertices of the elements of the mesh is denoted by $(\mathbf{y}_i)_{i \in \mathcal{N}}$. The set $\mathcal{N}$ is partitioned into $\mathcal{N} = \mathcal{N}_{\text{int}} \cup \mathcal{N}_{\text{ext}}$, where $(\mathbf{y}_i)_{i \in \mathcal{N}_{\text{ext}}}$ is the set of exterior nodes and $(\mathbf{y}_i)_{i \in \mathcal{N}_{\text{int}}}$ is the set of interior nodes.

* For any $K \in \mathcal{T}$, the family of the vertices of $K$ is denoted by $(\mathbf{y}_i)_{i \in \mathcal{N}_K}$ with $\mathcal{N}_K \subset \mathcal{N}$ contains $d+1$ elements. For any $K \in \mathcal{T}$, we denote by $|K|$ the measure in $\mathbb{R}^d$ of $K$, by $h_K$ the diameter of $K$ and by $\rho_K$ the diameter of the largest ball included in $K$. Then, we define the mesh size $h_\mathcal{T}$ and the mesh regularity $\theta_\mathcal{T}$ by

$$h_\mathcal{T} = \max_{K \in \mathcal{T}} h_K, \quad \theta_\mathcal{T} = \max_{K \in \mathcal{T}} \frac{h_K}{\rho_K}.$$

* For any $K \in \mathcal{T}$, the set of the faces of $K$ is denoted by $\mathcal{F}_K$. For any $F \in \mathcal{F}_K$, the family of the vertices of the simplex $F$ is denoted by $(\mathbf{y}_i)_{i \in \mathcal{N}_F}$ (then $\mathcal{N}_F \subset \mathcal{N}_K$ contains $d$ elements).
* The set $\mathcal{E}$ of the edges of the elements contains all pairs $\{i, j\}$ such that there exists $K \in \mathcal{T}$ with $\{i, j\} \subset \mathcal{N}_K$. For any $\sigma = \{i, j\} \in \mathcal{E}$, we denote by $\mathbf{y}_\sigma = \frac{1}{2}(\mathbf{y}_i + \mathbf{y}_j)$. For all $K \in \mathcal{T}$, we denote by $\mathcal{E}_K$ the subset of $\mathcal{E}$ containing all the edges of $K$.

For $i \in \mathcal{N}$ (resp. $\{i, j\} \in \mathcal{E}$), we define $\mathcal{T}_i$ (resp. $\mathcal{T}_{ij}$) as the set of all $K \in \mathcal{T}$ such that $i \in \mathcal{N}_K$ (resp. $\{i, j\} \in \mathcal{E}_K$) and we denote by

$$ \omega_i = \bigcup_{K \in \mathcal{T}_i} \bar{\mathcal{K}} \quad \text{and} \quad \omega_{ij} = \bigcup_{K \in \mathcal{T}_{ij}} \bar{\mathcal{K}}. $$

<figure><img src="image_5.png" /><figcaption>Fig. 2: An interior part of a 2D mesh. The degrees of freedom for the pressure are in pink, those for the components of the velocity are in black.</figcaption></figure>

In the whole section 4, we denote by $C_i$, with $i \in \mathbb{N}$, non negative real values which may depend on $d$ and on $\Omega$, and increasingly on $\theta_T$. The following lemma is a direct consequence of the definition of $\theta_T$.

**Lemma 4.1:** There exists $C_3$ such that

$$ \mathrm{card}\mathcal{T}_{ij} \leq C_3, $$

$$ \begin{gathered} \forall K \in \mathcal{T}_{ij}, \quad \frac{1}{C_3}h_K^d \leq |\omega_{ij}| \leq C_3h_K^d, \\ \forall K \in \mathcal{T}_{ij}, \quad \mathrm{diam}(\omega_{ij}) \leq C_3h_K. \end{gathered} \quad (48) $$

## 4.2 Discrete spaces and bases

We introduce the set $\mathbb{P}^k(K)$ of the polynomials over $K \in \mathcal{T}$ of degree less than or equal to $k$, regardless of the space dimension. We define

$$ \mathbb{P}^k(\mathcal{T}) = \{p \in C(\bar{\Omega}, \mathbb{R}) \text{ such that } p|_K \in \mathbb{P}^k(K) \text{ for any } K \in \mathcal{T}\}. $$

We denote by $(\varphi_i)_{i \in \mathcal{N}}$ the nodal basis of $\mathbb{P}^1(\mathcal{T})$, that is, for any $i \in \mathcal{N}$, the element $\varphi_i$ of $\mathbb{P}^1(\mathcal{T})$ such that $\varphi_i(\mathbf{y}_i) = 1$ and $\varphi_i(\mathbf{y}_j) = 0$ for all $j \in \mathcal{N} \setminus \{i\}$.

$$ \{\mathbf{x} \in \Omega \text{ such that } \varphi_i(\mathbf{x}) \neq 0\} = \omega_i \quad \text{for any } i \in \mathcal{N}. \quad (49) $$

Then the nodal basis of $\mathbb{P}^2(\mathcal{T})$ is given by the family $((\varphi_i(2\varphi_i - 1))_{i \in \mathcal{N}}, (4\varphi_i\varphi_j)_{\{i,j\} \in \mathcal{E}})$. Denoting this family by $(\phi_i)_{i \in \mathcal{N} \cup \mathcal{E}}$, for any $i \in \mathcal{N} \cup \mathcal{E}$, we have that $\phi_i(\mathbf{y}_i) = 1$ and $\phi_i(\mathbf{y}_j) = 0$ for all $j \in (\mathcal{N} \cup \mathcal{E}) \setminus \{i\}$. Using (49) and using [13, Proposition 11.6] we can state the following lemma.
**Lemma 4.2:** There exists $C_4$ such that

$$ \forall K \in \mathcal{T}, \max_{i \in \mathcal{N}_K} \|\nabla \varphi_i\|_{L^\infty(\Omega_i)^d} \le \frac{C_4}{h_K} \qquad (50) $$

which leads that there exists $C_5$

$$ \forall K \in \mathcal{T}, \forall \mathbf{x} \in K, \forall i \in \mathcal{N}_K \cup \mathcal{E}_K, |\nabla \phi_i(\mathbf{x})| \le \frac{C_5}{h_K}. \qquad (51) $$

The spaces of approximation for the velocity and the pressure are defined by

$$ \mathbf{X}_N = \mathbb{P}^2(\mathcal{T})^d \cap H_0^1(\Omega)^d \quad \text{and} \quad M_N = \mathbb{P}^1(\mathcal{T}) \cap L_0^2(\Omega). \qquad (52) $$

## 4.3 The $\mathbb{P}^2(\mathcal{T})$ Lagrange interpolator

The Lagrange interpolator, denoted by $\Pi_L : C_c^\infty(\Omega)^d \to \mathbf{X}_N$, is defined by

$$ \forall \mathbf{v} = (v_k)_{k=1,...,d} \in C_c^\infty(\Omega)^d, \quad (\Pi_L \mathbf{v})_k = \sum_{i \in \mathcal{N} \cup \mathcal{E}} v_k(\mathbf{y}_i)\phi_i, \quad k = 1, ..., d. \qquad (53) $$

**Lemma 4.3:** There exists $C_6$ such that for any $\mathbf{v} \in C_c^\infty(\Omega)^d$, we have

$$ \|\Pi_L \mathbf{v}\|_{L^\infty(\Omega)^d} \le C_6 \|\nabla \mathbf{v}\|_{L^\infty(\Omega)^{d \times d}}. \qquad (54) $$

**PROOF.** Let $\underline{x}_1, \bar{x}_1 \in \mathbb{R}$ such that, for any $\mathbf{x} = (x_1, ..., x_d) \in \Omega$, $\underline{x}_1 \le x_1 \le \bar{x}_1$. Writing, for $k = 1, ..., d$,

$$ v_k(\mathbf{y}_i) = \int_{\underline{x}_1}^{(y_i)_1} \partial_1 v_k(s, ((y_i)_\ell)_{\ell=2,...,d})ds, $$

we get

$$ |v_k(\mathbf{y}_i)| \le (\bar{x}_1 - \underline{x}_1) \|\nabla \mathbf{v}\|_{L^\infty(\Omega)^{d \times d}}. $$

Since $|\phi_i(\mathbf{x})| \le 1$ for all $\mathbf{x} \in \Omega$, and using that the maximum number of $i \in \mathcal{N} \cup \mathcal{E}$ such that $\phi_i(\mathbf{x}) \ne 0$ is $(d+1)(d+2)/2$, we get from (53) that

$$ |v_k(\mathbf{x})| \le \frac{(d+1)(d+2)}{2} (\bar{x}_1 - \underline{x}_1) \|\nabla \mathbf{v}\|_{L^\infty(\Omega)^{d \times d}}, $$

hence concluding the proof. □

We have the following results.

**Lemma 4.4:** There exists $C_7$ such that for any $\mathbf{v} = (v_k)_{k=1,...,d} \in C_c^\infty(\Omega)^d$, it holds, for any $K \in \mathcal{T}$,

$$ \|\mathbf{v} - \Pi_L \mathbf{v}\|_{L^\infty(K)^d} \le C_7 h_K \|\nabla \mathbf{v}\|_{L^\infty(\Omega)^{d \times d}}, \qquad (55) $$

$$ \|\nabla \Pi_L \mathbf{v}\|_{L^\infty(K)^{d \times d}} \le C_7 \|\nabla \mathbf{v}\|_{L^\infty(\Omega)^{d \times d}}, \qquad (56) $$

$$ \|\mathbf{v} - \Pi_L \mathbf{v}\|_{L^\infty(K)^d} \le C_7 h_K^2 \|\mathbf{v}\|_{W^{2,\infty}(\Omega)}, \qquad (57) $$

and

$$ \|\nabla \mathbf{v} - \nabla \Pi_L \mathbf{v}\|_{L^\infty(K)^{d \times d}} \le C_7 h_K \|\mathbf{v}\|_{W^{2,\infty}(\Omega)}. \qquad (58) $$
PROOF. Let $\mathbf{v} = (v_k)_{k=1,...,d} \in C_c^\infty(\Omega)^d$. For $k = 1, ..., d$, the mean value theorem implies that there exists $\mathbf{x}_i$ belonging to the segment $[\mathbf{y}_i, \mathbf{x}]$ such that

$$v_k(\mathbf{y}_i) - v_k(\mathbf{x}) = (\mathbf{y}_i - \mathbf{x}) \cdot \nabla v_k(\mathbf{x}_i),$$

which implies

$$|v_k(\mathbf{y}_i) - v_k(\mathbf{x})| \le h_K \|\nabla \mathbf{v}\|_{L^\infty(\Omega)^{d \times d}}.$$

Noticing that, for any $K \in \mathcal{T}$ and for all $\mathbf{x} \in K$, we have $\sum_{i \in \mathcal{N}_K \cup \mathcal{E}_K} \phi_i(\mathbf{x}) = 1$, we get

$$|\Pi_L v_k(\mathbf{x}) - v_k(\mathbf{x})| = \left| \sum_{i \in \mathcal{N}_K \cup \mathcal{E}_K} (v_k(\mathbf{y}_i) - v_k(\mathbf{x})) \phi_i(\mathbf{x}) \right| \le \frac{(d+1)(d+1)}{2} h_K \|\nabla \mathbf{v}\|_{L^\infty(\Omega)^{d \times d}}.$$

This yields (55) as well as

$$\nabla v_k(\mathbf{x}) = \sum_{i \in \mathcal{N}_K \cup \mathcal{E}_K} (v_k(\mathbf{y}_i) - v_k(\mathbf{x})) \nabla \phi_i(\mathbf{x}) = \sum_{i \in \mathcal{N}_K \cup \mathcal{E}_K} (\mathbf{y}_i - \mathbf{x}) \cdot \nabla v_k(\mathbf{x}_i) \nabla \phi_i(\mathbf{x}).$$

Applying (51), we conclude that (56) holds. Writing, for $k = 1, ..., d$ and any $j \in \mathcal{N}_K \cup \mathcal{E}_K$, the Taylor expansion

$$v_k(\mathbf{y}_j) = v_k(\mathbf{x}) + (\mathbf{y}_j - \mathbf{x}) \cdot \nabla v_k(\mathbf{x}) + \frac{1}{2}(\mathbf{y}_j - \mathbf{x})^t D_2 v_k(\mathbf{x}_j)(\mathbf{y}_j - \mathbf{x}),$$

for some point $\mathbf{x}_j$ in the segment $[\mathbf{y}_j, \mathbf{x}]$, and using (53) in addition to

$$\sum_{j \in \mathcal{N}_K \cup \mathcal{E}_K} \phi_j(\mathbf{x}) = 1 \text{ and } \sum_{j \in \mathcal{N}_K \cup \mathcal{E}_K} \phi_j(\mathbf{x})(\mathbf{y}_j - \mathbf{x}) = 0,$$

we obtain (57). We deduce, from the preceding Taylor expansion, that

$$\sum_{j \in \mathcal{N}_K \cup \mathcal{E}_K} v_k(\mathbf{y}_j) \nabla \phi_i(\mathbf{x}) = \sum_{j \in \mathcal{N}_K \cup \mathcal{E}_K} v_k(\mathbf{x}) \nabla \phi_i(\mathbf{x}) + \sum_{j \in \mathcal{N}_K \cup \mathcal{E}_K} (\mathbf{y}_j - \mathbf{x}) \cdot \nabla v_k(\mathbf{x}) \nabla \phi_i(\mathbf{x}) \\ + \sum_{j \in \mathcal{N}_K \cup \mathcal{E}_K} \frac{1}{2} (\mathbf{y}_j - \mathbf{x})^t D_2 v_k(\mathbf{x}_j) (\mathbf{y}_j - \mathbf{x}) \nabla \phi_i(\mathbf{x}).$$

We notice that, since the interpolation of first order polynomials is exact,

$$\sum_{j \in \mathcal{N}_K \cup \mathcal{E}_K} v_k(\mathbf{x}) \nabla \phi_i(\mathbf{x}) = 0 \text{ and } \sum_{j \in \mathcal{N}_K \cup \mathcal{E}_K} (\mathbf{y}_j - \mathbf{x}) \cdot \nabla v_k(\mathbf{x}) \nabla \phi_i(\mathbf{x}) = \nabla v_k(\mathbf{x}),$$

and

$$\left| \sum_{j \in \mathcal{N}_K \cup \mathcal{E}_K} \frac{1}{2} (\mathbf{y}_j - \mathbf{x})^t D_2 v_k(\mathbf{x}_j) (\mathbf{y}_j - \mathbf{x}) \nabla \phi_i(\mathbf{x}) \right| \le C_8 h_K \|\mathbf{v}\|_{W^{2,\infty}(\Omega)},$$

which allows to conclude (58). □

## 4.4 Divergence correcting operator

As stated in the introduction of Section 4, we adapt and simplify the edge-bubble operator defined in [11], without requesting that in the general case the divergence correcting operator vanishes at the boundary of the domain.
### 4.4.1 Tangential edge bubble functions

We define the family of normalized tangential edge bubble function

$$ \forall \mathbf{x} \in \Omega, \mathbf{b}_{i,j}(\mathbf{x}) = \frac{(d+2)(d+1)}{|w_{i,j}|} \varphi_i(\mathbf{x}) \varphi_j(\mathbf{x}) (\mathbf{y}_j - \mathbf{y}_i) \text{ for any } \{i,j\} \in \mathcal{E}. $$

Note that, using (49) we have

$$ \{\mathbf{x} \in \Omega \text{ such that } \varphi_i(\mathbf{x})\varphi_j(\mathbf{x}) \neq 0\} = \omega_{ij} \text{ for any } \{i,j\} \in \mathcal{E}. \quad (59) $$

We observe that, for any $\{i, j\} \in \mathcal{E}$, we have $\mathbf{b}_{i,j} \in \mathbb{P}^2(\mathcal{T})$, and, in the case where the segment $]y_i, y_j[ \subset \Omega$, we additionally have $\mathbf{b}_{i,j} \in \mathbf{X}_N$ (this property holds even if $\{i, j\} \in \mathcal{E}_{\text{ext}}$).

**Lemma 4.5:** For any $\{i, j\} \in \mathcal{E}$ and for any $k \in \mathcal{N}$ we have the following relation

$$ (\operatorname{div}\mathbf{b}_{i,j}, \varphi_k) = \delta_{i,k} - \delta_{j,k}. $$

**PROOF.** Let $\{i, j\} \in \mathcal{E}$ and let $k \in \mathcal{N}$. Using (59) We have

$$ (\operatorname{div}\mathbf{b}_{i,j}, \varphi_k) = \int_{\omega_{ij}} \operatorname{div}\mathbf{b}_{i,j}\varphi_k \, d\mathbf{x}. $$

* The element $k$ is not an element of $\cup_{K \in \mathcal{T}_{ij}} \mathcal{N}_K$. Using (49) we have $\varphi_k(\mathbf{x}) = 0$ for any $\mathbf{x} \in \omega_{ij}$ and hence the expression vanishes.

* The element $k$ is an element of $\cup_{K \in \mathcal{T}_{ij}} \mathcal{N}_K$. We have

$$ (\operatorname{div}\mathbf{b}_{i,j}, \varphi_k) = -\frac{(d+2)(d+1)}{|w_{i,j}|} \sum_{K \in \mathcal{T}_{ij}} \int_K \varphi_i \varphi_j (\mathbf{y}_j - \mathbf{y}_i) \cdot \nabla \varphi_k \, d\mathbf{x} \\ + \frac{(d+2)(d+1)}{|w_{i,j}|} \sum_{K \in \mathcal{T}_{ij}} \sum_{F \in \mathcal{F}_K} \int_F \varphi_k \varphi_i \varphi_j (\mathbf{y}_j - \mathbf{y}_i) \cdot \mathbf{n}_F \, d\gamma(\mathbf{x}) $$

For any $K \in \mathcal{T}_{ij}$ and $F \in \mathcal{F}_K$ we remark that $(\mathbf{y}_j - \mathbf{y}_i) \cdot \mathbf{n}_F = 0$ in the case $\{i, j\} \subset \mathcal{N}_F$ and $\varphi_i(\mathbf{x})\varphi_j(\mathbf{x}) = 0$ for any $\mathbf{x} \in F$ in the case $\{i, j\} \not\subset \mathcal{N}_F$.

We then obtain for any $K \in \mathcal{T}_{ij}$ and for any $F \in \mathcal{F}_K$

$$ \int_F \varphi_k \varphi_i \varphi_j (\mathbf{y}_j - \mathbf{y}_i) \cdot \mathbf{n}_F \, d\gamma(\mathbf{x}) = 0 $$

which gives the following identity

$$ (\operatorname{div}\mathbf{b}_{i,j}, \varphi_k) = -\frac{(d+2)(d+1)}{|w_{i,j}|} \sum_{K \in \mathcal{T}_{ij}} \int_K \varphi_i \varphi_j (\mathbf{y}_j - \mathbf{y}_i) \cdot \nabla \varphi_k \, d\mathbf{x}. $$

For $k \notin \{i, j\}$ using $\varphi_k(\mathbf{y}_i) = 0$ and $\varphi_k(\mathbf{y}_j) = 0$ we obtain that the piecewise constant function $\nabla \varphi_k$ is orthogonal to the edge vector $(\mathbf{y}_j - \mathbf{y}_i)$ on each $K \in \mathcal{T}_{ij}$. Thus, the integrals vanishes and the claim follows. For $k \in \{i, j\}$ using $\varphi_k(\mathbf{y}_k) = 1$ we obtain for any $K \in \mathcal{T}_{ij}$ the following identity

$$ \nabla(\varphi_k)|_K \cdot (\mathbf{y}_j - \mathbf{y}_i) = \begin{cases} -1, & \text{for } k = i, \\ 1, & \text{for } k = j. \end{cases} $$
For $k = i$ we obtain

$$ (\operatorname{div}\mathbf{b}_{i,j}, \varphi_k) = \frac{(d+2)(d+1)}{|w_{i,j}|} \sum_{K \in \mathcal{T}_{ij}} \int_K \varphi_i \varphi_j \, d\mathbf{x} = 1, $$

and for $k = j$ we obtain

$$ (\operatorname{div}\mathbf{b}_{i,j}, \varphi_k) = -\frac{(d+2)(d+1)}{|w_{i,j}|} \sum_{K \in \mathcal{T}_{ij}} \int_K \varphi_i \varphi_j \, d\mathbf{x} = -1, $$

which gives the expected result.

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

**Lemma 4.6:** There exists $C_9$ such that

$$ \forall K \in \mathcal{T}, \forall \{i, j\} \in \mathcal{E}_K, \|\mathbf{b}_{ij}\|_{L^{\infty}(K)^d} \le C_9 h_K^{1-d}, \quad (60) $$

and

$$ \forall K \in \mathcal{T}, \forall \{i, j\} \in \mathcal{E}_K, \|\nabla \mathbf{b}_{ij}\|_{L^{\infty}(K)^d} \le C_9 h_K^{-d}. \quad (61) $$

**PROOF.** We have, for $\{i, j\} \subset \mathcal{N}_K$,

$$ |\varphi_i(\mathbf{x})\varphi_j(\mathbf{x})| \le 1 $$

and

$$ |w_{i,j}| \ge C_{10} h_K^d, $$

which leads to (60) since $|\mathbf{y}_j - \mathbf{y}_i| \le h_K$. From (50), we get

$$ |\nabla(\varphi_i \varphi_j)(\mathbf{x})| \le 2 \frac{C_4}{h_K}. $$

Combining with the preceding inequalities, we deduce (61). □

### 4.4.2 Divergence correcting operator

We can now define the divergence correcting operator $\Pi_D : W_0^{1,\infty}(\Omega)^d \to \mathbb{P}^2(\mathcal{T})^d$ for any $\mathbf{v} \in W_0^{1,\infty}(\Omega)^d$ by

$$ \forall \mathbf{x} \in \Omega, \Pi_D \mathbf{v}(\mathbf{x}) = \sum_{i \in \mathcal{N}} \Pi_{D,i} \mathbf{v}(\mathbf{x}) \quad \text{with} \quad \Pi_{D,i} \mathbf{v}(\mathbf{x}) = \sum_{\substack{j \in \mathcal{N} \\ \{i,j\} \in \mathcal{E}}} (\operatorname{div}(\varphi_i \mathbf{v}), \varphi_j) \mathbf{b}_{j,i}(\mathbf{x}). \quad (62) $$

Using the fact that $\mathbf{b}_{i,j} = -\mathbf{b}_{j,i}$ for any $\{i, j\} \in \mathcal{E}$ we remark that this operator satisfies

$$ \forall \mathbf{v} \in W_0^{1,\infty}(\Omega)^d, \forall \mathbf{x} \in \Omega, \Pi_D \mathbf{v}(\mathbf{x}) = \sum_{\{i,j\} \in \mathcal{E}} (\mathbf{v}, -\varphi_i \nabla \varphi_j + \varphi_j \nabla \varphi_i) \mathbf{b}_{j,i}(\mathbf{x}). \quad (63) $$

**Lemma 4.7:** For any $\mathbf{v} \in W_0^{1,\infty}(\Omega)^d$ and for any $\mathbf{x} \in \Omega$, if $\Pi_D \mathbf{v}(\mathbf{x}) \ne 0$, then there exists $\mathbf{y} \in B(\mathbf{x}, C_3 h_{\mathcal{T}})$ such that $\mathbf{v}(\mathbf{y}) \ne 0$ (recall that $C_3$ is introduced in Lemma 4.1).
PROOF. If $\Pi_D v(x) \neq 0$, there exists $\{i, j\} \in \mathcal{E}$ such that $(v, -\varphi_i \nabla \varphi_j + \varphi_j \nabla \varphi_i) b_{j,i}(x) \neq 0$, which implies $x \in \omega_{ij}$, which yields that there exists $K \in T_{ij}$ such that $x \in \bar{K}$. This also implies that

$$ \int_{\omega_{ij}} |v(y)| dy > 0. $$

Therefore there exists $y \in \omega_{ij}$ such that $v(y) \neq 0$. Remarking that (48) implies that $\text{diam}(\omega_{ij}) \leq C_3 h_K$. Hence $|x - y| \leq C_3 h_K$. □

We can now prove that this operator preserves the divergence in the following sense.

Lemma 4.8: For any $v \in W_0^{1,\infty}(\Omega)^d$ and any $\varphi \in \mathbb{P}^1(\mathcal{T})$ we have

$$ (\text{div} \Pi_D v, \varphi) = (\text{div} v, \varphi). $$

PROOF. Using Lemma 4.5 we have, for any $k \in \mathcal{N}$,

$$ \forall i \in \mathcal{N}, (\text{div} \Pi_{D,i} v, \varphi_k) = \sum_{\substack{j \in \mathcal{N} \\ \{i,j\} \in \mathcal{E}}} (\text{div}(\varphi_i v), \varphi_j)(\delta_{j,k} - \delta_{i,k}). $$

Hence we get $(\text{div} \Pi_D v, \varphi_k) = T_1 + T_2$, with

$$ T_1 = \sum_{i \in \mathcal{N}} \sum_{\{j, i\} \in \mathcal{N}} (\text{div}(\varphi_i v), \varphi_j) \delta_{j,k} \quad \text{and} \quad T_2 = - \sum_{i \in \mathcal{N}} \sum_{\{j, i\} \in \mathcal{N}} (\text{div}(\varphi_i v), \varphi_j) \delta_{i,k}. $$

Remarking that, for any quantity $Q_{i,j}$,

$$ \sum_{i \in \mathcal{N}} \sum_{\{j, i\} \in \mathcal{N}} Q_{i,j} = \sum_{\{i, j\} \in \mathcal{E}} (Q_{i,j} + Q_{j,i}) = \sum_{j \in \mathcal{N}} \sum_{\{i, j\} \in \mathcal{E}} Q_{i,j}, $$

we obtain

$$ T_1 = \sum_{j \in \mathcal{N}} \sum_{\{i, j\} \in \mathcal{E}} (\text{div}(\varphi_i v), \varphi_j) \delta_{j,k} = (\text{div}(\sum_{\{i, k\} \in \mathcal{E}} \varphi_i v), \varphi_k). $$

For any $x \in \omega_k$, we have

$$ \sum_{\substack{i \in \mathcal{N} \\ \{i, k\} \in \mathcal{E}}} \varphi_i(x) = 1 - \varphi_k(x), $$

which leads to

$$ T_1 = (\text{div}((1 - \varphi_k)v), \varphi_k). $$

Next we have

$$ T_2 = - \sum_{\substack{j \in \mathcal{N} \\ \{k, j\} \in \mathcal{E}}} (\text{div}(\varphi_k v), \varphi_j) = -(\text{div}(\varphi_k v), \sum_{\substack{j \in \mathcal{N} \\ \{k, j\} \in \mathcal{E}}} \varphi_j). $$

For any $x \in \omega_k$, we have

$$ \sum_{\substack{j \in \mathcal{N} \\ \{k, j\} \in \mathcal{E}}} \varphi_j(x) = 1 - \varphi_k(x), $$

which implies, since $(\text{div}(\varphi_k v), 1) = 0$ (recall that $v$ vanishes on the boundary of $\Omega$),

$$ T_2 = -(\text{div}(\varphi_k v), -\varphi_k) = (\text{div}(\varphi_k v), \varphi_k). $$

This yields

$$ T_1 + T_2 = (\text{div}(v), \varphi_k), $$

which concludes the proof. □
**Lemma 4.9:** There exists $C_{11}$ such that, for any $v \in W_0^{1,\infty}(\Omega)^d$, we have

$$ \|\Pi_D v\|_{L^\infty(\Omega)^d} \le C_{11} \|v\|_{L^\infty(\Omega)^d}. $$

**PROOF.** Using (63) and the definition of the tangential bubble functions we have, for $x \in K$ with $K \in \mathcal{T}$,

$$ |\Pi_D v(x)| \le \sum_{\{i,j\} \in \mathcal{E}} |(v, -\varphi_i \nabla \varphi_j + \varphi_j \nabla \varphi_i)| \|b_{i,j}(x)\|. $$

If $K \subset \omega_{ij}$, then there exists $C_{12}$ such that $\omega_{ij} \subset B(x, C_{12}h_K)$. Since $m(\omega_{ij}) \ge C_{13}h_K^d$, we get the existence of $C_{14}$ such that the number of $\{i, j\} \in \mathcal{E}$ with $b_{i,j}(x) \ne 0$ is bounded by $C_{14}$. Using (60) and (50), we obtain

$$ |\Pi_D v(x)| \le \|v\|_{L^\infty(\Omega)^d} C_{14} C_9 h_K \frac{C_4}{h_K}, $$

which concludes the proof. □

**Lemma 4.10:** There exists $C_{15}$ such that

$$ \forall v \in W_0^{1,\infty}(\Omega)^d, \quad \|\nabla \Pi_D v\|_{L^\infty(\Omega)^d \times \mathcal{T}} \le C_{15} \max_{K \in \mathcal{T}} \frac{\|v\|_{L^\infty(K)^d}}{h_K}. \quad (64) $$

**PROOF.** Let us denote, for a given $K \in \mathcal{T}$,

$$ C_v = \max_{K \in \mathcal{T}} \frac{\|v\|_{L^\infty(K)^d}}{h_K}. \quad (65) $$

Using (63), we have, for $x \in K$,

$$ \nabla \Pi_D v(x) = \sum_{\{i,j\} \in \mathcal{E}} (v, -\varphi_i \nabla \varphi_j + \varphi_j \nabla \varphi_i) \nabla b_{i,j}(x). $$

We remark that, using (65) and (50) for $K \subset \omega_{ij}$ and $x \in K$,

$$ |v(x)(-\varphi_i(x)\nabla\varphi_j(x) + \varphi_j(x)\nabla\varphi_i(x))| \le \|v\|_{L^\infty(K)^d} 2\frac{C_4}{h_K}. $$

This yields

$$ |(v, -\varphi_i \nabla \varphi_j + \varphi_j \nabla \varphi_i)| \le 2|\omega_{ij}|C_v C_4, $$

and therefore, using (61), we conclude (64). □

## 4.5 The operator $\Pi_N$

We let $W = C_c^\infty(\Omega)^d \cap V(\Omega)$. We notice that, for any $v \in W$, then $\Pi_L v \in W_0^{1,\infty}(\Omega)^d$. One can then consider the function $\Pi_D(v - \Pi_L v) \in \mathbb{P}^2(\mathcal{T})^d$. We notice that this function is not necessarily null on the boundary of the domain, and in this case, it is not an element of $X_N$.

We then define $\Pi_N: W \to X_N$, for any $v \in W$, by:

$$ \text{if } \Pi_D(v - \Pi_L v) \in H_0^1(\Omega)^d, \text{ then } \Pi_N v = \Pi_L v + \Pi_D(v - \Pi_L v) \text{ else } \Pi_N v = 0. \quad (66) $$

Notice that the operator $\Pi_N$ defined by (66) is not linear.

**Lemma 4.11:** For any $v \in C_c^\infty(\Omega)^d \cap V(\Omega)$ we have $\Pi_N v \in E_N$ that is

$$ \forall \varphi \in \mathbb{P}^1(\mathcal{T}), (\text{div}\Pi_N v, \varphi) = 0. \quad (67) $$
PROOF. If $\Pi_D(v - \Pi_L v) \in H_0^1(\Omega)^d$, we get from Lemma 4.8

$$(\operatorname{div}\Pi_D(v - \Pi_L v), \varphi) = (\operatorname{div}(v - \Pi_L v), \varphi) = -(\operatorname{div}\Pi_L v, \varphi),$$

which provides (67). Otherwise, we get (67) since $\Pi_N v = 0$. □

**Lemma 4.12:** There exists $C_{16}$ such that for any $v \in W$, we have

$$\|\Pi_N v\|_{L^\infty(\Omega)^d} \le C_{16} \|\nabla v\|_{L^\infty(\Omega)^{d \times d}}. \qquad (68)$$

and

$$\|\nabla \Pi_N v\|_{L^\infty(\Omega)^d} \le C_{16} \|\nabla v\|_{L^\infty(\Omega)^{d \times d}}. \qquad (69)$$

PROOF. It suffices to consider the case $\Pi_D(v - \Pi_L v) \in H_0^1(\Omega)^d$ (otherwise the conclusion is straightforward). We apply Lemmas 4.3 and 4.9. We then get (68).

Now applying Lemmas (55) in Lemma 4.4 and Lemma 4.10, we obtain (69). □

**Lemma 4.13:** There exists $C_{17}$ such that, for any $v \in W$ such that $\Pi_D(v - \Pi_L v) \in H_0^1(\Omega)$, then

$$\|\Pi_N v - v\|_{W^{1,\infty}(\Omega)^d} \le C_{17} h_{\mathcal{T}} \|v\|_{W^{2,\infty}(\Omega)}. \qquad (70)$$

PROOF. We first remark that, from (57) and (58) in Lemma 4.4, we get

$$\|\Pi_L v - v\|_{W^{1,\infty}(\Omega)^d} \le C_{18} h_{\mathcal{T}} \|v\|_{W^{2,\infty}(\Omega)}. $$

Again applying (57) in the same Lemma, and using Lemma 4.10, we get that

$$\|\Pi_D(v - \Pi_L v)\|_{W^{1,\infty}(\Omega)^d} \le C_{15} C_7 h_{\mathcal{T}} \|v\|_{W^{2,\infty}(\Omega)}, $$

which concludes the proof of (70). □

**Theorem 4.14:** Let $W = C_c^\infty(\Omega) \cap V(\Omega)$. Let $(\mathcal{T}_N)_{N \ge 1}$ be a sequence of simplicial meshes in the sense of this section, such that $(h_{\mathcal{T}_N})_{N \ge 1}$ converges to zero while $(\theta_{\mathcal{T}_N})_{N \ge 1}$ remains bounded. Let $(X_N, M_N, \Pi_N)_{N \ge 1}$ be such that $X_N$ and $M_N$ are the Taylor-Hood spaces defined by (52), $\Pi_N$ is defined by (66), in the case where $\mathcal{T} = \mathcal{T}_N$.

Then the family $(W, (X_N, M_N, \Pi_N)_{N \ge 1})$ satisfies the convergence assumptions (21)-(25) in Section 3.1.

PROOF. The density of $C_c^\infty(\Omega)^d \cap V(\Omega)$ in $H_0^1(\Omega)^d \cap V(\Omega)$ for the $H_0^1(\Omega)^d$ norm is proved in particular in [4, Lemma IV.3.4 p. 249].

Let $v \in W$. Applying (66), we have that $\Pi_N v \in E_N$.

Indeed, if $\Pi_D(v - \Pi_L v) \in H_0^1(\Omega)$, then $\Pi_D(v - \Pi_L v) \in X_N$ as well as $\Pi_L v \in X_N$. In this case, the application of Lemma 4.11 yields (21).

Otherwise, if $\Pi_D(v - \Pi_L v) \notin H_0^1(\Omega)$, then $\Pi_N v = 0 \in E_N$ and (21) is also satisfied.

Let $v \in C_c^\infty(\Omega) \cap V(\Omega)$. Since the support of $v$ is compact, there exists $a > 0$ such that, for any $x \in \partial\Omega$, $v = 0$ on the ball $B(x, a)$. Let $N_0 \ge 1$ be such that, for all $N \ge N_0$, $h_{\mathcal{T}_N} \le a/(C_3 + 2)$. Let $N \ge N_0$, and let us consider $\Pi_L$ defined by (53) in the case where $\mathcal{T} = \mathcal{T}_N$. Then, for any $x \in \partial\Omega$, $v - \Pi_L v = 0$ on the ball $B(x, a - h_{\mathcal{T}_N})$.

We then deduce from Lemma 4.7 that, for any $x \in \Omega$ with $B(x, h_{\mathcal{T}_N}) \cap \partial\Omega \ne \emptyset$, then $v - \Pi_L v = 0$ on the ball $B(x, C_3 h_{\mathcal{T}_N})$, which leads to $\Pi_D(v - \Pi_L v) \in H_0^1(\Omega)$ for all $N \ge N_0$ in the case where $\mathcal{T} = \mathcal{T}_N$. Therefore, $\Pi_D(v - \Pi_L v) \in H_0^1(\Omega)$ for all $N \ge N_0$ in the case where $\mathcal{T} = \mathcal{T}_N$. This yields that, for any $N \ge N_0$, the conclusion of Lemma 4.13 holds. This proves (24).

Finally, (25) is a consequence of the convergence of $\mathbb{P}^1$ finite elements for the Neumann problem with null average. □
**Remark 4.1 (A pair of finite element which does not fulfill the inf-sup condition):**

Let $\omega > 1$ be a given real and let $(\mathcal{T}_N)_{N \ge 1}$ be a sequence of simplicial meshes in the sense of this section, such that $(h_{\mathcal{T}_N})_{N \ge 1}$ converges to zero while $(\theta_{\mathcal{T}_N})_{N \ge 1}$ remains bounded. For any $N \ge 1$, let $\tilde{\mathcal{T}}_N$ be the set of the elements of $\mathcal{T}_N$ whose distance to the boundary is greater than $\omega h_{\mathcal{T}_N}$. We then define $\Omega_N = \bigcup_{K \in \tilde{\mathcal{T}}_N} \bar{K}$. We consider $\mathbb{P}^1(\mathcal{T})$ finite element for the pressure, $\mathbb{P}^2(\mathcal{T})$ finite element for the velocity on

any $K \in \tilde{\mathcal{T}}_N$, and $\mathbb{P}^1(\mathcal{T})$ or transition finite element as in Figure 3 for any $K \in \mathcal{T}_N \setminus \tilde{\mathcal{T}}_N$. We then define $\Pi_D$ only for functions whose support is included in $\Omega_N$ and we define $\Pi_N : \mathbf{W} \to \mathbf{E}_N = \mathbf{X}_N \cap \mathbf{V}_N \cap L^\infty(\Omega)^d$, for any $\mathbf{v} \in \mathbf{W}$, by:

if $\text{support}(\mathbf{v} - \Pi_L \mathbf{v}) \subset \Omega_N$, then $\Pi_N \mathbf{v} = \Pi_L \mathbf{v} + \Pi_D(\mathbf{v} - \Pi_L \mathbf{v})$ else $\Pi_N \mathbf{v} = 0$.

Following the proof of Theorem 4.14, we can prove that this operator $\Pi_N$ fulfills the convergence assumptions (21)-(25) in Section 3.1.

<figure><img src="image_7.png" /><figcaption>Fig. 3: Finite element pair which does not fulfill the inf-sup condition. The degrees of freedom for the components of the velocity are in black, those for the pressure are in pink. The two hatched elements include 4 degrees of freedom for each of the components of the velocity.</figcaption></figure>

# 5 Has the function $\bar{u}$ provided by Theorem 3.8 some additional regularity?

Theorem 3.8 shows the convergence of the incremental projection scheme to a weak solution $\bar{u}$ in the sense of Definition 1.1. The question that $\bar{u}$ be suitable in the sense introduced by Caffarelli-Kohn-Nirenberg in [7] arises. In [21] and later in [15], the authors proved that the weak solution, obtained as the limit of a semi-discrete or fully discrete numerical scheme coupling the velocity and the pressure under the existence of an inf-sup condition, is a suitable weak solution. Further works have to be done in order to determine whether this property also holds for $\bar{u}$ obtained as the limit of the projection method.

# References

[1] M. Bercovier and O. Pironneau. Error estimates for finite element method solution of the Stokes problem in the primitive variables. *Numerische Mathematik*, 33:211–224, 1979.

[2] D. Boffi. Stability of higher order triangular hood-taylor methods for the stationary stokes equations. *Mathematical Models and Methods in Applied Sciences*, 04:223–235, 1994.
[3] J. Boland and R. A. Nicolaides. Stable and semistable low order finite elements for viscous flows. *SIAM Journal on Numerical Analysis*, 22:474–492, 1985.

[4] F. Boyer and P. Fabrie. *Mathematical Tools for the Study of the Incompressible Navier-Stokes Equations and Related Models*. Applied Mathematical Sciences 183. Springer, 2013.

[5] F. Brezzi and R. S. Falk. Stability of higher-order Hood-Taylor methods. *SIAM Journal on Numerical Analysis*, 28:581–590, 1991.

[6] F. Brezzi and M. Fortin. *Mixed and hybrid finite element methods*, volume 15 of Springer Series in Computational Mathematics. Springer-Verlag, New York, 1991.

[7] L. Caffarelli, R. Kohn, and L. Nirenberg. Partial regularity of suitable weak solutions of the Navier-Stokes equations. *Comm. Pure Appl. Math.*, 35(6):771–831, 1982.

[8] L. Chen. A simple construction of a Fortin operator for the two dimensional Taylor-Hood element. *Comput. Math. Appl.*, 68:1368–1373, 2014.

[9] A. J. Chorin. On the convergence of discrete approximations to the Navier-Stokes equations. 1969.

[10] P. G. Ciarlet. The finite element method for elliptic problems. In *Studies in Mathematics and its Applications*, Vol. 4, pages xix+530. North-Holland Publishing Co., Amsterdam-New York-Oxford, 1978.

[11] L. Diening, J. Storn, and T. Tscherpel. Fortin operator for the Taylor-Hood element. *Numerische Mathematik*, 150:671–689, 2021.

[12] A. Ern and J.-L. Guermond. Finite element quasi-interpolation and best approximation. *ESAIM Math. Model. Numer. Anal.*, 51(4):1367–1385, 2017.

[13] A. Ern and J.-L. Guermond. Finite elements i. *Texts in Applied Mathematics*, 2021.

[14] R. S. Falk. A Fortin operator for two-dimensional Taylor-Hood elements. *Mathematical Modelling and Numerical Analysis*, 42:411–424, 2008.

[15] E. Fernández-Cara and I. Marín-Gayte. A new proof of the existence of suitable weak solutions and other remarks for the Navier-Stokes equations. *Applied Mathematics*, 9:383–402, 2018.

[16] T. Gallouët and R. Herbin. Equations aux dérivées partielles. *HAL SCIENCE*, (cel-01196782), 2015.

[17] T. Gallouët, R. Herbin, and J.-C. Latché. $W^{1,q}$ stability of the Fortin operator for the MAC scheme. *Calcolo*, 49(1):63–71, 2012.

[18] T. Gallouët, R. Herbin, J.-C. Latché, and D. Maltese. Convergence of the fully discrete incremental projection scheme for incompressible flows. *Journal of Mathematical Fluid Mechanics*, 25(63), 2023.

[19] K. Goda. A multistep technique with implicit difference schemes for calculating two- or three-dimensional cavity flows. *Journal of Computational Physics*, 30:76–95, 1979.

[20] J.-L. Guermond. Some implementations of projection methods for Navier-Stokes equations. *Mathematical Modelling and Numerical Analysis*, 30:637–667, 1996.

[21] J.-L. Guermond. Faedo-Galerkin weak solutions of the Navier-Stokes equations with Dirichlet boundary conditions are suitable. *J. Math. Pures Appl. (9)*, 88(1):87–106, 2007.

[22] J.-L. Guermond, P. D. Minev, and J. Shen. An overview of projection methods for incompressible flows. *Computer Methods in Applied Mechanics and Engineering*, 195:6011–6045, 2006.
[23] J.-L. Guermond and L. Quartapelle. Calculation of incompressible viscous flows by an unconditionally stable projection fem. *Journal of Computational Physics*, 132:12–33, 1997.

[24] J.-L. Guermond and L. Quartapelle. On stability and convergence of projection methods based on pressure Poisson equation. *International Journal for Numerical Methods in Fluids*, 26:1039–1053, 1998.

[25] J.-L. Guermond and L. Quartapelle. On the approximation of the unsteady Navier-Stokes equations by finite element projection methods. *Numerische Mathematik*, 80:207–238, 1998.

[26] F. H. Harlow and J. E. Welch. Numerical calculation of time-dependent viscous incompressible flow of fluid with free surface. *Physics of Fluids*, 8:2182–2189, 1965.

[27] H. Kuroki and K. Soga. On convergence of Chorin's projection method to a Leray-Hopf weak solution. *Numerische Mathematik*, 146:401–433, 2020.

[28] J. Leray. Essai sur les mouvements plans d'un fluide visqueux que limitent des parois. *Journal de Mathématiques Pures et Appliquées*, 13:331–418, 1934. URL: http://eudml.org/doc/235388.

[29] J.-L. Lions. *Quelques méthodes de résolution des problèmes aux limites non linéaires*. Dunod, Paris; Gauthier-Villars, Paris, 1969.

[30] K. Mardal, J. Schöberl, and R. Winther. A uniformly stable Fortin operator for the Taylor-Hood element. *Numerische Mathematik*, 123:537–551, 2013.

[31] R. Rannacher. On Chorin's projection method for the incompressible Navier-Stokes equations. 1992.

[32] L. R. Scott and S. Zhang. Finite element interpolation of nonsmooth functions satisfying boundary conditions. *Mathematics of Computation*, 54:483–493, 1990.

[33] J. Shen. On error estimates of projection methods for Navier-Stokes equations: first-order schemes. *SIAM Journal on Numerical Analysis*, 29:57–77, 1992.

[34] J. Shen. Remarks on the pressure error estimates for the projection methods. *Numerische Mathematik*, 67:513–520, 1994.

[35] H. Sohr. *The Navier-Stokes equations*. Modern Birkhäuser Classics. Birkhäuser/Springer Basel AG, Basel, 2001. An elementary functional analytic approach, [2013 reprint of the 2001 original] [MR1928881].

[36] R. Stenberg. Analysis of mixed finite elements methods for the Stokes problem: a unified approach. *Mathematics of Computation*, 42:9–23, 1984.

[37] R. Temam. Sur l'approximation de la solution des équations de Navier-Stokes par la méthode des pas fractionnaires (ii). *Archive for Rational Mechanics and Analysis*, 33:377–385, 1969.

[38] J. van Kan. A second-order accurate pressure correction scheme for viscous incompressible flow. *Siam Journal on Scientific and Statistical Computing*, 7:870–891, 1986.