The Optimal Error Estimate of the Fully Discrete Locally Stabilized Finite Volume Method for the Non-Stationary Navier-Stokes Problem

This paper proves the optimal estimations of a low-order spatial-temporal fully discrete method for the non-stationary Navier-Stokes Problem. In this paper, the semi-implicit scheme based on Euler method is adopted for time discretization, while the special finite volume scheme is adopted for space discretization. Specifically, the spatial discretization adopts the traditional triangle P1−P0 trial function pair, combined with macro element form to ensure local stability. The theoretical analysis results show that under certain conditions, the full discretization proposed here has the characteristics of local stability, and we can indeed obtain the optimal theoretic and numerical order error estimation of velocity and pressure. This helps to enrich the corresponding theoretical results.


Introduction
Recently, due to the characteristics of simple implementation, the finite volume method has been widely used in many scientific research and engineering fields. It has obtained many ideal numerical simulation and calculation results, and often used to solve the complex engineering calculation problems well. Nevertheless, compared with a wide range of application scenarios, its theoretical analysis, such as stability and convergence analysis, is far behind, which inevitably shadow benefits of the finite volume method, so it needs to be studied continuously. Among them, the theoretical analysis of Navier Stokes process is one of the important field.
Finite volume method is an effective method to solve differential equations. In the past several decades, the calculation methods used to solve Navier-Stokes problems have developed rapidly. The results are richer and richer, but to our dismay, the theoretical analysis of the algorithm is still insufficient [1,2]. As we all know, based on the current advanced computing equipment, simple numerical methods are easy to distribute and suitable for large-scale computing, which makes them the hope of solving complex problems, such as incompressible Navier-Stokes equations. Among them, using simple coordinated low-order elements and local stabilization method is a good choice [3,4].
However, it is well known that low order coordinated finite volume function pairs, such as P 1 − P 0 , are unstable for numerical solution of Navier-Stokes equations. The common way to overcome this shortage is to use local stabilization technique, that is, add "macro element condition" to improve the stability of the algorithm. This kind of low order method has been widely analyzed and applied, and has been proved to be effective in practice [5][6][7][8][9]. The basic idea of this method was first proposed by Boland and Nicolaides, and has been vigorously developed since then. The recent work of Wen [10] and Li [11] paves the way for the numerical analysis of Stokes and Navier Stokes problems.
In addition, He [12] and Li [13] have given the locally stable finite volume method for partial spatial discretization of Navier-Stokes problems, which has a good effect except some fully discrete results.
This paper continues to analyze the convergence results of using FVM to solve twodimensional time-dependent Navier-Stokes equations, so as to enrich the relevant theories. Here we will still focus on P 1 − P 0 pairs. For this purpose, let's assume T h is the uniform and regular triangulation ofΩ. It should be reminded that the finite element space here does not have the inf-sup condition for (x h , m h ), so a similar skill in the paper [11][12][13] is required. Because these papers mainly discuss the spatial discrete case, this paper studies a approximation based on time-discretization is Euler semi-implicit and space-discretization is P 1 − P 0 Locally stable FVM.
For brevity, this paper assumes T h is a regular triangulation that satisfies the general regular condition [14,15]. Let the mesh size is h and the time step is 0 < ∆t < 1, the theoretical results of the optimal order convergence of the fully discrete FVM based on the low order coordinated finite element local stabilization are as follows: For such a finite volume solution (u n h , p n h ) obtained by the fully discrete locally stabilized FVM, we derive in this paper the following order error estimates: where N = T ∆t , t n = n∆t ∈ (0, T], τ(t) = min{1, t}. The rest of the paper is organized as follows. Section 2 introduces some basic concepts and function definitions related to Navier-Stokes problem and the locally stable FVM. Some basic results are prepared in Section 3. Section 4 mainly analyzes the error estimation of time semi-discretization based on Euler semi-implicit scheme. Section 5 proves the error estimations of time and spatial fully discretization. Section 6 contains some numerical results and a summary of the article is included in Section 7.

Foundation of Finite Volume Method for the Navier-Stokes Problem
This article consider the follow non-stationary Navier-Stokes equations where Ω be a bounded domain in R 2 assumed to have a common continuous boundary ∂Ω (stronger than Lipschitz continuity) [14,16]; u = u(x, t) = (u 1 (x, t), u 2 (x, t)) are the velocity vector and p = p(x, t) is the pressure, f = f (x, t) is the body force, u 0 (x) is the initial velocity, and ν > 0 is the viscosity. For the convenience of analysis of problem (3), we introduce the following common abbreviations where V is the closed subset of X and H is denote the closed subset of Y, respectively. The spaces L 2 (Ω) m , m = 1, 2, 4, are endowed with the common norm denoted by (·, ·) and · 0 , respectively. The norms of the Hilbert space H 1 0 (Ω) and X are ((u, v)) = (∇u, ∇v), u = ((u, u)) 1/2 .
For more information about the above marks, we refer the reader to [14,16,17]. We also need to denote A = −∆ as the Laplace operator andĀ = −P∆ as the Stokes operator, where P is the L 2 -orthogonal projection of Y onto H.
It is known [17] that where D(A) = H 2 (Ω) 2 ∩ X, γ 0 is positive constant depending only on Ω. C, like the quantities C i , i = 1, 2, · · · , appear subsequently, is a positive constant depending on Ω. This paper uses following kind of continuous bilinear forms a(·, ·) and d(·, ·) on X × X and X × M, respectively, and a trilinear form on With the above notations, the general variational formulation of problem (3) for all (v, q) ∈ (X, M).
In order to get the fully dispersed error estimates, we need the following smoothness.
Theorem 1. Assume some continuity of u 0 and f (x) are valid [14,17]. Then problem (6) admits a unique solution (u, p) satisfying the following estimates: for all t ∈ [0, T], where κ is a positive constant. Now, we consider the fully discrete locally stabilized FVM for two-dimensional timedependent incompressible Navier-Stokes Equation (3). For convenience, let τ h = τ h (Ω) a partitioning ofΩ into triangles satisfied the regular in the usual sense (see [5,14]). (X h , M h ) are the corresponding finite element subspace of (X, M). Γ h is the set of all interelement boundaries.
In order to define the finite volume method for Equation (3), We need introduce a popular configuration of dual partition T * h for T h : the interior point p i is chosen to the barycenter of element K i ∈ T h , and the midpoint m ij on side of v i v j . See Figure 1. This type of dual partition is locally regular if T h is locally regular. Corresponding to Hilbert space V, H, we define the following finite element velocity subspaces and the pressure subspace The L 2 dimensional reduction P h : X→X h and J h : M→M h are defined as follows: The finite volume forms of velocity a(·, ·) on X h × X h is, where n is the unit outnormal vector. The finite volume form d(·, ·) of pressure on X h × M h is defined as To facilitate the analysis, we need the following two trilinear forms.
The last time difference part is The finite volume form of the right side is For the convenience of reading, we introduce the following generalized form In this paper, the norms are defined as following: where S Q is the area of P i P j P k (see Figure 2). To describe the locally stabilized formulation of the non-stationary Navier-Stokes problem, we use the classic not overlap macroelement partitioning Λ h [18]. For every macroelement K in Λ h , the set of interelement(small finite element) edges is denoted by Γ K , and the length of an edge e ∈ Γ K is denoted by h e .
With the above definitions, a locally stabilized formulation of the non-stationary Navier-Stokes problem (3) can be stated as follows.

Definition 1.
Locally stabilized finite volume formulation for non-stationary Navier-Stokes: [p] e [q] e ds, for all p, q in the algebraic sum H 1 (Ω) + M h , and [·] e is the jump operator across e ∈ Γ K and β > 0 is the local stabilization parameter.
In order get the regularity of the above definitions, we need the following stability results [5,8,9,12].

Technical Preliminaries
The main task of this section is to prepare many basic estimates which will help the error analyses for the finite volume solution (u h , p h ).
Since the bilinear forms Moreover, we also need the discrete gradient operators: Firstly, we have the following classical properties (see [8,19]) where For the P 1 − P 0 triangular element, It follows from (4), (5) and (13) that [20,21] As for the trilinear forms b and b,we can deduce the following results.
Proof. Since the following discrete analogue of the Sobolev inequality holds [17] combining the above formula with (14) and (16), we can get the first formula in (15).
To prove the others in (15), we need the follow discrete results [17,22], namely for any For Due to Theorem 2, we know that ( R h , Q h ) is well defined and have the properties [12]: Beside, we need the specific result in He et al. [12] .
Since our error analysis for the time discretization depends heavily on these regularity estimates, we then provide some smoothness estimates of (u h , p h ). The main idea is similar to the work in He et al. [9,23].
Proof. Note from (13) and (14) we can deduce the following estimates: and from Theorems 1 and 3, we have for all t ∈ [0, T]. Then using the similar method in [9], we can get these estimates (20).
Finally, in order to get the upper bounders of velocity and pressure in the time related case, we state the classical Gronwall lemma used in [24]. Lemma 2. Let C 0 and a k , b k , c k , d k , for integers k ≥ 0, be nonnegative numbers such that a n + ∆t Then, a n + ∆t The following is dual Gronwall lemma.

Lemma 3.
Given integer m > 0 and let C and a k , b k , c k , d k , for integers 0 ≤ k ≤ m, be nonnegative numbers such that a n + ∆t Then, a n + ∆t

Error Estimates for Semi-Discrete for Time Depended Navier-Stokes Equations
In this section we consider the time discretization of the locally stabilized and get some useful estimates. Let T time to stop calculation and N the time corresponding step.
For the first part, we need to analysis the errors of finite element Original case. It's well known that the common Euler semi-implicit scheme applied to the spatially discrete problem (9) can be described as: , we integrate and differentiate (9), respectively, to get for all (v h , q h ) ∈ (X h , M h ). Subtracting (25) from (26) and using (27) and the relation: Under the assumptions of Theorem 3 the error E n satisfies the following bounds: Proof. By using (14), (15) and (29), we derive Applying Theorem 4 in (31), we obtain Utilizing Theorem 4, if we sum (32) from n = 1 to n = m, we can derive the first inequality in (30) directly. Next, we deduce from (25) and (29) that Because We can deduce from (33) and Theorem 4 to get for i = 0, 1. Summing (35) from n = 1 to n = m , we derive the second inequality in (30).
As for the last one, Deriving from (14) and (29), we have Hence, Formula (36) and Theorem 4 imply that Similarly, summing (37) from n = 1 to n = m and using Theorem 4 deduces for i = 0,1,2, which yields the third one in (30). Now, let's discuss the second part: the error of time discrete duality argument corresponding to (25). Firstly, the dual problem corresponding to (25) usually describes as: find where Φ m h = 0. For the dual part, we also need the existence, uniqueness and regularity of problem (38), so we introduce the following results: With the similar method in He [9], we have the following two lemmas.

Lemma 5. Assume that the assumptions of Theorem 3 is valid and ∆t satisfies
Then,

Error Analysis for Time and Spacial Discrete
In this section we proof the upper bounds for the error (e n h , µ n h ) = (u h (t n ) − u n h ,p h (t n ) − p n h ) in L 2 and H 1 norms, and deduce the last optimal order estimates. Firstly, we have Proof. Setting v h = e n h , q h = µ n h in (28) and using (14) we obtain From (14) and (15) Summing (44) from n = 1 to n = m, we obtain Applying Lemmas 3 and 6 in (45) and by (14) yields (43).
Proof. We derive from (28) that (47), we obtain From (5) and Lemma 1, we get that Combining (48) with the above estimates and using Lemma 7 yields Summing (48) from n = 1 to n = m and applying Lemmas 5 and 7, we obtain Combining the above estimate with (14) which yields (46).

Lemma 9.
Under the assumptions of Lemma 7, we have that for all 1 ≤ m ≤ N.

Theorem 5.
Under the assumptions of Lemma 7, the error (u(t n ) − u n h , p(t n ) − p n h ) satisfies the following optimal bounds: for all t m ∈ (0, T].

Single Numerical Example
In this section, some numerical results are computed to test the rationality of the theoretical analysis. Because it is difficult to obtain the analytical solution of the general problem governed by the Navier-Stokes equation, we show the relevant numerical results through an example with analytical solution for simplicity. So we consider the following model problem in the unit square area [0, 1] × [0, 1]. Here this example might as well takes ν as 0.01. Only the velocities and pressure are given here. The right term f of the equations can be obtained by bringing the relationship between p(x) and u(x, y) = (u 1 (x, y), u 2 (x, y)) into the NS equations, and the initial values of u 1 (x, y), u 2 (x, y) and p(x, y) can be obtained by bringing t = 0 into the calculation. Now, consider a unit square domain with an exact solution given by u(x, y) = (u 1 (x, y), u 2 (x, y)), p(x, y) = e −16π 2 ν 10(2x − 1)(2y − 1). (79) f is determined by (3). It can be verified that such u 1 (x, y), u 2 (x, y) satisfy the non divergence condition. For simplicity, we can record the time and spatial discretization of the problem as follows: where the matrices in (80) correspond to the differential operators: A − diag(v∆ + 1/∆t), N(u n ) ∼ u n · ∇, B ∼ ∇, C ∼ J h div, C ∼ C h (·, ·), and I is the identity matrix. The right-hand side f n+1 contains the source term.
To make the next iterations less complex, here are a few new notations. Let W = u n+1 , q = p n+1 , then we can further record the above equation as follows: Besides, in order to improve the calculation efficiency, we can generally adopt Newton iterative method to solve the above nonlinear problems. The typical calculation steps are as follows: It is worth noting that since Newton iterative method requires high initial values, we need to use the following Picard method to obtain the initial values of Newton iterative method: In this way, we can finally transform the time-dependent Navier Stokes problem into a large-scale linear system of equations and solve it through such links: Euler time discretization → finite volume space discretization → Newton iterative transformation → Picard format transformation → large linear equations → solving equations.
Based on the above description, we can get some simulation results (β = 6). The following Figure 3 is the result of one iteration based on the initial value. The time step here is ∆t = 0.001, and the spatial grid is divided into two congruent triangles on the basis of 100 × 100 equidistant rectangular grid. We only pay attention to that only 1/5 of the data in both X and Y directions are selected. The following Figure 4 is the result after 1000 iterations with time step ∆t = 0.001. Here, the reference value of streamline is the same as that when t = 0.001. Compared Figure 2 with Figure 3, it can be seen from Figure 3 that the streamline and flow field are weakened accordingly. It is easy to understand that since the interpretation constructed here is decaying, both velocity and pressure show a decaying trend. This also shows that our numerical method maintains strong robustness, and the calculation results are more intuitive.
The following table shows the L 2 error order of velocity and pressure after 1000 iterations with 20 × 20, 40 × 40, 80 × 80 grids selected under the condition of ∆t = 0.001.
From the above Table 1, we can see the optimal L 2 convergence rate, almost 2 for velocities and 1 for pressure are really obtained, which confirm the numerical analysis above. The following Figure 5 shows the error curve calculated based on 100 × 100 spatial grid with different time steps: ∆t = 0.01, 0.005, 0.0025, 0.00125 ('dt' in Figure 5 is ∆t). It can be seen from here that the initial error tends to increase, but the error also decreases with the decrease of flow field energy, which shows that the time iteration is stable. The following Table 2 shows the error ratio of different time steps: dt = 0.01, 0.005, 0.0025, 0.00125: The error ratio curve in the Figure 6 below shows the convergence of one section of time and is relatively stable.   Figure 6 tell us the optimal convergence rate of time is 1 which is consistent with the theoretical analysis.
Due to time constraints, our numerical results only show these. It is also worth noting that if the solution does not decay but increases, and if the growth rate is fast, the error of numerical results will increase with the increase of numerical calculation time. At that time, the method may not converge or inefficient. This requires a little attention in specific applications.

Conclusions
After detailed theoretical analysis, this article finally proves that if we use the finite volume method based on P 1 − P 0 element to approximate the non-stationary Navier-Stokes equation, we can achieve the follow optimal numerical error estimation: The optimal error estimate (83) shows that the time discretization of Euler method is 1 order and the space discretization is 2 order in this space-time full discretization finite volume method, which is consistent with the theoretical optimal order error estimation of P 1 − P 0 element in solving differential equations. Although the proof process is challenging and cumbersome, the optimal result is also obvious and certain. However, this work is still helpful to reveal some special aspects of the finite volume method which is different from finite element and other methods in solving complex differential equations. Therefore, it is helpful to improve the corresponding numerical analysis theory.
In addition, with the continuous research and exploration of using neural network to solve differential equations recently [25][26][27], although very gratifying numerical results have been obtained, but the effectiveness and convergence of neural network-based methods for solving differential equations are still unclear, and there is a lack of theory. However, we know gradually that the key point of neural network in the calculation of differential equations is to use low-order functions for numerical approximation, which is similar to the basic principle of using low-order continuous functions to discretize and approximate differential equations, except that the former is global approximation while the latter is piecewise approximation. Therefore, improving and enriching the low-order function approximation theory also provides some important reference in further understanding the neural network in solving differential equations efficiently, It is worth studying.