{'fraction_non_alphanumeric': 0.06465001495662578, 'fraction_numerical': 0.024977565061322166, 'mean_word_length': 3.6088230225745304, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 11, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 13, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Convergence of the fully discrete incremental projection scheme for incompressible flows

The present paper addresses the convergence of a first order in time incremental projection scheme for the time-dependent incompressible Navier-Stokes equations to a weak solution, without any assumption of existence or regularity assumptions on the exact solution. We prove the convergence of the approximate solutions obtained by the semi-discrete scheme and a fully discrete scheme using a staggered finite volume scheme on non uniform rectangular meshes. Some first a priori estimates on the approximate solutions yield the existence. Compactness arguments, relying on these estimates, together with some estimates on the translates of the discrete time derivatives, are then developed to obtain convergence (up to the extraction of a subsequence), when the time step tends to zero in the semi-discrete scheme and when the space and time steps tend to zero in the fully discrete scheme; the approximate solutions are thus shown to converge to a limit function which is then shown to be a weak solution to the continuous problem by passing to the limit in these schemes.


Introduction
The incompressible Navier-Stokes equations for a homogeneous fluid read: where the density and the viscosity are set to one for the sake of simplicity, and where (2) T > 0, and Ω is a connected, open and bounded subset of R 3 , with a Lipschitz boundary ∂Ω.
Note that we only consider the three dimensional setting in this work, but the analysis may be carried out in a similar (and often somewhat simpler) manner in the one or two dimensional setting. The variables u and p are respectively the velocity and the pressure in the flow, and Eqns. (1a) and (1b) respectively enforce the momentum conservation and the mass conservation and incompressibility of the flow. This system is supplemented with the boundary condition  The function u 0 is the initial datum for the velocity and the function f is the source term. Throughout the paper, we shall assume that (5) f ∈ L 2 ((0, T ) × Ω) 3 and u 0 ∈ E(Ω), where E(Ω) is the subset of H 1 0 (Ω) 3 of divergence-free functions, defined by E(Ω) = {u ∈ H 1 0 (Ω) 3 such that divu = 0}. Note that in fact, the initial condition is assumed to be in E(Ω) for the sake of simplicity. It could be considered in L 2 (Ω) 3 only, see Remark 2.3.
Let us define the weak solutions of Problem (1)- (4) in the sense of Leray [17].
. Under the assumptions (2) and (5), a function u ∈ L 2 (0, T ; E(Ω)) ∩ L ∞ (0, T ; L 2 (Ω) 3 ) is a weak solution of the problem (1)-(4) if in Ω × (0, T ) . The first projection method to solve the system (1) was designed over 50 years ago, and is known as the Chorin-Temam algorithm [4,21,22]. It consists in a prediction step based on a linearized momentum equation without the pressure gradient, and a pressure correction step that enforces the divergence-free constraint. This method and its variants are now often referred to (following [13]) as non incremental projection schemes, in opposition to the incremental projection schemes that were obtained by adding the old pressure gradient in the prediction step (see [12] for a first-order time scheme and [23] of a second order time scheme). These latter schemes are indeed incremental in the sense that the correction step may now be seen as solving an equation on the time increment of the pressure. They seem to be much more efficient from a computational point of view [13] and have been the object of several error analysis, under some regularity assumptions on the solution of the continuous problem, in the semi discrete setting, see [13] and references therein.
The non incremental schemes have been the object of some analyses in the fully discrete setting. In [1] some error estimates are derived for a non incremental scheme with a discretization by the finite element, under some regularity assumptions on the exact solution, In a recent paper the approximate solutions of a fully discrete non incremental scheme with a uniform staggered discretization [15] are shown to converge to a weak solution (and so without any regularity assumption on the solution of (1)) under the condition that h ≤ δt 3−α where h and δt are respectively the mesh size and the time step and with 0 < α ≤ 2.
However, to our knowledge, up to now, no proof of convergence exists for the fully discrete incremental projection schemes, even though they are the most used in practice. The purpose of the present work is therefore to fill this gap and to show the convergence of the incremental projection method with a discretization by a staggered finite volume scheme based on a (non uniform) MAC grid, without any regularity assumption on the exact solution.
The Marker-And-Cell (MAC) scheme, introduced in the middle of the sixties (see [14]), is one of the most popular methods (see e.g. [18] and [24]) for the approximation of the Navier -Stokes equations in the engineering framework, because of its simplicity, its efficiency and its remarkable mathematical properties. Although originally presented as a finite difference scheme on uniform meshes, the MAC scheme is in fact a finite volume scheme and as such can be used on non uniform meshes. The convergence analysis of the staggered finite volume scheme on the MAC mesh using a fully implicit time scheme may be found in [11], and we shall use several of the tools developed therein. We also refer to this latter paper for some more references on studies of the MAC scheme.
The paper is organized as follows. Section 2 deals with the convergence analysis for the semi-discrete projection algorithm. The fully discrete scheme is analysed in Section 3; we only give the main ingredients of the staggered space discretization that we use, and which is often referred to as the MAC scheme. To avoid a lengthy description, the precise definitions of the now classical discrete MAC operators are to be found in [11].
Before starting the analysis of the semi-discrete and fully discrete schemes, we wish to recall, for the sake of clarity, that: • In a Banach space E equipped with a norm · E , a sequence (u n ) n∈N ⊂ E is said to converge to u ∈ E if u n − u E → 0 as n → +∞, while it is said to weakly converge to u ∈ E if for any continuous linear form T ∈ E ′ , one has T (u n ) → T (u) as n → +∞. • A sequence (T n ) n∈N ⊂ E ′ is said to ⋆-weakly converge to T ∈ E ′ if for any u ∈ E, one has T n (u) → T (u) as n → +∞. • If E = L p (Ω), where 1 ≤ p < +∞ and Ω is an open set of R 3 , the space E ′ is identified to L q (Ω), q = p/(p − 1). • For T > 0 and E = L 1 ((0, T ), L 2 (Ω)), the space L ∞ ((0, T ), L 2 (Ω)) is identified with E ′ .
In the appendix, we give some useful technical lemmas.

Analysis of the time semi-discrete incremental projection scheme
We consider a partition of the time interval [0, T ], which we suppose uniform to alleviate the notations, so that the assumptions read: N ≥ 1, δt N = T N , t n N = n δt N for n ∈ 0, N .

2.1.
The time semi-discrete scheme. Under the assumptions (7), the usual first order time semi-discrete incremental projection scheme (see [20]) reads: Initialization: Solve for 0 ≤ n ≤ N − 1 : where n stands for the outward normal unit vector to the boundary ∂Ω and f n+1 Let us briefly account for the existence of a solution at each step of this algorithm.
Prediction step -A weak form of Eqns. (8b)-(8c) reads The existence of the predicted velocity is then a consequence of Lemma A.1.
The existence of (u n+1 N , p n+1 N ) ∈ V (Ω) × H 1 (Ω) satisfying (10) is a consequence of the decomposition result of Lemma A.2 given in the appendix. Indeed, this correction step is the decomposition stated in Lemma A.2 applied to the predicted velocityũ n+1 . Note that p n+1 N is uniquely defined thanks to (8f).
Remark 2.1. Summing (8b) at step n and (8d) at step n − 1, we obtain for n ∈ 1, N − 1 We may thus state the following existence result and define the approximate solutions obtained by the projection scheme (8).
Definition 2.1 (Approximate solutions, semi-discrete case). Under the assumptions (2),(5) and (7), there exists (ũ n N , u n N , p n N ) n∈ 1,N ⊂ H 1 0 (Ω) 3 × V (Ω) ×H 1 (Ω) satisfying (8). We then define the functions u N : where (ũ n N ) n∈ 1,N and (u n N ) n∈ 1,N are a solution to (8), where 1 A denotes the indicator function of a given set A. (1) is imposed on the weak solution through the functional space H 1 0 (Ω) 3 . Note that this condition is only imposed on the predicted velocity in the algorithm (8). Indeed, the corrected velocity does not satisfy the full Dirichlet condition (3) but only the no slip condition imposed by (8e). The compactness of the sequence of predicted velocitiesũ together with the convergence of u −ũ towards zero in L 2 as the time step tends to zero will be the mean to prove that the Dirichlet boundary condition is finally satisfied on the limit of the numerical approximations. Note also that there is no need for a boundary condition on the pressure in the correction step. In fact, it can be inferred from the correction step (10) that the incremental pressure ψ n+1 = p n+1 − p n satisfies a Poisson equation on Ω with a Neumann boundary condition on the boundary, but this is a redundant information that does not need to be implemented. We refer to [19] for an interesting discussion on these boundary conditions. Remark 2.3 (On the initial condition). In fact, the existence of a solution (see Lemma A.1) only requires the initial velocity u 0 N to be in V (Ω), so that we could relax the assumption on the initial condition u 0 ∈ E(Ω) to u 0 ∈ L 2 (Ω) 3 and take u 0 = P V (Ω) u 0 as the orthogonal projection of u 0 onto the closed subspace V (Ω) of L 2 (Ω) 3 , also known as the Leray projection. In this case, u 0 N can be computed as u 0 = u 0 − ∇ψ where ψ ∈ H 1 (Ω) is a solution (unique, up to a constant) of the following problem (see Lemma A.2) ψ ∈ H 1 (Ω), Ω ∇ψ · ∇ϕ dx = Ω u 0 · ∇ϕ dx, for any ϕ ∈ H 1 (Ω).

Remark 2.2 (On the boundary conditions). The original homogeneous Dirichlet boundary conditions (3) of the strong formulation
Theorem 2.1 (Convergence of the semi-discrete in time projection algorithm). Under the assumptions (2) and (5), consider for N ≥ 1, the time discretization defined by (7), and the approximate solutions u N andũ N of the projection algorithm (8) as given in Definition 2.1. Then there existsū ∈ L 2 (0, T ; E(Ω)) ∩ L ∞ (0, T ; L 2 (Ω) 3 ) such that up to a subsequence, • the sequence (ũ N ) N ≥1 converges toū in L 2 (0, T ; L 2 (Ω) 3 ) and weakly in . Moreover the functionū is a weak solution to (1) in the sense of Definition 1.1.
Proof. Here are the main steps of the proof; each step is detailed in one of the following paragraphs.
• Step 2: compactness and convergence in L 2 (detailed in section 2.3) This is the tricky part of the proof. Since the sequence (ũ N ) N ≥1 converges weakly in L 2 (0, T ; H 1 0 (Ω) 3 ), some estimate on the discrete time derivative would be sufficient to obtain the convergence in L 2 (0, T ; H 1 0 (Ω) 3 ) by a Kolmogorovlike theorem. A difficulty to obtain this estimate arises from the presence of the pressure gradient in Equation (8b), which needs to be "killed" by multiplying this latter equation by a divergence-free function. This function ϕ should also be regular enough so that the nonlinear divergence term makes sense: hence we choose ϕ ∈ L 2 (0, T ; W 1,3 0 (Ω) 3 ) such that divϕ = 0, and define the following semi-norm on (L 2 (Ω)) 3 : with W (Ω) = {ϕ ∈ W 1,3 0 (Ω) 3 : Estimates on the L 2 (|·| * ,1 ) semi-norm of the time translates of the predicted velocity u N are then obtained from the semi-discrete momentum equation (8b): see Lemma 2.3. Note that this is only an intermediate result; indeed, in order to gain compactness, we need an estimate on the time translates of the predicted velocity in the L 2 (L 2 ) norm. The idea is then to first introduce the following semi-norm on (L 2 (Ω)) 3 .
where V (Ω) is the space of L 2 divergence-free functions defined by (11). Note that where P V (Ω) is the orthogonal projection operator onto V (Ω). Then, thanks to a Lions-like lemma (Lemma 2.4 below), we get that for any ε > 0, there exists C ε ∈ R + such that By (15), we have an L 2 (0, T ; H 1 0 (Ω) 3 ) bound on the predicted velocities; we have also seen that the time translates of u N for the L 2 (| · | * ,1 ) semi-norm tend to 0 as N → +∞ (Lemma 2.3 below). Therefore, by (20), the time translates of u N for the L 2 (| · | * ,0 ) semi-norm also tend to 0 as N → +∞. In order to show that the L 2 (L 2 ) norm of the time translates ofũ N tend to 0, we remark that if v ∈ V (Ω), then |v| * ,0 = ||v|| L 2 (Ω) and conclude thanks to (16), see Lemma 2.5).

• Step 3: convergence towards the weak solution (detailed in section 2.4)
Owing to a Kolmogorov-type theorem (see e.g. [9,Corollary 4.41]), the estimates of steps 1 and 2 yield that there exist subsequences, still denoted (u N ) N ≥1 and (ũ N ) N ≥1 , that converge toū in L 2 (0, T ; L 2 (Ω) 3 ). In section 2.4, we pass to the limit in the scheme to obtain thatū satisfies (6); thereforeū is a weak solution to (1) in the sense of Definition 1.1.
Remark 2.4 (Uniqueness and convergence of the whole sequence). In the case where uniqueness of the solution is known, then the whole sequence converges ; this is for instance the case in the 2D case [17], see e.g. [3, Chapter 5, Section 1.3] for more on this subject.   (7), the functions u N andũ N defined by (13) satisfy (15) and (16), with C 1 depending only on |Ω|, u 0 L 2 (Ω) 3 and f L 2 (Ω) 3 .
2.3. Proof of step 2: compactness and L 2 congergence. Following Step 2 of the sketch of proof of Theorem 2.1, we start by the following lemma.

Passing to the limit in this inequality yields that
Owing to Lemma A.3, this in turn implies that there exists ξ ∈ H 1 (Ω) such that P V (Ω) w = ∇ξ. Using the fact that P V (Ω) w ∈ V (Ω) we have and is therefore relatively compact in L 2 (0, T ; L 2 (Ω) 3 ).
Then the functionū is a weak solution to (1) in the sense of Definition 1.1.
The regularity of f and ϕ implies that: Multiplying (12) by δt N ϕ n N , integrating over Ω and summing over n ∈ 1, N − 1 yields Using the fact that ϕ N N = 0 in Ω the first term of the left hand side reads By the triangle inequality, Since the sequence (ũ N ) N ≥1 converges toū in L 2 (0, T ; L 2 (Ω) 3 ), we obtain The second term in the left hand-side reads T δtN Ωũ The third term in the left hand-side may be written The weak convergence of the sequence (∇ũ N ) N ≥1 in L 2 (0, T ; L 2 (Ω) 3 ) and the con- The right hand-side satisfies The convergence of the sequence (f N ) N ≥1 in L 2 (0, T ; L 2 (Ω) 3 ) and the convergence of the sequence Using (28)-(31) and passing to the limit in (27) gives the expected result.

Analysis of the fully discrete projection scheme
Our purpose is now to adapt the proof of convergence of the semi-discrete case to the fully discrete case. We choose as an example of space discretization a staggered discretization on a (possibly non uniform) rectangular grid of R 3 . The resulting scheme, often referred to as a MAC scheme, was analysed in [11] for an timeimplicit scheme. The idea here is to prove its convergence for the incremental projection scheme. We consider the following assumptions on Ω and on the timespace discretization, indexed by N (in the convergence analysis, the time and space steps will tend to 0 as N tends to +∞).
T > 0, Ω is an open rectangular subset of R 3 , with boundary faces that are orthogonal to one of the vectors of the canonical basis of (32a) Note that at this step, we are only considering one time step δt N = T N and one discretization mesh D N , which is also indexed by N . This might seem strange, but it is in view of the convergence analysis for which a sequence (D N , δt N ) N ≥1 will be considered, with h N , δt N → 0 as N → +∞.
The regularity of the mesh is measured by the following parameter: with | · | the Lebesgue measure (this notation is used in the following for either the R 3 or the R 2 measure). We refer to [11] for the precise definition ot the discrete spaces and operators. The approximate pressure belongs to the set L N (Ω) of functions that are piecewise constant on the so called primal cells K of the (primal) mesh M N : where x K denotes the mass center of K and x K,∂Ω the orthogonal projection of x K on ∂Ω We thus define three dual meshes of Ω.
3.1. The fully discrete scheme. The space discretization of the time-discrete scheme (8) reads: Initialization: In this algorithm, the terms C N (ũ n+1 N )u n N , -∆ Nũ n+1 N , ∇ N p n N and div N u n+1 N are the MAC discretization of the terms div(ũ n+1 N ⊗ u n N ), ∆ũ n+1 N , ∇p n N and divu n+1 N in the algorithm (8) and are defined in [11,Section 2]. In (34c), the vector function f n+1 N is defined by its components (f n+1 N,i , i = 1, 2, 3) where f n+1 N,i is the piecewise constant function from Ω × (0, T ) to R 3 defined by Let us briefly account for the existence of a solution at each step of this algorithm.
First remark that the discrete no slip boundary condition (34c) and (34f) are equivalent to requiring that the i-th component u i of the approximate predicted and corrected velocities belongs to the space H  Note that since u 0 ∈ E(Ω), we also have u 0 N ∈ E N (Ω). Prediction step -The existence of the predicted velocity follows from Lemma B.1.

Correction step -A weak form of the correction step (34d) which computes a divergence-free velocity and an associated pressure reads
Note that if u n+1 N satisfies (35), then Ω u n+1 N · ∇ N q dx = 0 for any q ∈ L N (Ω), so that u n+1 ∈ E N (Ω). The existence of (u n+1 , p n+1 ) ∈ E N (Ω) × L N (Ω) satisfying (35) is a consequence of the decomposition result of Lemma B.2 (given in the appendix).
We may then define the approximate solutions as follows.
For a given N ≥ 1 and the associated (uniform) time discretization  (10) in [11,Section 2]. As in the semi-discrete case, it is not imposed in the correction step (34g)-(34d). Note also that, as in the semi-discrete case, there is no need for a boundary condition on the pressure in the correction step. In fact, it can be inferred from the correction step that the incremental pressure ψ n+1 = p n+1 − p n satisfies a discrete Poisson equation on Ω with a Neumann boundary condition on the boundary.

Remark 3.2 (On the initial condition).
If the initial condition u 0 ∈ E(Ω) is relaxed to u 0 ∈ L 2 (Ω) 3 as in Remark 2.3, the discrete initial condition should be taken as the orthogonal projection onto E N (Ω) of the function u 0 defined by (34a). Remark 3.3. Summing (34b) and (34d), we get the discrete equivalent of (12): Let us now state the convergence of the algorithm (34) as the time step δt N and the mesh step h N tend to 0 (or N = T δtN → +∞); the proof of this result is the object of the following sections. Proof. We give here the main steps of the proof, which follows that of the semidiscrete case; these steps are detailed in the following paragraphs.
• Step 1: first estimates and weak convergence (detailed in Section 3.2). Let us define, for q ∈ N * , a discrete W 1,q 0 (Ω) 3 -norm for the discrete velocity fields.
There remains to show thatū is a weak solution in the sense of Definition 1.1 and in particular thatū satisfies (6). The weak convergence is not sufficient to pass to the limit in the scheme, because of the nonlinear convection term, so that we first need to get some compactness on one of the subsequences (ũ N ) N ∈N or (u N ) N ∈N (since, by (16), their difference tends to 0 in the L 2 norm).
• Step 2: compactness and convergence in L 2 (detailed in section 3.3) We adapt Step 2 of the convergence proof of the semi-discrete case. Using the bound (40) on the sequence (ũ N ) N ≥1 , some estimate on the discrete time derivative would be sufficient to obtain the convergence in L 2 (0, T ; H 1 0 (Ω) 3 ) by a Kolmogorov-like theorem. As in the semi-discrete case, a difficulty arises from the presence of the (discrete) pressure gradient in Equation (12); we get rid of it by multiplying this latter equation by a discrete divergence-free function, chosen as the interpolate of a regular function ϕ ∈ L 2 (0, T ; (W 1,3 0 (Ω)) 3 ) such that divϕ = 0. Let us then define the discrete equivalent of the semi-norm (17) on H N,0 (Ω) by: Estimates on the L 2 (| · | * ,1,N ) semi-norm of the time translates of the predicted velocity u N are then obtained from the discrete momentum equation (12): see Lemma 3.3. Again, this is only an intermediate result since we seek an estimate on the time translates of the predicted velocity in the L 2 (L 2 ) norm. So next, as in the semi-discrete case, we introduce the discrete equivalent of the semi-norm | · | * ,0,N .
(46) |w| * ,0,N = P EN (Ω) w L 2 (Ω) 3 , for any w ∈ H N,0 (Ω), where P EN (Ω) is the orthogonal projection operator onto E N (Ω). Then, thanks to a discrete equivalent of the Lions-like 2.4 lemma (Lemma 3.4 below), we get that for any ε > 0, there exists C ε ∈ R + such that From this latter inequality, using Lemma 3.3 on the time translates of u N for the L 2 (| · | * ,1 ) semi-norm and the bound (40), we get that the time translates of u N for the L 2 (| · | * ,0,N ) semi-norm also tend to 0 as N → +∞.
Remark 3.4 (Uniqueness and convergence of the whole sequence). If the solution of the continuous problem is unique, then again the whole sequence converges.
Proof. By Lemma B.1 with α = 1 δtN , we have for n ∈ 0, N − 1 Squaring the relation (34d), integrating over Ω, multiplying by δtN 2 and owing to (34e) and to the discrete duality property of the MAC scheme [11, Lemma 2.4], we get Summing this latter relation with the previous relation yields for n ∈ 0, N − 1 We then get the relation (48) using the Cauchy-Schwarz inequality and the discrete Poincaré estimate [8, Lemma 9.1] after summing over the time steps.
defined by (22). Using (38), we thus get that where ∇ N is the gradient operator of the velocity defined on each dual rectangular grid, see [11, Section 2]), we get that Let us reproduce at the fully discrete level the computations done for each of these terms in the proof of Lemma 2.3.
Therefore, thanks to the boundedness assumptions on u N and u N , Again invoking the discrete Sobolev inequality, there exists C 7 ∈ R + only depending on |Ω| such that v L 3 (Ω) 3 ≤ C 7 v 1,3,N , for any v ∈ H N,0 (Ω). Consequently, Thanks to the fact that ϕ ∈ E N (Ω) and to the discrete duality property stated in [11,Lemma 2.4], A p (t) = 0.
For v = (v 1 , v 2 , v 3 ) ∈ (L 2 (Ω) 3 , we define P N v as the vector function with piecewise constant components: the i-th component of P N v is constant on each dual cell D σ , σ ∈ E, and equal to the mean value of v i on the face σ. By [11,Lemma 3.7], P N is a Fortin operator in the sense that it preserves the divergence; in particular, Lemma 3.4 (Lions-like, fully discrete version). Consider a rectangular domain Ω of R 3 and a sequence of MAC grids (D N ) N ≥1 of Ω satisfying (32c) such that h N → 0 as N → +∞ and such that the mesh regularity parameter θ N defined by (33) remains bounded. Then, for any ε > 0, there exists C ε > 0 and N ε ≥ 1 depending on ε such that for any N ≥ N ε and for any w ∈ H N,0 (Ω), (47) is satisfied.
By the W 1,q stability of the operator P N stated in [10,Theorem 1], there exists C 8 only depending on |Ω| and on θ N in a nondecreasing way, such that 3 , for any ϕ ∈ W 1,3 0 (Ω) 3 . Passing to the limit in this inequality yields that Ω P V (Ω) w · ϕ dx = 0, for any ϕ ∈ W (Ω).
Proof. Using the fact that (v N ) N ≥1 is bounded in L 2 (Ω) 3 we obtain that the sequence (P EN (Ω) v N ) N ≥1 is bounded in L 2 (Ω) 3 . Hence there exists a subsequence still denoted by (P EN (Ω) v N ) N ≥1 that converges to a functionṽ weakly in L 2 (Ω) 3 . Thanks to the discrete duality property stated in [11,Lemma 2.4], we have, for any for all x ∈ K, K ∈ M N . The discrete gradient ∇ N is consistent in the sense of [11,Lemme 2.3] and therefore there exists C 9 ∈ R + depending only on Ω and on θ N in a nondecreasing way, such that Passing to the limit in the previous identity gives Ωṽ · ∇ϕ dx = 0, for any ϕ ∈ C ∞ c (R 3 ).
Passing to the limit in the previous identity gives Ω v · ϕ dx = Ωṽ · ϕ dx for any ϕ ∈ V (Ω).
We then obtain thatṽ = P V (Ω) v and the sequence (P EN (Ω) v N ) N ≥1 converges to P V (Ω) v weakly in L 2 (Ω) 3 . We can write Using the convergence of the sequence (v N ) N ≥1 to v in L 2 (Ω) 3 and the weak convergence of the sequence ( The weak convergence of the sequence (P EN (Ω) v N ) N ≥1 to P V (Ω) in L 2 (Ω) 3 and convergence of the sequence ( P EN (Ω) v N L 2 (Ω) 3 ) N ≥1 to P V (Ω) v L 2 (Ω) 3 gives the expected result.
Lemma 3.6. Under the assumptions of Theorem 3.1 the sequence (ũ N ) N ≥1 satisfies uniformly with respect to N , and is therefore relatively compact in L 2 (0, T ; L 2 (Ω) 3 ).
Proof. We follow the proof of Lemma 2.5. By the triangle inequality, For any N , A N (τ ) → 0 as τ → 0, but owing to (43), we get that A N (τ ) → 0 as τ → 0, uniformly with respect to N . Let us prove that this is also the case for B N (τ ) → 0. Since u N (t) ∈ E N (Ω) for any t ∈ (0, T ) we have for any t ∈ (0, T − τ ) Now thanks to Lemma 3.4, for any ε > 0, there exists C ε > 0 and N ε ≥ 1 such that for any N ≥ N ε and for any t ∈ (0, T − τ ) In particular for any N ≥ N ε and for any τ ∈ (0, T ) we have Therefore, owing to lemmas 3.2 and 3.3, for any N ≥ N ε and for any τ ∈ (0, T ) Hence, for any N ≥ N ε and for any τ ∈ (0, T ), . Let ζ > 0 be given, and let: • τ 0 > 0 such that for any τ ∈ (0, τ 0 ), 2A N (τ ) ≤ ζ for any N ≥ 1.
Let us now give the decomposition result which was used for the proof of existence of a solution to the correction step (10).
The following lemma gives a characterisation of the gradient which is used in the proof of Lemma 2.4. Its proof is a simple consequence of a result of M. E. Bogovskii [2] and refer to the very clear presentation of [5] for more on this subject.
A consequence of this lemma is the following interesting per se density result. Proof. Equipped with the L 2 (Ω) 3 -norm, the space V (Ω) is a Hilbert space. In order to prove this density result, we prove that, in this Hibert space, V(Ω) ⊥ = {0}.

Appendix B. Some discrete technical lemmas
We assume that Ω is an open rectangular parallelepiped. In addition, we assume that the edges (respectively the faces) of Ω are orthogonal to one vector of the canonical basis of R 3 .
The following lemma is the discrete version of lemma A.2 which was used for the proof of existence of a solution to the correction step (35). Proof. Let ψ be a solution (unique, up to a constant) of ψ ∈ L N (Ω), Ω ∇ N ψ · ∇ N ξ dx = Ω w · ∇ N ξ dx, for any ξ ∈ L N (Ω).